help-minimization and equilibration of Ubiquitin in a water box

From: Zhi Qi (zqi_at_purdue.edu)
Date: Thu Jan 12 2006 - 18:38:00 CST

Dear NAMD,

Thank you for your help last time. And sorry to trouble you again, I
meet another little weird problem during minimization and equilibration
of Ubiquitin in a water box.

Here are the steps I did (by reading the NAMD Tutorial):
----------------------------------------------------------------------------
Step1: Take the crystallographic structure of ubiquitin 1UBQ, and
protonate it, including crystallographic water

Step2: Immersed the structures in a 10.0 A shell of water molecules (TIP3P).

Step3: Add Ions. (for Ubiquitin, total charge is 0, so no this step)

Step4: Minimized the structure with 100 steps of restrained Powell
minimization on all atoms.

Step5: Then, ~Sfreeze~T all heavy atoms in protein, including
crystallographic water --Simulated Annealing:
(5.1) A heating stage from 100K to 500K in steps of 50K with 200 steps
of MD at each temperature.
(5.2) A refinement stage with 200 steps at 500K.
(5.3) A cooling stage from 500K to 25K in steps of 95K with 200 steps of
MD at each temperature.

Step6: Minimized the structure with 200 steps of restrained Powell
minimization on all atoms.
-----------------------------------------------------------------------------

The problem is at Step5, Here is the conf. file.

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

# Minimization and Equilibration of
# proteinL in a Water Box

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

structure ubq_step3_Ions_wb.psf
coordinates Minimization_step5_MD_V4.coor

set temperature 308
# 35 degree + 273

set outputname equil_water_Box_V5_MD_V4

firsttimestep 0

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

# Input
paraTypeCharmm on
parameters par_all22_prot.inp
temperature $temperature

# Force-Field Parameters
exclude scaled1-4
1-4scaling 1.0
cutoff 9.5

#dielectric 2.0
# for water, 1.0; for protein, 2 to 20

switching on
switchdist 8
pairlistdist 11.5
#margin 2.5

# Integrator Parameters
timestep 2.0 ;# 2fs/step

#SHAKE use
rigidBonds all ;# needed for 2fs steps
nonbondedFreq 1
fullElectFrequency 2 ;# no more than 4fs
stepspercycle 10
#rigidTolerance 0.00001

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

# Periodic Boundary Conditions
cellBasisVector1 50. 0. 0.
cellBasisVector2 0. 49. 0.
cellBasisVector3 0. 0 60.
cellOrigin -0.4015398 -18.0065 -42.967

wrapAll on

# PME (for full-system periodic electrostatics)
PME yes
PMEGridSizeX 50
#5^2*2
PMEGridSizeY 50
#5^2*2
PMEGridSizeZ 64
#5*3*2^2

# Constant Pressure Control (variable volume)
  useGroupPressure yes ;# needed for rigidBonds
  useFlexibleCell no
  useConstantArea no

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

# Output
outputName $outputname
binaryoutput off

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

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

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

# Minimization

fixedAtoms on
fixedAtomsForces on
fixedAtomsFile ubq_step53_Ions_heavyfreeze_wb.pdb
fixedAtomsCol B

#minimize 100

#################################
reinitvels $temperature
#################################

#############################################################
## Simulated Annealing ##
#############################################################

#reassignFreq 100
#reassignTemp 1000
#reassignIncr -10
#reassignHold 300
#reassignFreq 100

for { set TEMP 100 } { $TEMP <= 500 } { incr TEMP 50 } {

   run 200
   langevinTemp $TEMP
   langevinPistonTemp $TEMP
   }

run 1000

for { set TEMP 500 } { $TEMP >= 25 } { incr TEMP -95 } {

   langevinTemp $TEMP
   langevinPistonTemp $TEMP
   run 200

}

Questions:

(1) for Simulated Annealing, is the script right?
#############################################################
## Simulated Annealing ##
#############################################################

#reassignFreq 100
#reassignTemp 1000
#reassignIncr -10
#reassignHold 300
#reassignFreq 100

for { set TEMP 100 } { $TEMP <= 500 } { incr TEMP 50 } {

   run 200
   langevinTemp $TEMP
   langevinPistonTemp $TEMP

}

run 1000

for { set TEMP 500 } { $TEMP >= 25 } { incr TEMP -95 } {

   langevinTemp $TEMP
   langevinPistonTemp $TEMP

   run 200

}

(2) when I use this script, if I change temperature from 160K -> 500k ->
160K, there is no problem, but if I change from 100K -> 500K -> 25K,
then it give some mistake information:

ERROR: Exceeded RATTLE iteration limit for atom 1250!
ERROR: Constraint failure; simulation has become unstable.
ERROR: Exceeded RATTLE iteration limit for atom 1391!
ERROR: Constraint failure; simulation has become unstable.
ERROR: Exiting prematurely.

Or sometimes the information:
FATAL ERROR: Periodic cell has become too small for original patch grid!
Possible solutions are to restart from a recent checkpoint,
increase margin, or disable useFlexibleCell for liquid simulation.
------------- Processor 0 Exiting: Called CmiAbort ------------
Reason: FATAL ERROR: Periodic cell has become too small for original
patch grid!Possible solutions are to restart from a recent checkpoint,
increase margin, or disable useFlexibleCell for liquid simulation.

(3) Then, I make some modification:

  (A) margin 2.5
set margin as 2.5 (at the beginning, I use default value)
  (B) timestep 1.0
change timestep from 2.0 to 1.0

Then, I can change temperature from 30K -> 500k -> 30K, but I can not
decrease to 25K

Could you tell me whether I can modify margin and timestep? Is that
reasonable? How to deal with this problem: I mean, in Simulated
Annealing, decrease the temperature to 25K.

BTW, I found that if I do not fix the heavy atoms in protein, there is
no any problem, and do you think the problem came from this?

In my feeling, I do not think temperature will give so great influence,
anyway, hope get your help and suggestion.

Sorry for this so long letter, and have a great day!
Zhi Qi
Physics Department
Purdue University

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:43:14 CST