[gmx-users] Not able to continue with Equilibration
Justin A. Lemkul
jalemkul at vt.edu
Thu Mar 29 17:37:51 CEST 2012
francesca vitalini wrote:
> My minimization was done without any constraints. I'm starting using
> constraints only in the nvt equilibration as I'm trying to equilibrate
> my solvent, so, maybe here I'm being extremely naïve, but, shouldn't I
> avoid fixing the water in the nvt, as equilibrating it is exactly the
> reason I'm running the nvt equilibration? Or am I misunderstanding
> something?
You should likely be restraining your protein, not the water. Bad initial
geometry (sometimes stuff survives EM) can give your structure a nasty kick that
causes things to spin out of control. Restraints on the protein mitigate these
effects while the water relaxes.
-Justin
> Thanks
> Francesca
>
>
>
> 2012/3/29 <chris.neale at utoronto.ca>:
>> I didn't follow this whole thread, but I sometimes need to turn off all
>> constraints when doing minimization. In fact, for that reason I entirely
>> stopped ever using restraints during energy minimization. In extreem cases,
>> I have had success also by forcing the water to be flexible with a
>> -DFLEXIBLE (or whatever is appropriate for your water model; check the .itp)
>> statement in the .mdp file.
>>
>> Once EM is done, use rigid water and restraints and everything works out ok.
>>
>> Chris.
>>
>> -- original message --
>>
>> Hallo Felix,
>> thank you for your answer. I tried the constraints = h-bonds but no
>> change in the output. If I look at the step.pdb file that is produced
>> after the running I have some strange outcome. For example some of my
>> atoms are not recognized as part of my protein any more and my
>> structure is destroied. Am I using some strange parameters for nvt
>> without realizing it? I've started from the mdp file provided by the
>> lysozyme tutorial non the Gromacs webpage.
>> if anyone has any ideas it is more than welcome.
>> Francesca
>>
>> integrator = md ; leap-frog integrator
>> nsteps = 10000 ; 2 * 10000 = 20 ps
>> dt = 0.002 ; 2 fs
>> ; Output control
>> nstxout = 1 ; save coordinates every 0.2 ps
>> nstvout = 1 ; save velocities every 0.2 ps
>> nstenergy = 50 ; save energies every 0.2 ps
>> nstlog = 50 ; update log file every 0.2 ps
>> ; Bond parameters
>> unconstrained_start = no
>> ;continuation = no ; first dynamics run
>> constraint_algorithm = lincs ; holonomic constraints
>> constraints = h-bonds ; all bonds (even heavy atom-H bonds)
>> constrained
>> lincs_iter = 1 ; accuracy of LINCS
>> lincs_order = 4 ; also related to accuracy
>> ; Neighborsearching
>> ns_type = grid ; search neighboring grid cells
>> nstlist = 5 ; 10 fs
>> rlist = 3.0 ; short-range neighborlist cutoff (in nm)
>> rcoulomb = 3.0 ; short-range electrostatic cutoff (in nm)
>> rvdw = 3.0 ; short-range van der Waals cutoff (in nm)
>> ; Electrostatics
>> coulombtype = PME ; Particle Mesh Ewald for long-range
>> electrostatics
>> pme_order = 4 ; cubic interpolation
>> fourierspacing = 0.16 ; grid spacing for FFT
>> ; Temperature coupling is on
>> tcoupl = berendsen ; modified Berendsen thermostat
>> tc-grps = Protein Non-Protein ; two coupling groups - more
>> accurate
>> tau_t = 0.1 0.1 ; time constant, in ps
>> ref_t = 300 300 ; reference temperature, one
>> for each group, in K
>> ; Pressure coupling is off
>> pcoupl = no ; no pressure coupling in NVT
>> ; Periodic boundary conditions
>> pbc = xyz ; 3-D PBC
>> ; Dispersion correction
>> DispCorr = EnerPres ; account for cut-off vdW scheme
>> ; Velocity generation
>> gen_vel = yes ; assign velocities from Maxwell
>> distribution
>> gen_temp = 300 ; temperature for Maxwell distribution
>> gen_seed = -1 ; generate a random seed
>>
>>
>>
>> 2012/3/29 Rausch, Felix <frausch at ipb-halle.de>:
>>
>> [Hide Quoted Text]
>> Hello again,
>>
>> Well, I once had problems with simulations crashing randomly during
>> production runs (sometimes after tens of nanoseconds) with the LINCS
>> warnings you described. Switching LINCS from "all-bonds" to only "h-bonds"
>> did the trick for me, although I never exactly figured out why.
>> Maybe it's worth a try for you, too?
>> Cheers,
>> Felix
>>
>> -----Ursprüngliche Nachricht-----
>> Von: gmx-users-bounces at gromacs.org [mailto:gmx-users-bounces at gromacs.org] Im
>> Auftrag von francesca vitalini
>> Gesendet: Donnerstag, 29. März 2012 15:15
>> An: Discussion list for GROMACS users
>> Betreff: Re: [gmx-users] Not able to continue with Equilibration
>>
>> Hi!
>> I'm having a similar problem. I have a dimer solvated in a big box of water
>> plus ions that I have managed to minimize correctly (see output of
>> minimization at the end) but when I try to run NVT equilibration (see later)
>> I get LINCS warnings(see below) refearred to atoms which are not in a
>> cluster or in a strange position. I have added dihedral restraints on them
>> but still the same type of error. I'm using GROMACS 3.3.1. I have tried to
>> switch to a newer version of GROMACS but still the same error.
>> Does anyone have any suggestions?
>> Thanks
>> Francesca
>>
>>
>> MINIMIZATION OUTPUT
>>
>> Steepest Descents converged to Fmax < 1000 in 681 steps Potential Energy =
>> -1.9597512e+07
>> Maximum force = 2.4159846e+02 on atom 13087
>> Norm of force = 2.1520395e+04
>>
>>
>> MINIMIZATION.MDP
>>
>> define = -DEFLEXIBLE ; flexible water
>> integrator = steep ; Algorithm (steep = steepest descent
>> minimization)
>> emtol = 1000.0 ; Stop minimization when the maximum
>> force < 1000.0 kJ/mol/nm
>> emstep = 0.001 ; Energy step size
>> nsteps = 50000 ; Maximum number of (minimization)
>> steps to perform
>>
>> ; Parameters describing how to find the neighbors of each atom and how to
>> calculate the interactions
>> nstlist = 1 ; Frequency to update the neighbor
>> list and long range forces
>> ns_type = grid ; Method to determine neighbor list
>> (simple, grid)
>> rlist = 1.0 ; Cut-off for making neighbor list
>> (short range forces)
>> coulombtype = PME ; Treatment of long range
>> electrostatic interactions
>> rcoulomb = 1.0 ; Short-range electrostatic cut-off
>> rvdw = 1.0 ; Short-range Van der Waals cut-off
>> pbc = xyz ; Periodic Boundary Conditions (yes/no)
>>
>>
>> NVT.MDP
>>
>> define = -DPOSRES ; position restrain the protein
>> ; Run parameters
>> integrator = md ; leap-frog integrator
>> nsteps = 10000 ; 2 * 10000 = 20 ps
>> dt = 0.002 ; 2 fs
>> ; Output control
>> nstxout = 50 ; save coordinates every 0.2 ps
>> nstvout = 50 ; save velocities every 0.2 ps
>> nstenergy = 50 ; save energies every 0.2 ps
>> nstlog = 50 ; update log file every 0.2 ps
>> ; Bond parameters
>> unconstrained_start = no
>> ;continuation = no ; first dynamics run
>> constraint_algorithm = lincs ; holonomic constraints
>> constraints = all-bonds ; all bonds (even heavy atom-H bonds)
>> constrained
>> lincs_iter = 1 ; accuracy of LINCS
>> lincs_order = 4 ; also related to accuracy
>> ; Neighborsearching
>> ns_type = grid ; search neighboring grid cells
>> nstlist = 5 ; 10 fs
>> rlist = 3.0 ; short-range neighborlist cutoff (in nm)
>> rcoulomb = 3.0 ; short-range electrostatic cutoff (in nm)
>> rvdw = 3.0 ; short-range van der Waals cutoff (in nm)
>> ; Electrostatics
>> coulombtype = PME ; Particle Mesh Ewald for long-range
>> electrostatics
>> pme_order = 4 ; cubic interpolation
>> fourierspacing = 0.16 ; grid spacing for FFT
>> ; Temperature coupling is on
>> tcoupl = berendsen ; modified Berendsen thermostat
>> tc-grps = Protein Non-Protein ; two coupling groups - more
>> accurate
>> tau_t = 0.1 0.1 ; time constant, in ps
>> ref_t = 300 300 ; reference temperature, one
>> for each group, in K
>> ; Pressure coupling is off
>> pcoupl = no ; no pressure coupling in NVT
>> ; Periodic boundary conditions
>> pbc = xyz ; 3-D PBC
>> ; Dispersion correction
>> DispCorr = EnerPres ; account for cut-off vdW scheme
>> ; Velocity generation
>> gen_vel = yes ; assign velocities from Maxwell
>> distribution
>> gen_temp = 300 ; temperature for Maxwell distribution
>> gen_seed = -1 ; generate a random seed
>>
>> ; Dihedral restraints
>>
>> dihre = simple ; Some dihedrals are restrained for
>> instance peptide bonds are set to trans conformation.
>> dihre_fc = 500
>> dihre_tau = 0.0
>> nstdihreout = 50
>>
>> NVT LINCS ERROR
>>
>> Step 0, time 0 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 0.002851 (between atoms 10721 and 10723) rms 0.000161 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 672 673 33.7 0.1000 0.1000 0.1000
>> 1081 1082 33.8 0.1000 0.1000 0.1000
>> 1647 1648 34.4 0.1000 0.1000 0.1000
>> 2819 2820 37.1 0.1000 0.1000 0.1000
>> 2920 2921 39.6 0.1000 0.1000 0.1000
>> 4090 4091 35.2 0.1000 0.1000 0.1000
>> 4727 4728 37.2 0.1000 0.1000 0.1000
>> 5160 5161 31.7 0.1000 0.1000 0.1000
>> 6824 6825 33.3 0.1000 0.1000 0.1000
>> ...
>> step 0
>> Step 1, time 0.002 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 0.017976 (between atoms 14692 and 14693) rms 0.000520 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 672 673 67.8 0.1000 0.0999 0.1000
>> ....
>>
>> Step 2, time 0.004 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 0.060681 (between atoms 6906 and 6907) rms 0.001874 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 672 673 80.8 0.0999 0.0999 0.1000
>> 1081 1082 90.0 0.1000 0.1008 0.1000
>> ...
>>
>> Step 3, time 0.006 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 0.104238 (between atoms 10603 and 10604) rms 0.003281 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 672 673 67.1 0.0999 0.1000 0.1000
>> 1081 1082 69.6 0.1008 0.1000 0.1000
>> ...
>>
>> Step 4, time 0.008 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 0.569872 (between atoms 8709 and 8710) rms 0.013041 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 672 673 40.3 0.1000 0.1000 0.1000
>> 1081 1082 39.6 0.1000 0.1000 0.1000
>> 1249 1250 37.6 0.1000 0.1000 0.1000
>> ...
>>
>> Step 5, time 0.01 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 2.761976 (between atoms 8202 and 8204) rms 0.056483 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 3320 3321 43.0 0.1090 0.1090 0.1090
>> 3875 3876 32.0 0.1090 0.1090 0.1090
>> 4013 4014 34.3 0.1090 0.1082 0.1090
>> 4070 4071 42.1 0.1092 0.1092 0.1090
>> 4072 4073 45.2 0.1090 0.1086 0.1090
>> 4727 4729 34.2 0.1007 0.1002 0.1000
>> 4961 4962 31.4 0.1091 0.1089 0.1090
>> 4992 4993 46.0 0.1091 0.1089 0.1090
>> 8159 8160 30.9 0.1251 0.1465 0.1230
>> ...
>>
>> Step 6, time 0.012 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 22.009092 (between atoms 8202 and 8204) rms 0.488487 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 2920 2921 44.7 0.1000 0.1001 0.1000
>> 3320 3321 36.6 0.1090 0.1087 0.1090
>> 3592 3593 32.1 0.1090 0.1089 0.1090
>> ....
>> Warning: 1-4 interaction between 8159 and 8168 at distance 1.580 which is
>> larger than the 1-4 table size 1.000 nm These are ignored for the rest of
>> the simulation This usually means your system is exploding, if not, you
>> should increase table-extension in your mdp file
>>
>> Step 7, time 0.014 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 403.393890 (between atoms 8202 and 8204) rms 8.430617 bonds that rotated
>> more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 302 303 31.7 0.1003 0.1001 0.1000
>> 672 673 66.4 0.1001 0.1001 0.1000
>> 1050 1051 33.1 0.1003 0.1003 0.1000
>> 1081 1082 46.7 0.1000 0.1000 0.1000
>> 1225 1226 35.0 0.1003 0.1003 0.1000
>> ....
>>
>> Wrote pdb files with previous and current coordinates
>>
>> Step 8, time 0.016 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 7568.298828 (between atoms 8202 and 8204) rms 155.435715 bonds that
>> rotated more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length
>> 205 206 31.6 0.1003 0.1002 0.1000
>> 302 303 36.9 0.1001 0.1002 0.1000
>> 672 673 76.0 0.1001 0.0997 0.1000
>> 1050 1051 34.3 0.1003 0.1003 0.1000
>> 1081 1082 63.0 0.1000 0.0999 0.1000
>> 1225 1226 36.6 0.1003 0.1003 0.1000
>> 1249 1250 38.7 0.1001 0.1001 0.1000
>> 1647 1648 62.8 0.1000 0.0996 0.1000
>> 1679 1680 42.5 0.1000 0.1000 0.1000
>> 2214 2215 31.6 0.1002 0.1001 0.1000
>> 2577 2578 31.1 0.1003 0.1002 0.1000
>> 2819 2820 68.0 0.1000 0.0997 0.1000
>> 2920 2921 71.5 0.1023 0.1000 0.1000
>> 3468 3469 34.7 0.1092 0.1090 0.1090
>> 3676 3677 34.7 0.1090 0.1090 0.1090
>> 3714 3717 34.7 0.0998 0.0999 0.1000
>> 3746 3747 70.8 0.1086 0.1069 0.1090
>> ....
>>
>> Back Off! I just backed up step7.pdb to ./#step7.pdb.1# Wrote pdb files with
>> previous and current coordinates
>>
>> Step 9, time 0.018 (ps) LINCS WARNING
>> relative constraint deviation after LINCS:
>> max 148599.687500 (between atoms 8202 and 8204) rms 3002.328125 bonds that
>> rotated more than 30 degrees:
>> atom 1 atom 2 angle previous, current, constraint length ....
>>
>> Back Off! I just backed up step8.pdb to ./#step8.pdb.1# Wrote pdb files with
>> previous and current coordinates Segmentation fault
>>
>>
>> 2012/3/29 Mark Abraham <Mark.Abraham at anu.edu.au>:
>> On 29/03/2012 8:22 PM, Hendry wrote:
>>
>> Hi,
>>
>>
>>
>> I am using Gromacs 4.5.4. After successful minimization by SD, I
>> continued with equilibration step but I got the below errors. I tried
>> many times with different parameters but the problem still persists. I
>> have given errors and md parameters of equilibration step below. I
>> have also provided my minimization output at the end. Could you
>> provide some suggestions what went wrong?.
>>
>>
>> You are using PR-pressure coupling for equilibration, which is an
>> unstable combination. You are coupling ions to their own thermostat,
>> which is unstable. Do check out the GROMACS manual for discussion of
>> how these algorithms work, and also the website for some practical
>> observations.
>>
>> Mark
>>
>> -
>>
>>
>> --
>> gmx-users mailing list gmx-users at gromacs.org
>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>> Please search the archive at
>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>> Please don't post (un)subscribe requests to the list. Use thewww interface
>> or send it to gmx-users-request at gromacs.org.
>> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
>
>
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
More information about the gromacs.org_gmx-users
mailing list