% This is to illustrate diffusion of particles (i.e., Brownian motion)
% in nd dimensions, and exit from a region of radius L.
clear
 set(0,                           ...
   'defaultaxesfontsize', 20,   ...
   'defaultaxeslinewidth', 1.0, ...
   'defaultlinelinewidth', 1.2, ...
   'defaultpatchlinewidth', 0.7);


nd = 2;  % set the dimension

N = 2000;  %number of particles
 D = 1;
dt = .0005; 
L = 1; % radius of the bounding circle
Rsq = L^2;
diff = sqrt(2*D*dt);
 
for j = 1:N %number of particles
   clear X 
  
   x = zeros(nd,1); %the particle starts at the origin
X = x;

t = 0;
absx = 0;
    while(absx<Rsq)
    x = x + diff*randn(nd,1);
    X = [X,x];
     t = t+dt;
     absx = x'*x;
    end
if (j==1)% make a plot showing one of the particles
    th = 2*[0:.01:1]*pi;
    figure(1)
    plot(X(1,:),X(2,:),L*sin(th),L*cos(th))
    axis([-L L -L L])
    xlabel('x','fontsize',20)
    ylabel('y','fontsize',20)
    

end

T(j) = t;
end
mean(T)  % this is the mean exit time
%compare this with the theoretical answer
L^2/(2*nd*D)
figure(2)
hist(T,40)
 xlabel('t','fontsize',20)
  ylabel('First Exit Time Histogram','fontsize',20) 
