View source for Molecular Dynamics Algorithms

Navigation: [[Documentation]] / [[Technical Details]] / [[Molecular Dynamics Algorithms]]
----

=AMBER=
Both sander and pmemd from the [https://ambermd.org AMBER] package use the leap-frog integration algorithm. The algorithm is demonstrated for the SANDER driver:

<ol>
<li>'''Initialization''' (sander.F90)</li>
<pre>
call pmf_sander_init_preinit(mdin,natom,nres,ntb,nstlim,dt,temp0,a,b,c,alpha,beta,gamma)
! topology population
do i=1,nres
    call pmf_sander_set_residue(i,ih(m02-1+i),ix(i02-1+i))
end do
do i=1,natom
    call pmf_sander_set_atom(i,ih(m04-1+i),ih(m06-1+i))
end do
call pmf_sander_finalize_preinit(natom,x(lmass),x(lcrd))
call pmf_sander_cst_init_collisions(ntc,nbonh,ix(iifstwt),ix(iibh),ix(ijbh),x(l50))
call pmf_sander_init(natom,x(lmass),x(lcrd))
</pre>
<li>'''MD Loop''' (runmd.F90)
<ol>
<li>'''Calculate Forces</li>
<pre>
call force(x,f,ener)  ! AMBER MM potential and forces
! x - in t
! v - in t-dt/2
! f - in t
! ener%pot%tot - in t
! ener%kin%tot - in t-dt
if( ntb .ne. 0 ) then
    call pmf_sander_update_box(a,b,c,alpha,beta,gamma)
end if
call pmf_sander_force(natom,x,v,f,ener%pot%tot,ener%kin%tot,pmfene)
ener%pot%tot = ener%pot%tot + pmfene
</pre>
<li>'''Update Velocities</li>
<pre>
! backup old velocities
vold(t-dt/2) = v(t-dt/2)
! update velocities
v(t+dt/2) = v(t-dt/2) + a(t)*dt
</pre>
<li>'''Update Positions</li>
<pre>
! backup positions
f(t) = x(t)
! update positions
x(t+dt) = x(t) + v(t+dt/2)*dt
</pre>
<li>'''Apply Constraints''' for SHAKE or PMFLib constraints</li>
<pre>
call pmf_sander_constraints(natom,x,con_modified)
call shake
! f - old positions in t
! x - updated positions in t+dt
! re-evaluate velocities
v(t+dt/2) = [x(t+dt) - f(t)]/dt
</pre>
<li>'''Calculate Kinetic Energy'''</li>
<pre>
v(t) = [v(t+dt/2) + vold(t-dt/2)] / 2
ener%kin%tot <- v(t)
</pre>
</ol>
</li>
<li>'''Finalization''' (sander.F90)</li>
<pre>
call pmf_sander_finalize()
</pre>
</ol>