Molecular Dynamics Algorithms

Navigation: Documentation / Technical Details / Molecular Dynamics Algorithms

---

AMBER

Both sander and pmemd from the AMBER package use the leap-frog integration algorithm. The algorithm is demonstrated for the SANDER driver:

1. Initialization (sander.F90)

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))

2. MD Loop (runmd.F90)
   1. Calculate Forces

   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
   

   2. Update Velocities

   ! backup old velocities
   vold(t-dt/2) = v(t-dt/2)
   ! update velocities
   v(t+dt/2) = v(t-dt/2) + a(t)*dt
   

   3. Update Positions

   ! backup positions
   f(t) = x(t)
   ! update positions
   x(t+dt) = x(t) + v(t+dt/2)*dt
   

   4. Apply Constraints for SHAKE or PMFLib constraints

   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
   

   5. Calculate Kinetic Energy

   v(t) = [v(t+dt/2) + vold(t-dt/2)] / 2
   ener%kin%tot <- v(t)
   

3. Finalization (sander.F90)

call pmf_sander_finalize()