[gmx-users] Error in MD simulation

Mark Abraham mark.abraham at anu.edu.au
Mon Jun 28 00:54:11 CEST 2010



----- Original Message -----
From: teklebrh at ualberta.ca
Date: Monday, June 28, 2010 6:01
Subject: [gmx-users] Error in MD simulation
To: jalemkul at vt.edu, Discussion list for GROMACS users <gmx-users at gromacs.org>

> Dear Gromacs users,
> 
> First I have to mention few points. I read all the LINCS error 
> and tried a number of options.
> 
> 1, Minimized well the system. and the Epot and Max Forces are 
> very resonable.
> 
> Steepest Descents converged to Fmax < 100 in 379 steps
> Potential Energy  = -1.0573750e+05
> Maximum force     =  8.0876366e+01 on 
> atom 1117
> Norm of force     =  8.7360611e+00
> 
> 2, Equlibrate well the system. It works for 100ps
> 
> Once I run NPT simulation using the following .mdp parameter set 
> it shows the following error.
> 
> ; The etire lines starting with ';' ARE COMMENTS
> 
> title        	= NPT 
> equilibration using the berendsen thermostat and barostat	; 
> Title of run
> cpp             = /usr/bin/cpp ; location of cpp on linux
> 
> ; The following lines tell the program the standard locations 
> where to find certain files
> 
> ; VARIOUS PREPROCESSING OPTIONS
> ; Preprocessor information: use cpp syntax.
> ; e.g.: -I/home/joe/doe -I/home/mary/hoe
> include                  =
> 
> ; e.g.: -DI_Want_Cookies -DMe_Too
> 
> ; RUN CONTROL PARAMETERS
> 
> integrator               = md ; leap-frog integrator
> ; Start time and timestep in ns
> tinit                    = 0
> dt                       = 0.002     ; 2fs
> nsteps                   = 500000   ; 2*500 = 1ns
> 
> ; THOSE OPTIONS REMOVE MOTION OF THE ASPHALTENE MODEL COMPUNDS 
> RELATIVE TO THE SOLVENT/IONS
> 
> ; mode for center of mass motion removal
> comm-
> mode                = Linear
> ; number of steps for center of mass motion removal
> nstcomm                  = 1
> ; group(s) for center of mass motion removal
> comm-
> grps                = System
> 
> ; LANGEVIN DYNAMICS OPTIONS
> ; Friction coefficient (amu/ps) and random seed
> bd-
> fric                  = 0
> ld-
> seed                  = 1993
> 
> ; ENERGY MINIMIZATION OPTIONS
> ; Force tolerance and initial step-size
> emtol                    = 200
> emstep                   = 0.01
> 
> ; Frequency of steepest descents steps when doing CG
> nstcgsteep               = 1000
> nbfgscorr                = 10
> 
> ; OUTPUT CONTROL OPTIONS
> 
> ; Output frequency for coords (x), velocities (v) and forces (f)
> nstxout                  = 2000
> nstvout                  = 2000
> nstfout                  = 0
> ; Output frequency for energies to log file and energy file
> nstlog                   = 1000
> nstenergy                = 1000
> ; Output frequency and precision for xtc file
> nstxtcout                = 1000
> xtc-
> precision            = 1000
> 
> 
> ; This selects the subset of atoms for the xtc file. You can
> ; select multiple groups. By default all atoms will be written.
> xtc-
> grps                 = BIS  HEP
> ; Selection of energy groups
> energygrps               = BIS  HEP
> 
> ; NEIGHBORSEARCHING PARAMETERS
> 
> ; nblist update frequency
> nstlist                  = 5
> ; ns algorithm (simple or grid)
> ns-
> type                  = grid
> ; Periodic boundary conditions: xyz, no, xy
> pbc                      = xyz
> ; nblist cut-off
> rlist                    = 1.2
> 
> ; OPTIONS FOR ELECTROSTATICS AND VDW
> 
> ; Method for doing electrostatics
> coulombtype              = PME
> rcoulomb                 = 1.2
> 
> 
> ; Method for doing Van der Waals
> vdw-
> type                 = Cut-off
> ; cut-off lengths
> rvdw                     = 1.2
> 
> ; Apply long range dispersion corrections for Energy and Pressure
> DispCorr                 = EnerPres ; account for the cut-off vdw scheme
> 
> ; Spacing for the PME/PPPM FFT grid
> fourierspacing           = 0.16
> 
> ; EWALD/PME/PPPM parameters
> pme_order                = 4
> ewald_rtol               = 1e-05
> ewald_geometry           = 3d
> 
> ; IMPLICIT SOLVENT ALGORITHM
> implicit_solvent         
> = No
> 
> ; OPTIONS FOR WEAK COUPLING ALGORITHMS
> ; Temperature coupling is on
> tcoupl                   = V-rescale; modified berendsen thermostat
> ; Groups to couple separately
> tc-
> grps                  = BIS  HEP ; three coupling gropus
> ; Time constant (ps) and reference temperature (K)
> tau-
> t                    = 0.1  0.1   ; time constant
> ref-
> t                    = 298  298   ; reference temprature
> ; Pressure coupling is on
> Pcoupl                   = berendsen; Pressure coupling is on
> Pcoupltype               = Isotropic
> ; Time constant (ps), compressibility (1/bar) and reference P (bar)
> tau-
> p                    = 3  3
> compressibility          = 1.47e-4  1.47e-4  ; isothermal compressibility, bar^-1
> ref-
> p                    = 1  1
> ; Random seed for Andersen thermostat
> andersen_seed            = 815131
> 
> ; GENERATE VELOCITIES FOR STARTUP RUN
> gen-
> vel                  = yes ; assign velocities from Maxwell distribution (MXD)
> gen-
> temp                 = 298 ; temperature for MXD
> gen-
> seed                 = 173529    ; used to initialize random generator for random velocities
> 
> ; OPTIONS FOR BONDS
> constraints              = all-bonds
> ; Type of constraint algorithm
> constraint-algorithm     = Lincs
> ; Do not constrain the start configuration
> ; Highest order in the expansion of the constraint coupling matrix
> lincs-
> order              = 4
> ; Number of iterations in the final step of LINCS. 1 is fine for
> ; normal simulations, but use 2 to conserve energy in NVE runs.
> ; For energy minimization with constraints it should be 4 to 8.
> lincs-
> iter               = 2
> ; Lincs will write a warning to the stderr if in one step a bond
> ; rotates over more degrees than
> lincs-
> warnangle          = 30
> ; Convert harmonic bonds to morse potentials
> morse                    = no
> 
> ==========
> 
> It run the simulation for quite sometime but it gives me an 
> error of LINCS some how in the middle.
> 
> I try changing some parameter sets like step size,tau-p , tau-t, 
> lincs-warnangle but did not overcame the problem.

You are trying to equilibrate straight into NPT, which can be problematic if (for example) the density isn't right. Consider the advice in point 8 of http://www.gromacs.org/Documentation/How-tos/Steps_to_Perform_a_Simulation

Mark

> Step 197685, LINCS WARNING
> >relative constraint deviation after LINCS:
> >rms 0.000000, max 0.000001 (between atoms 1856 and 15982)
> >bonds that rotated more than 30 degrees:
> > atom 1 atom 2  angle  previous, current, constraint length
> >    97     98   
> 35.3    0.1000   
> 0.1000      0.1000
> 
> can some one give me an advice or help?
> 
> For my other 10 MD set up it works fine but with diffrent solvents.
> 
> rob
> -- 
> 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/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/mailing_lists/users.php
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20100628/87deeb54/attachment.html>


More information about the gromacs.org_gmx-users mailing list