# heat up MD system through Berendsen thermostat

#

msf="step6_heatup.msf"

mds="step9_equil_npt_more.mds"

rst="step9_equil_npt_more.rst"

nstep=2500  #number steps for tempering

ttar=500    #target temperature

#load md system with minimized state

loadmsf file="$msf"

#load mds and state

loadmds file="$mds"

loadrst file="$rst"

#change to Berendsen thermostat for tempering

ttcoup=0.1*($nstep/1000)

ensemble nvt tcoup=Berendsen ttcoup=$ttcoup tref=$ttar

mdcp trjfile="step10_heat_to400k.trj" frqtrj=100 trjfmt=modar \

     rstfile="step10_heat_to400k.rst" frqrst=1000 \

     enefile="step10_heat_to400k.ene" frqene=10 \

     resetstartstep=true

runmd nstep=$nstep

#change to Nose-Hoover thermostat

ensemble nvt tcoup=Nose tmass=1000

mdcp trjfile="step10_equil_nvt.trj" frqtrj=100 trjfmt=modar \

     rstfile="step10_equil_nvt.rst" frqrst=10000 \

     enefile="step10_equil_nvt.ene" frqene=10

runmd nstepmore=2500

#change to NPT

tell tmass #print total mass, and put total mass number to $totmass

set pmass = int($totmass/50.0)

ensemble npt pmass=$pmass pgamma=20.0 \

             pxx=1.0 pyy=1.0 pzz=1.0 \

             tmass=2000

mdcp trjfile="step10_equil_npt.trj" frqtrj=100 trjfmt=modar \

     rstfile="step10_equil_npt.rst" frqrst=10000 \

     enefile="step10_equil_npt.ene" frqene=10

runmd nstepmore=2500

#save mds and pdb

savemds file="step10_equil_npt.mds"

savepdb file="step10_equil_npt.pdb"