Re: FATAL ERROR: Periodic cell has become too small for original patch grid!

From: Joshua Adelman (jadelman_at_berkeley.edu)
Date: Wed Jan 30 2008 - 11:47:16 CST

Hi Peter,

Can you give me an example of how to implement the heating protocol
using langevinTemp? Your thermostat would have something like:
# Constant Temperature Control
langevin on ;# do langevin dynamics
langevinDamping 5 ;# damping coefficient (gamma) of 5/ps
langevinTemp $temperature
langevinHydrogen no ;# don't couple langevin bath to hydrogens

Would you just tack the following on the end of your input script?
set temperature 15
langevinTemp $temperature
run 1000

set temperature 30
langevinTemp $temperature
run 1000

set temperature 45
langevinTemp $temperature
run 1000

etc.?

Josh

On Jan 30, 2008, at 9:34 AM, Peter Freddolino wrote:

> Hi Josh,
> actually, what I'd recommend instead is doing your initial heating
> in NPT; *more* NV(E) equilibration is unlikely to help because it is
> actually exacerbating the problem that because solvation protocols
> are imperfect, we usually end up with a box that needs to shrink
> somewhat; if you then feed the barostat a box that has shrunk,
> there's a lot of vacuum and bad things start to happen.
>
> If you have constraints on your solute doing the heating in NPT is
> unlikely to hurt the system; you can change the langevin dynamics
> temperature with langevinTemp $NEWTEMP and follow the same heating
> protocol as before. You may also want to increase your piston period
> and piston decay by a factor of 2; 100/50 and 200/100 are both
> reasonable values, but the latter averages over a longer amount of
> time and thus is more likely to tolerate rapid initial fluctuations.
>
> BTW, there's no need for it now, but it's generally helpful to
> continue these discussions on-list so that they're searchable in the
> future.
> Best,
> Peter
>
> Joshua Adelman wrote:
>> Hi Peter,
>>
>> Good call. I printed out the DCD and XST files for every step and
>> the system is definitely blowing up:
>> # 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
>> 35000 119.6 0 0 0 87.2 0 0 0 112 0 0 0 0 0 0 0 0 0
>> 35001 119.6 0 0 0 87.2 0 0 0 112 0 0 0 -0.000284665 -0.000284665
>> -0.000284665 0 0 0
>> 35002 119.397 0 0 0 87.0524 0 0 0 111.81 0 0 0 0.000385176
>> 0.000385176 0.000385176 0 0 0
>> 35003 119.397 0 0 0 87.0524 0 0 0 111.81 0 0 0 0.00202167
>> 0.00202167 0.00202167 0 0 0
>> 35004 120.763 0 0 0 88.0481 0 0 0 113.089 0 0 0 -0.00277211
>> -0.00277211 -0.00277211 0 0 0
>> 35005 120.763 0 0 0 88.0481 0 0 0 113.089 0 0 0 -0.014042 -0.014042
>> -0.014042 0 0 0
>> 35006 111.611 0 0 0 81.3755 0 0 0 104.519 0 0 0 0.0167054 0.0167054
>> 0.0167054 0 0 0
>> 35007 111.611 0 0 0 81.3755 0 0 0 104.519 0 0 0 0.0898706 0.0898706
>> 0.0898706 0 0 0
>> 35008 185.206 0 0 0 135.033 0 0 0 173.437 0 0 0 -0.324938 -0.324938
>> -0.324938 0 0 0
>> 35009 185.206 0 0 0 135.033 0 0 0 173.437 0 0 0 -1.23019 -1.23019
>> -1.23019 0 0 0
>>
>> Although it starts to expand a little then shrinks a lot and then
>> explodes, I'm guessing because it got compressed and is reacting to
>> that. I never get any errors about atoms moving too fast though.
>>
>> Would the solution be to equilibrate longer in NVT using the
>> restraints before turning on the barostat? Additionally do you
>> think it would help to reduce the restraints on non-backbone atoms
>> while still in NVT, before switching to NPT (and then I would relax
>> the backbone restraints in the NPT ensemble)?
>>
>> Thanks again for your help.
>>
>> Josh
>>
>> On Jan 30, 2008, at 8:52 AM, Peter Freddolino wrote:
>>
>>> Hi Josh,
>>> you may want to check how immediate 'immediately' really is. To do
>>> this, try the same NPT run as before, but set your restart,
>>> output, and dcd frequencies. See if it's really crashing on step
>>> 0, or if the periodic cell rapidly shrinks before the crash (have
>>> a look at your xst file). I've seen similar behavior in the past
>>> when the water box contracts too much, too quickly; you may also
>>> want to look at your velocity rescaling trajectories and see
>>> whether there are any voids opening up during that portion of the
>>> simulation, which would leave the barostat with a lot to fix very
>>> quickly once you turn on pressure control.
>>> Best,
>>> Peter
>>>
>>> Joshua Adelman wrote:
>>>> I am try to equilibrate a system, by first using harmonic
>>>> constraints (restraints) on the protein atoms. I run a short
>>>> simulation in which I hold the protein atoms and let the water
>>>> and ions go as I slowly heat the system using progressive
>>>> temperature rescaling without any thermo- or barostats. Once I
>>>> reach 300K, I try to run another simulation where I relax the
>>>> restraints on the non-backbone atoms of the protein in the NPT
>>>> ensemble, but the simulation immediately crashes out giving the
>>>> error:
>>>>
>>>> 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.
>>>>
>>>> I tried increasing the margin to 2.0, but it does the same thing.
>>>> useFlexibleCell is off already. If I turn off the barostat (NVT),
>>>> the simulation runs fine and if I remove the constraints
>>>> completely, NPT runs fine. I'm not sure how the barostat and the
>>>> constraints are interacting to cause this problem. I'm including
>>>> the input script below.
>>>>
>>>> Any suggestions would be appreciated.
>>>>
>>>> Josh
>>>>
>>>> #############################################################
>>>> ## ADJUSTABLE PARAMETERS ##
>>>> #############################################################
>>>>
>>>>
>>>> structure ../2sub2.psf
>>>> coordinates ../2sub2.pdb
>>>> outputName 2s2-equil4
>>>>
>>>> set temperature 300.0
>>>>
>>>> # Continuing a job from the restart files
>>>> set inputname 2s2-equil3
>>>> binCoordinates $inputname.restart.coor
>>>> binVelocities $inputname.restart.vel ;# remove the
>>>> "temperature" entry if you use this!
>>>> extendedSystem $inputname.restart.xsc
>>>>
>>>> set xscfile $inputname.restart.xsc
>>>>
>>>> #########################################################
>>>> proc get_first_ts { xscfile } {
>>>> set fd [open $xscfile r]
>>>> gets $fd
>>>> gets $fd
>>>> gets $fd line
>>>> set ts [lindex $line 0]
>>>> close $fd
>>>> return $ts
>>>> }
>>>>
>>>> set firstts [get_first_ts $xscfile]
>>>> ########################################################
>>>>
>>>> firsttimestep $firstts
>>>>
>>>>
>>>> #############################################################
>>>> ## SIMULATION PARAMETERS ##
>>>> #############################################################
>>>>
>>>> # Input
>>>> paraTypeCharmm on
>>>> parameters ../par_all27_prot_na.prm
>>>>
>>>>
>>>> # NOTE: Do not set the initial velocity temperature if you
>>>> # have also specified a .vel restart file!
>>>>
>>>>
>>>> # Periodic Boundary conditions
>>>> # NOTE: Do not set the periodic cell basis if you have also
>>>> # specified an .xsc restart file!
>>>> #cellBasisVector1 119.5 0. 0.
>>>> #cellBasisVector2 0. 87. 0.
>>>> #cellBasisVector3 0. 0 111.9
>>>> cellOrigin 0. 0. 0.
>>>>
>>>> wrapWater on
>>>> wrapAll off
>>>>
>>>>
>>>> # Force-Field Parameters
>>>> exclude scaled1-4
>>>> 1-4scaling 1.0
>>>> cutoff 12.
>>>> switching on
>>>> switchdist 10.
>>>> pairlistdist 13.5
>>>> COMmotion yes
>>>> margin 2.0
>>>>
>>>> # Integrator Parameters
>>>> timestep 2.0 ;# 2fs/step
>>>> rigidBonds all
>>>> nonbondedFreq 1
>>>> fullElectFrequency 2
>>>> stepspercycle 10
>>>>
>>>> #PME (for full-system periodic electrostatics)
>>>> PME yes
>>>> PMEGridSizeX 125
>>>> PMEGridSizeY 125
>>>> PMEGridSizeZ 125
>>>>
>>>> # Restraints
>>>> constraints on
>>>> consref 1VYM_2sub2_restrained.pdb
>>>> conskfile 1VYM_2sub2_restrained.pdb
>>>> conskcol B
>>>>
>>>> # Constant Temperature Control
>>>> langevin on ;# do langevin dynamics
>>>> langevinDamping 5 ;# damping coefficient (gamma) of 5/ps
>>>> langevinTemp $temperature
>>>> langevinHydrogen no ;# don't couple langevin bath to
>>>> hydrogens
>>>>
>>>>
>>>> # Constant Pressure Control (variable volume)
>>>> useGroupPressure yes ;# needed for 2fs steps
>>>> useFlexibleCell no ;# no for water box, yes for membrane
>>>> useConstantArea no ;# no for water box, yes for membrane
>>>>
>>>> langevinPiston on
>>>> langevinPistonTarget 1.01325 ;# in bar -> 1 atm
>>>> langevinPistonPeriod 100.
>>>> langevinPistonDecay 50.
>>>> langevinPistonTemp $temperature
>>>>
>>>> restartfreq 500 ;# 1000steps = every 1ps
>>>> dcdfreq 500
>>>> xstFreq 500
>>>> outputEnergies 100
>>>> outputPressure 100
>>>>
>>>> #protocol
>>>> run 5000
>>>>
>>>>
>>>>
>>>>
>>>> ------------------------------------------------------------------------------------------------------
>>>> Joshua L. Adelman
>>>> Biophysics Graduate Group Lab: 510.643.2159
>>>> 218 Wellman Hall Fax: 510.642.7428
>>>> University of California, Berkeley http://nature.berkeley.edu/~jadelman
>>>> Berkeley, CA 94720 USA jadelman_at_berkeley.edu
>>>> ------------------------------------------------------------------------------------------------------
>>>>
>>>>
>>
>> ------------------------------------------------------------------------------------------------------
>> Joshua L. Adelman
>> Biophysics Graduate Group Lab: 510.643.2159
>> 218 Wellman Hall Fax: 510.642.7428
>> University of California, Berkeley http://nature.berkeley.edu/~jadelman
>> Berkeley, CA 94720 USA jadelman_at_berkeley.edu
>> ------------------------------------------------------------------------------------------------------
>>
>>

------------------------------------------------------------------------------------------------------
Joshua L. Adelman
Biophysics Graduate Group Lab: 510.643.2159
218 Wellman Hall Fax: 510.642.7428
University of California, Berkeley http://nature.berkeley.edu/~jadelman
Berkeley, CA 94720 USA jadelman_at_berkeley.edu
------------------------------------------------------------------------------------------------------

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