[gmx-users] Problem with energy minimization parameters

[Justin Lemkul]

On 7/16/13 11:05 AM, Sainitin Donakonda wrote:
> Hi all,
>
> I have 3 protein-drug complexes (same protein and different drugs)  i
> performed 20 ns molecular dynamics simulation..using same parameters (i
> mean energy minimzation parameters) ..2 of them was successful for MD
> simulation..(im using Charmm 27 Force field)
>
> But 1 drug protein complex is failed in Energy minimization step itself..
> showing following error
>
> *Error : whole protein drug complex
> *
> *
> *
> *708 particles communicated to PME node 2 are more than 2/3 times the
> cut-off out of the domain decomposition cell of their charge group in
> dimension x.*
> *This usually means that your system is not well equilibrated
>
> *
> Then i looked up http://www.gromacs.org/Documentation/Terminology/Blowing_Up
>
> In this i followed up suggestions which are mentioned in diagnosing
> unstable system as follows.
>
> First i did energy minimization for protein alone which worked perfectly
> with out any error.
>
> Second i did energy minimization for ligand alone (with out protein)..this
> showed again following error..
>
> Error 2 : Ligand only
>
> Fatal error:
> *184 particles communicated to PME node 2 are more than 2/3 times the
> cut-off out of the domain decomposition cell of their charge group in
> dimension x.
> This usually means that your system is not well equilibrated.*
>
> So there is problem with small molecule...not with the protein ..here i
> think i should tune parameters in .MDP file..which i m just not getting how
> to do..
>
> here is the EM.mdp file..
>
> ; LINES STARTING WITH ';' ARE COMMENTS
> title        = Minimization    ; Title of run
> define                  = -DFLEX_TI3P    ; defines to pass to the
> preprocessor
> constraints             = none
>
> ; 7.3.3 Run Control
> integrator              = steep         ; steepest descents energy
> minimization
> nsteps                  = 2000          ; maximum number of steps to
> integrate
> energygrps              = LIG        ; group(s) to write to energy file
>
> ; 7.3.5 Energy Minimization
> emtol                   = 1000          ; [kJ/mol/nm] minimization is
> converged when max force is < emtol
> emstep                  = 0.01          ; [nm] initial step-size
> dt                      = 0.01
>
>
> ; 7.3.9 Neighbor Searching
> nstlist                 = 1            ; [steps] freq to update neighbor
> list
> nstxout                 = 1
> ns_type                 = grid          ; method of updating neighbor list
> pbc                     = xyz           ; periodic boundary conditions in
> all directions
> rlist                   = 1.2           ; [nm] cut-off distance for the
> short-range neighbor list
> rlistlong               = 1.4           ; [nm] Cut-off distance from the
> long-range neighbor list
>
> ; 7.3.10 Electrostatics
> coulombtype             = PME           ; Particle-Mesh Ewald electrostatics
> rcoulomb                = 1.2           ; [nm] distance for Coulomb cut-off
>
> ; 7.3.11 VdW
> vdwtype                 = switch       ; twin-range cut-off with rlist
> where rvdw >= rlist
> rvdw                    = 1.2           ; [nm] distance for LJ cut-off
> rvdw_switch             = 0.2           ; Start switching th LJ potential
> DispCorr                = Ener          ; apply long range dispersion
> corrections for energy
>
> ; 7.3.13 Ewald
> fourierspacing          = 0.12          ; [nm] grid spacing for FFT grid
> when using PME
> fourier_nx              = 0
> fourier_ny              = 0
> fourier_nz              = 0
> pme_order               = 4             ; interpolation order for PME, 4 =
> cubic
> ewald_rtol              = 1e-5          ; relative strength of
> Ewald-shifted potential at rcoulomb
> optimize_fft            = yes
>
>
> Can any body please tell me how to optimize these parameters
>

The value of rvdw_switch should be 1.0 for correct usage of the force field. 
One should not seek to "optimize" the .mdp settings when something is crashing. 
  Fixing errors (like rvdw_switch here) is one thing, but beyond that you 
shouldn't do any adjusting.  More often than not, in a case like yours where a 
molecule in vacuo blows up, the topology is the problem.

-Justin

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

Justin A. Lemkul, Ph.D.
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

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