Chemistry Reference
In-Depth Information
eout
¼
zeros(1,nstep
þ
1);
x(:,1)
¼
x0;
v(:,1)
¼
v0;
[e,f,ffast,fslow]
¼
feval(fun,x0,0);
eout(1)
¼
e
þ
.5*sum(M.*v0.*v0);
for i
¼
1:nstep
xx
¼
x(:,i);
vv
¼
v(:,i)
þ
(tau*h/2)*fslow./M;
for j
¼
1:tau
vv
¼
vv
þ
(h/2)*ffast./M;
xx
¼
xx
þ
h*vv;
[e,f,ffast,fslow]
¼
feval(fun,xx,1);
vv
¼
vv
þ
(h/2)*ffast./M;
end
x(:,i
þ
1)
¼
xx;
[e,f,ffast,fslow]
¼
feval(fun,xx,0);
v(:,i
þ
1)
¼
vv
þ
(tau*h/2)*fslow./M;
eout(i
þ
1)
¼
e
þ
.5*sum(M.*v(:,i
þ
1).*v(:,i
þ
1));
end
%%%%%%%%%%end of function %%%%%%%%%%%%%%%%%%%%%%
function [t,x,v,eout]
¼
bbk1(fun,h,nstep,x0,v0,M,
gamma,temp);
%[t,x,v,eout]
¼
bbk1(fun,h,nstep,x0,v0,M,gamma,
temp);
% Langevin dynamics integrator
%
%INPUT ARGUMENTS
% fun string containing file name of energy and force
routine
% h timestep
% nstep number of steps
% x0 vector of initial atom positions
% v0 vector of initial atom velocities
% M vector of atomic masses
% gamma Langevin friction parameter
% temp target temperature for Langevin dynamics
%OUTPUT ARGUMENTS
% t vector of time values
% x trajectory of positions
% v trajectory of velocities
% eout energy along computed trajectory
%
%Eric Barth
%Kalamazoo College, 2007
