cleanup of temperature and pressure control for MD integrators

Author: jayponder

source/baoab.f

@@ -156,10 +156,9 @@ subroutine baoab (istep,dt)
 c 
       call gradient (epot,derivs)
 c
-c     make half-step correction to control the pressure
+c     compute the kinetic energy from half-step velocities
 c
       call kinetic (eksum,ekin,temp)
-      call pressure2 (epot,temp)
 c
 c     second B step for accelerations and full-step velocities
 c
@@ -188,11 +187,12 @@ subroutine baoab (istep,dt)
          end do
       end if
 c
-c     compute the kinetic energy and control the pressure;
-c     half-step kinetic energy gives better temperature control
+c     compute full-step kinetic energy and pressure correction;
+c     half-step kinetic energy gives better pressure control
 c
 c     call kinetic (eksum,ekin,temp)
       call pressure (dt,epot,ekin,temp,pres,stress)
+      call pressure2 (epot,temp)
 c
 c     final constraint step to enforce position convergence
 c

source/ghmcstep.f

@@ -146,7 +146,7 @@ subroutine ghmcstep (istep,dt)
 c
       if (use_rattle)  call rattle2 (dt)
 c
-c     determine the kinetic energy, temperature and total energy
+c     compute the kinetic energy, temperature and total energy
 c
       call kinetic (eksum,ekin,temp)
       etot = eksum + epot
@@ -175,7 +175,7 @@ subroutine ghmcstep (istep,dt)
          end do
       end if
 c
-c     evolve velocities according to midpoint Euler for half-step
+c     update velocities using midpoint Euler for half-step
 c
       call ghmcterm (istep,dt,alpha,beta)
       do i = 1, nuse
@@ -225,6 +225,7 @@ subroutine ghmcstep (istep,dt)
 c
       call mdstat (istep,dt,etot,epot,eksum,temp,pres)
       call mdsave (istep,dt,epot,eksum)
+      call mdrest (istep)
       return
       end
 c

source/mdinit.f

@@ -70,14 +70,14 @@ subroutine mdinit (dt)
       nfree = 0
       ndummy = 0
       irest = 100
+      iprint = 100
       use_wrap = .true.
       velsave = .false.
       frcsave = .false.
       uindsave = .false.
       friction = 1.0d0
       if (use_solv)  friction = 91.0d0
       use_sdarea = .false.
-      iprint = 100
 c
 c     set default values for temperature and pressure control
 c