Pressure in equilibration stage of MD

From: Mariana Graterol (marianagraterol_at_gmail.com)
Date: Thu Oct 18 2012 - 16:52:05 CDT

Dear NAMD user's:
I have a problem to balance the pressure in my system to the
conditions of 300 K and 1 atm.
First, I am modeling a protein using periodic boundary conditions. I
built the box and neutralized (considering an optimal size between
images). I have done the equilibration in stages:
The first stage corresponds to the relaxation of the water added with
minimization steps: 1) fix the protein (relaxing the water), 2)
restricting the backbone of it. approximately during 30 ps.
Second stage: It consisted in gradual heating to 100 K (in total 200
ps, 0,0005 K heating step each step of 1 fs).
Third stage: NPT run to adjust volume;at this point I obtained my box
at 100 K with a suitable volume to obtain the aqueous solution density
1.024 g/ml.

The issue arises when I want to reach 300 K, when I do it slow enough
(I think) in a NVT run for 800 ps in total, and 0,00025 K each step,
the system reaches the desired temperature, but rather the pressure
increases (to approx. PRESSAVG 5000 bar).
Even after warming up, I simulated 300 ps additional unheated, and the
pressure does not relax in the NVT system using Langevin Control ON,
yielding values of PRESSAVG ranging between 3000 and 6000 Bar.

In trying to balance the pressure again, I made anothe NPT run, this
time at 300 K, the pressure decreases at the expense of an increase of
10 times the initial volume. Which obviously does not correspond to my
equilibrate system.

I hesitated about pressure values in NAMD, are they nominal? after
completing the protocol: Minimize-> Heat 100 K-> balance volume (NPT)
-> Heat 300 K, it is necessary to lower these values of pressure,
which may affect my later production run?

I have read about calculating NAMD pressure, which depends on the
virial, and therefore its fluctuations are about 1000 Bar, but I do
not know if this increase is to be attributed to my system failure or
disruption by increasing temperature.

Thanking you for any comments, is the first time I dynamics.
Below I include the keywords in the protocols used:

Second stage (heat 100 K)

# Input files
coordinates C:\\Users\\serranom\\Documents\\Mariana\\eq100K\\br100p.pdb
structure C:\\Users\\serranom\\Documents\\Mariana\\eq100K\\br100p.psf
binvelocities C:\\Users\\serranom\\Documents\\Mariana\\eq100K\\br100p2.vel
bincoordinates
C:\\Users\\serranom\\Documents\\Mariana\\eq100K\\br100p2.coor
parameters C:\\Users\\serranom\\Desktop\\parm.prm
paraTypeCharmm on

# Output files
binaryoutput yes
outputname C:\\Users\\serranom\\Documents\\Mariana\\eq100K\\br100
restartfreq 1000
restartsave yes
binaryrestart yes
DCDfile C:\\Users\\serranom\\Documents\\Mariana\\eq100K\\br100.dcd
DCDfreq 100
DCDUnitCell no
outputEnergies 100
mergeCrossterms yes
outputMomenta 0
outputPressure 100
outputTiming 10000

# Basic dynamics
exclude scaled1-4
COMmotion no
zeroMomentum no
dielectric 1.000000
nonbondedScaling 1.000000
1-4scaling 1.000000
vdwGeometricSigma no
seed 12345
rigidBonds none

# PME parameters
PME on
PMETolerance 1.000000e-06
PMEInterpOrder 4
PMEGridSpacing 1.000000
PMEGridSizeX 120
PMEGridSizeY 120
PMEGridSizeZ 90
FFTWEstimate no
FFTWUseWisdom yes

# Full direct parameters
FullDirect no

# Multiple timestep parameters
fullElectFrequency 4
MTSAlgorithm impulse
longSplitting c1

# Harmonic constraints
constraints on
consexp 2
conskcol B
constraintScaling 1.000000
conskfile C:\\Users\\MS\\Documents\\Mariana\\cal300\\br100p.pdb
consref C:\\Users\\MS\\Documents\\Mariana\\cal300\\br100p.pdb

# Periodic boundary conditions
cellBasisVector1 108.447000 0.000000 0.000000
cellBasisVector2 0.000000 111.005000 0.000000
cellBasisVector3 0.000000 0.000000 85.008000
cellOrigin 0.000000 0.000000 0.000000
XSTfreq 1000
wrapWater on
wrapAll on
wrapNearest off

# Temperature reassignment
reassignFreq 30
reassignTemp 0.000
reassignIncr 0.0015
reassignHold 100.00

Third stage (NPT run)

..
# Timestep parameters
numsteps 700000
timestep 1.000000
firsttimestep 200000
stepspercycle 20

# Basic dynamics
exclude scaled1-4
COMmotion no
zeroMomentum no
dielectric 1.000000
nonbondedScaling 1.000000
1-4scaling 1.000000
vdwGeometricSigma no
seed 12345
rigidBonds none
..

# Langevin dynamics
langevin on
langevinTemp 100.000000
langevinHydrogen off
langevinDamping 5.000000

# Constant Pressure Control

useGroupPressure no

useFlexibleCell no

useConstantArea no

langevinPiston on

langevinPistonTarget 1.01325

langevinPistonPeriod 1000

langevinPistonDecay 500

langevinPistonTemp 100

Heating 300 K

# Temperature reassignment
reassignFreq 20
reassignTemp 100.000000
reassignIncr 0.005000
reassignHold 300.000000

NVT Run at 300 K

# Langevin dynamics
langevin on
langevinTemp 300.000000
langevinHydrogen off
langevinDamping 5.000000

--
* mari *

This archive was generated by hypermail 2.1.6 : Tue Dec 31 2013 - 23:22:40 CST