From: Arneh Babakhani (ababakha_at_mccammon.ucsd.edu)
Date: Wed Feb 22 2006 - 14:09:23 CST
Hi All,
Just a quick follow-up to this message that I posted last week, in case
any future membrane simulators run into the same problem!
First, thank you all for your suggestions. I learned a lot about MD and
NAMD while trouble-shooting this. It turns out that the problem was in
my original membrane construction; I used a single (poorly structured)
DMPC residue, to make my membrane. The result was a system with a lot
of steric clashes and bad contacts (I did not bother to quantitate or
specifically identify which interactions were causing the problem. You
could tell it was a bad structure just by looking at it).
These bad interactions could not be resolved, no matter how long I ran
the minimization. So eventually, when I tried to heat up the system and
apply constant pressure, the simulation would become unstable and abort.
Upon using a better template DMPC (generously donated by my labmate,
Justin Gullingsrud) , in which there we no such clashes, I built a more
stable membrane. The minimization proceeded more smoothly, I was able
to heat up, apply constant pressure, equilibrate, with no "RATTLE"
error or aborting.
So the morale of the story is a classic one in computational science:
Garbage in equals garbage out! Be mindful of what you use as your
building block for your membrane system.
Thanks!
Arneh
Arneh Babakhani wrote:
> Dear NAMD community,
>
> I constructed and solvated a lipid bilayer membrane. This consists of
> 200 DMPC lipids (100 per monolayer), about 5100 water molecules,
> for a total of ~40,000 atoms.
>
> I conducted 10,000 steps of minimization, ultimately achieving a
> gradient tolerance of ~2.91. At the end of this minimization, I
> zeroed the velocities of all atoms. No problem here . . .
>
> I then gradually heated up my system, over 2500 steps (2fs timestep),
> to 330K. Looking at the plot of the temperature, the system does heat
> up gradually, eventually plateauing out at 330K, towards the end of
> the run. No problem here . . .
>
> I then turn on constant pressure and attempt to run for another 2500
> steps. The temperature spikes from 330K to ~375K, by the second step!
> (I am implementing langevin temp control). The temp. gradually
> decreases through the course of the simulation, but it never quite
> gets back down to 330K. NAMD prematurely aborts halfway through the
> run, at about 350K, giving the following error:
>
> ERROR: Constraint failure in RATTLE algorithm for atom 13561!
>
> ERROR: Constraint failure; simulation has become unstable.
>
> ERROR: Exiting prematurely.
>
>
> I found through some trial & error that by increasing the langevin
> Damping parameter, I can decrease the magnitude of this initial
> temperature spike (although I can't get rid of it completely, very
> strange). By decreasing this spike, the system comes back down to
> 330K much faster, and the simulation proceeds to completion (with
> langevin Damping = 7).
>
> I then try to run a 1 ns equilibration of my membrane, again using
> lang. Damping = 7 (since this seemed to work for the constant pressure
> run). This runs for about 850 ps and then aborts, with the same
> RATTLE algorithm error. However, during this equilibration, the temp.
> stays constant, as does the pressure (no spikes), and nothing else is
> observably wrong. I repeated the attempt for lang. Damping =5, and
> 9, no joy, still aborts.
>
> Would greatly appreciate any ideas? (I've pasted my conf file below
> for reference)
>
> Regards,
>
> Arneh
>
> #############################################################
> # 1NS equilibration dmpc lipid bilayer in a Water Box
> #############################################################
> structure ../common/solvated.psf
> coordinates ../common/solvated.pdb
> bincoordinates lipwaterconstpress2.coor
> binvelocities lipwaterconstpress2.vel
> extendedsystem lipwaterconstpress2.xsc
>
> set temperature 330
> set outputname lipwaterequil1
>
> firsttimestep 0
>
> ## SIMULATION PARAMETERS ##
>
> paraTypeCharmm on
> parameters ../common/par_all27_prot_lipid.inp
>
> # Force-Field Parameters
> exclude scaled1-4
> 1-4scaling 1.0
> cutoff 10.
> switching on
> switchdist 9.
> pairlistdist 11.5
>
>
> # Integrator Parameters
> timestep 2.0 ;# 2fs/step
> rigidBonds all ;# needed for 2fs steps
> nonbondedFreq 1
> fullElectFrequency 2
> stepspercycle 10
>
>
> # Constant Temperature Control
> langevin on ;# do langevin dynamics
> langevinDamping 7 ;# damping coefficient (gamma) of 9/ps
> langevinTemp $temperature
> langevinHydrogen off ;# don't couple langevin bath to hydrogens
>
> wrapAll on
>
> # PME (for full-system periodic electrostatics)
> PME yes
> PMEGridSizeX 64
> PMEGridSizeY 64
> PMEGridSizeZ 64
>
> # 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.
> langevinPistonDecay 200.
> langevinPistonTemp $temperature
>
> # Output
> outputName $outputname
> restartfreq 500 ;# 500steps = every 1ps
> dcdfreq 500
> outputEnergies 100
> outputPressure 100
>
> run 500000 ; # 1ns
>
>
>
>
>
-- ---------------------------- Arneh Babakhani University of California at San Diego Physical Chemistry / Department of Chemistry & Biochemistry Laboratory of Prof. J. A. McCammon 9500 Gilman Drive MC 0365 La Jolla, CA 92093-0365 (619)895-6540 (858)534-4974 (FAX) http://mccammon.ucsd.edu/~ababakha/ ababakha_at_mccammon.ucsd.edu
This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:41:39 CST