energy jump running a simulation after minimization

From: Eva Gonzalez Noya (evanoya_at_uchicago.edu)
Date: Tue Feb 15 2011 - 16:37:59 CST

Hi:

I am learning how to use NAMD and found that when running a MD
simulation after a minimization the total energy suffers a jump when
switching from the minimization to the simulation run (the conf file I
am using is copied at the end of this message).

The most surprising thing is that in the log file there are two lines
with the same step number (the last step of minimization and the first
on of the MD run) and the energies are different for those lines:

ENERGY: 10 9830.6905 6728.5509 1233.5493
56.2839 -306280.7690 33101.0201 1709.8377
0.0000 0.0000 -253620.8366 0.0000
-253620.8366 -253620.8366 0.0000 -4869.3857
-4869.1463 639048.6136 -4869.3857 -4869.1463

TCL: Running for 10 steps
ETITLE: TS BOND ANGLE DIHED
IMPRP ELECT VDW BOUNDARY
MISC KINETIC TOTAL TEMP
POTENTIAL TOTAL3 TEMPAVG PRESSURE
GPRESSURE VOLUME PRESSAVG GPRESSAVG

ENERGY: 10 1330.6813 1237.8926 1233.3710
56.2073 -277911.3863 32478.0121 1710.0486
0.0000 0.0000 -239865.1735 0.0000
-239865.1735 -239868.4111 0.0000 -2159.4309
-2095.4452 639048.6136 -2159.4309 -2095.4452

As I minimized previously the structure for 16,000 steps I would expect
it to be close to the minimum so that when I do the MD simulation, atoms
shouldn't move much and the energy should be similar to the one of the
last step of the minimization.

Also I wonder whether the line labeled as 10 in the simulation run
corresponds to the evaluation of the energy before starting the
simulation or it corresponds to the energy after the first Molecular
Dynamics step. If it is the energy before starting the MD simulation I
think that I should exactly match the last step of the minimization,
doesn't it?

I copied below the conf file I am using (I am restarting it from a
previous minimization run of 16000 steps, so that I expect to be close
to a local minimum).

I tried already using the velocity restart files instead of assigning
the temperature a very small value and does not fix the problem.

Anyone found a similar problem before? Why is this jump in energy
happening?

Than you very much in advance for you help

Eva

#############################################################
## JOB DESCRIPTION ##
#############################################################

# Minimization and Equilibration of
# Ubiquitin in a Water Box

#############################################################
## ADJUSTABLE PARAMETERS ##
#############################################################

structure ./apo-tubz_wb_ionized.psf
coordinates ./apo-tubz_wb_ionized.pdb

set temperature 0.000000001
set outputname apo-tubz_wb_ionized_eq

if {1} {
set inputname apo-tubz_wb_ionized_eq
binCoordinates $inputname.restart.coor
#binVelocities $inputname.restart.vel ;# remove the "temperature"
entry if you use this!
extendedSystem $inputname.xsc
}

firsttimestep 0

#############################################################
## SIMULATION PARAMETERS ##
#############################################################

# Input
paraTypeCharmm on
parameters
/scratch/noya/apo-monomer/common/par_all27_prot_lipid.inp
temperature $temperature

# Force-Field Parameters
exclude scaled1-4
1-4scaling 1.0
cutoff 12.0
switching on
switchdist 10.0
pairlistdist 13.5

# Integrator Parameters
timestep 2.0 ;# 2fs/step
rigidBonds all ;# needed for 2fs steps
nonbondedFreq 1
fullElectFrequency 2
stepspercycle 10

# Constant Temperature Control
#if {0} {
#langevin on ;# do langevin dynamics
#langevinDamping 0.5 ;# damping coefficient (gamma) of 5/ps
#langevinTemp $temperature
#langevinHydrogen off ;# don't couple langevin bath to hydrogens
#}

# Periodic Boundary Conditions
#cellBasisVector1 81.7999997138977 0. 0.
#cellBasisVector2 0. 83.28000068664551 0.
#cellBasisVector3 0. 0 93.80800342559814
#cellOrigin 47.18317794799805 59.78273010253906 48.607540130615234

wrapAll on

# PME (for full-system periodic electrostatics)
PME yes
PMEGridSpacing 1.0

#manual grid definition
#PMEGridSizeX 45
#PMEGridSizeY 45
#PMEGridSizeZ 48

# Constant Pressure Control (variable volume)
#useGroupPressure yes ;# needed for rigidBonds
#useFlexibleCell no
#useConstantArea no
#
#langevinPiston on
#langevinPistonTarget 1.01325 ;# in bar -> 1 atm
#langevinPistonPeriod 100.0
#langevinPistonDecay 50.0
#langevinPistonTemp $temperature

# Output
outputName $outputname

restartfreq 500 ;# 500steps = every 1ps
dcdfreq 250
xstFreq 250
outputEnergies 1
outputPressure 100

#############################################################
## EXTRA PARAMETERS ##
#############################################################

# bound the Calpha of TubZ to harmonic springs with strength 100kcal/mol/A

constraints on
consexp 2
consref apo_tubz_fixed.pdb
conskfile apo_tubz_fixed.pdb
conskcol O
constraintScaling 100.0

#############################################################
## EXECUTION SCRIPT ##
#############################################################

# Minimization
#reassignFreq 100 # 0.2ps
#reassignTemp 1
#reassignIncr 6.2 # this means heating at a reta 31K/ps
#reassignHold 310
#
minimize 10
run 10

This archive was generated by hypermail 2.1.6 : Mon Dec 31 2012 - 23:19:49 CST