Restart, discontinuity of energy

From: juan du (dujuan709_at_gmail.com)
Date: Wed Dec 01 2010 - 09:56:51 CST

Hi, everyone,
I am continuing an NAMD equilibration runs using restart files.
I obtained an odd results. I give the last step of previous run as the
firsttimestep for the new run. I noticed the discontinuity of energy
in restart calculations. What's more, I calculated
the RMSD of the backbone atoms using the initial structure as
reference. The RMSD values increase significantly after each restart
and then reach stability. I wonder what's the problem with the
restart. Your kind help will be greatly appreciated.

Juan Du
Last step
ENERGY: 250000 43292.9642 29263.0362 3945.3575
441.5854 -506610.3555 47059.5022 1164.0174
0.0000 117463.1912 -263980.7016 311.2331 -381443.8928
 -262914.2342 309.6929 -239.3309 -267.4445
1220176.8696 17.9805 19.0874

new run first step
ENERGY: 0 43292.9642 29263.0362 3945.3575
441.5854 -506613.0293 47061.8423 0.0000
0.0000 117462.0335 -265146.2102 311.2300 -382608.2437
 -264084.8769 311.2300 186.2934 158.2625
1220172.1958 186.2934 158.2625

Laststep
# NAMD extended system configuration output file
#$LABELS step a_x a_y a_z b_x b_y b_z c_x c_y c_z o_x o_y o_z s_x s_y
s_z s_u s_v s_w
250000 106.858 0 0 0 106.858 0 0 0 106.858 0.037 0.054 -0.008
2.27649e-05 2.27649e-05 2.27649e-05 0 0 0

new run first
# NAMD extended system trajectory file
#$LABELS step a_x a_y a_z b_x b_y b_z c_x c_y c_z o_x o_y o_z s_x s_y
s_z s_u s_v s_w
0 106.858 0 0 0 106.858 0 0 0 106.858 0.037 0.054 -0.008 2.27649e-05
2.27649e-05 2.27649e-05 0 0 0
1000 106.885 0 0 0 106.885 0 0 0 106.885 0.037 0.054 -0.008
-2.07681e-05 -2.07681e-05 -2.07681e-05 0 0 0

.conf file

#set temperature 310.0
set targetemp 310.0
set coorfile ../hv57-334.pdb
set laststep myeqrst1

#--Input File
coordinates $coorfile
structure ../hv57-334.psf
bincoordinates ${laststep}.coor
binvelocities ${laststep}.vel
extendedSystem ${laststep}.xsc

#--Output File
outputname myeqout2
restartname myeqrst2
outputEnergies 100
binaryoutput yes
xstFreq 1000
DCDFreq 1000
restartFreq 1000

#--Initialization
COMmotion no
timestep 1 ;# 2fs/step
#rigidBonds all ;# needed for 2fs/step
nonbondedFreq 1
fullElectFrequency 2

#--Temperature ReScaling
rescaleFreq 100
#rescaleTemp $temperature
rescaleTemp $targetemp

#--Constant T
langevin on
langevinDamping 1.
langevinTemp $targetemp
langevinHydrogen off ;# do not couple bath to hydrogens

#--Constant P
useGroupPressure yes
useFlexibleCell no
useConstantArea no
useConstantRatio no

langevinPiston on
langevinPistonTarget 1.01325 ;# in bar --> 1 atm
langevinPistonPeriod 100.
langevinPistonDecay 50.
langevinPistonTemp $targetemp

#--PME
PME on
PMETolerance 0.000001
PMEGridSpacing 1.0

#--Constraint
constraints on
consexp 2
consref $coorfile
conskfile ../prot-fix2.cnst ;# force constant value
conskcol B
constraintScaling 20.0 ;# the FC is 20 kcal/mol Angstrom
margin 10

#--parameter file
paraTypeCharmm on
parameters ../par_all27_prot_lipid.prm

exclude scaled1-4
1-4scaling 1.0
dielectric 1.0
switching on
switchdist 10.0
cutoff 12.0
pairlistdist 14.0 ;# cutoff + 2
stepspercycle 10 ;# follow NAMD configuration file

wrapWater on
wrapAll on

#--simulation time
run 250000

RMSD
mol new hv57-334.pdb
mol new hv57-334.psf type psf first 0 last -1 step 1 filebonds 1
autobonds 1 waitfor all
mol addfile myeqout1.dcd type dcd first 0 last -1 step 1 filebonds 1
autobonds 1 waitfor all
mol addfile myeqout2.dcd type dcd first 0 last -1 step 1 filebonds 1
autobonds 1 waitfor all

proc print_rmsd_through_time {} {
     set outfile [open rmsd.dat w]
     set reference [atomselect 0 "protein and backbone"]
     set compare [atomselect 1 "protein and backbone"]
     set num_steps [molinfo 1 get numframes]

     for {set frame 0} {$frame < $num_steps} {incr frame} {
             $compare frame $frame
             set trans_mat [measure fit $compare $reference]
             $compare move $trans_mat
             set rmsd [measure rmsd $compare $reference]
             puts $outfile "RMSD of $frame is $rmsd"
     }
}

print_rmsd_through_time

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:54:49 CST