clear % clear all variables
D=10; % diffusion coefficient (um^2/s)
tau=0.1; % step time (s)
Nsteps=1000; % number of steps
Nmolec=1; % number of molecules
stepsize=sqrt(2*D*tau);% rms displacement (um)
t=[0:tau:(Nsteps-1)*tau]; % set array of times (s)
x(1)=0; % initial position (um)
for j=1:Nmolec % simulate for each molecule
    for i=1:(Nsteps-1) % for each molecule take Nsteps-1 steps
        tmp=rand; % generate a random number
        if tmp<0.5 % if rand#<0.5, then step to right
            dx=stepsize; % step increment is positive
        else
            dx=-stepsize; % otherwise step increment is negative
        end
        x(i+1)=x(i)+dx; % update the position
    end
   xfinal(j)=x(Nsteps); % record the final position
   x(1)=0; % reset position to zero before simulating the next molec
end
figure(1) 
plot(t,x) % plot the last molecule's history
axis([0 Nsteps*tau -40 40])% adjust the axes on the plot
xlabel('TIME (s)','FontSize',16)
ylabel('POSITION (um)','FontSize',16)
%figure(2)
%hist(xfinal) % plot the histogram of final positions
%axis([-1000 1000 0 400]) % adjust the axes on the plot
%xlabel('POSITION (um)','FontSize',16)
%ylabel('NUMBER OF MOLECULES','FontSize',16)