[gmx-users] A question about "blowing up"

[Justin Lemkul]

On 3/1/14, 7:22 AM, sujith ks wrote:
> Hello Justin,
>
>    I performed energy minimization , the way you performed for the free
> energy calculation for the methane-water system in the tutorial. In the
> first step steepest descent was used followed by l-bfgs. The emtol (force)
> was set to 10 kj/mol/nm.
>
> The final results can be summerized as:
>
> Low-Memory BFGS Minimizer converged to Fmax < 10 in 7929 steps
> Potential Energy  = -9.4460383e+04
> Maximum force     =  7.4954224e+00 on atom 1760
> Norm of force     =  5.8917433e-01
>
> The potential vs time graph looks good (converged).
>
>   The .mdp is as follows:
>
> ; Run control
> integrator               = l-bfgs
> nsteps                   = 100000
> define                    = -DFLEXIBLE
> ; EM criteria and other stuff
> emtol                      = 10
> emstep                   = 0.01
> niter                        = 20
> nbfgscorr                 = 10
> ; Output control
> nstlog                      = 1
> nstenergy                = 1
> ; Neighborsearching and short-range nonbonded interactions
> nstlist                     = 1
> ns_type                  = grid
> pbc                        = xyz
> rlist                        = 0.9
> ; Electrostatics
> coulombtype              = pme
> rcoulomb                   = 0.9
> ; van der Waals
> vdw-type                 = Cut-off
> rvdw                       = 1.4
> ; Apply long range dispersion corrections for Energy and Pressure
> DispCorr                   = EnerPres
> ; Spacing for the PME/PPPM FFT grid
> fourierspacing            = 0.12
> ; EWALD/PME/PPPM parameters
> pme_order                 = 6
> ewald_rtol                 = 1e-06
> epsilon_surface          = 0
> optimize_fft               = no
> ; Temperature and pressure coupling are off during EM
> tcoupl                    = no
> pcoupl                   = no
> ; Generate velocities to start
> gen_vel                  = no
> ; options for bonds
> constraints              = none
> ; Type of constraint algorithm
> constraint-algorithm     = lincs
> ; Do not constrain the starting configuration
> continuation             = no
> ; Highest order in the expansion of the constraint coupling matrix
> lincs-order               = 12

This explains why turning off constraints solved the problem.  You minimized 
without any constraints on anything (-DFLEXIBLE and constraints = none), but 
then turned them on during NVT.  What can happen is that a configuration that is 
unconstrained minimizes happily, but the geometry is distorted.  You need to run 
EM under the same conditions that you will use for the MD (in terms of 
constraints).  Minimizing without constraints initially is fine if there are 
problems in the structure, but if you intend to constraint anything later, you 
need to minimize with constraints before you can do so.

-Justin

-- 
==================================================

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 601
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

==================================================