clear all;
close all;
s = 0.01; % an individual displacement of a molecule as a result of a single collision
n = input('Enter the number of molecular collisions, n = ');
t = 1 : n; % the vector of time instances
Nvect = 100 : 100 : 1000; % the vector of number of molecules in an ensemble
for k = 1 : length(Nvect)
N = Nvect(k); % the given ensemble has N molecules
A = rand(N,n) - 0.5; % A is an N-by-n matrix of random numbers equally distributed in the interval [-0.5, 0.5]
A = s*sign(A);
X = zeros(N,n); % X is an N-by-n matrix of total displacements of molecules
% columns of X correspond to positions of all N molecules at a time instance
for j = 1:n
if ( j ~= 1)
X(:,j) = sum(A(:,1:j)')';
% the sum is taken for each row of A between 1 and j columns
% the result is assigned to the j column of matrix X
% X contains trajectories of all N molecules simultaneously
else
X(:,1) = A(:,1);
% no summations are performed for the first collision
end
end
Xav = sum(X)/N; % Xav is a row vector of n entries which give the mean at different time instances
XX = zeros(N,n);
for j = 1 : n
XX(:,j) = Xav(j)*ones(N,1); % the j-th column of XX has the value of Xav(j)
end
Xav2 = sum((X - XX).^2)/N; % Xav2 is a row vector of n entries which give the variance at different time instances
figure(1);
plot(Xav); % plot the mean as function of time at figure(1)
hold on;
figure(2);
plot(Xav2); % plot the variance as function of time at figure(2)
hold on;
Dtime = Xav2./t; % the value of the coefficient diffusion D depends on the time t, when we compute it
D(k) = sum(Dtime)/n; % this is a naive and rough way to approximate the theoretical value of D
% NOTE: a correct solution of the approximation problem should be done with
% the use of MATLAB function "polyfit"
end
figure(1)
title('Evolution of the mean \mu in time');
xlabel('Time'); ylabel('\mu'); hold off;
figure(2)
title('Evolution of the variance \sigma^2 in time');
xlabel('Time'); ylabel('\sigma^2'); hold off;
figure(3)
plot(Nvect,D,'*r',Nvect,D,'-b');
xlabel('N'); ylabel('D');
title('Coefficient of diffusion versus the number of molecules');