[gmx-users] Fwd: Help with Gromacs 5.1 mdp options for CHARMM27 force field

[Justin Lemkul]

On 12/6/15 8:57 AM, Rakesh Ramachandran wrote:
> Hi Justin,
>
>      Thanks for the reply. Sorry if I am asking too many questions. I wanted
> to clarify some more doubts just to make sure I am performing the
> simulation the right way.
>
> 1) What is the ideal Tcoupl option for CHARMM36 Nose-hoover or V-rescale.
> Moreover I get the following note and warning when I use Nose-Hoover
> thermostat.
>

Use either; there is no force field dependence here.

> NOTE 1 [file npt.mdp]:
>    leapfrog does not yet support Nose-Hoover chains, nhchainlength reset to 1
>

Notes are for information only; they do not indicate a problem.

> WARNING 1 [file npt.mdp]:
>    For proper integration of the Nose-Hoover thermostat, tau-t (1) should be
>    at least 20 times larger than nsttcouple*dt (0.08)
>

Check your settings here.  Warnings are important.

> So is V-rescale more optimal for Verlet cutoff and GPU based Gromacs.
>
> 2) Moreover I see in some publications that 0.9 nm to 1.0 nm for minimum
> protein-edge distance is being used for Gromacs. Even in your tutorial you
> mention using 1.0 nm. Taking into account minimum image convention should
> it not be greater than 1.2 nm or is it highly dependent on the protein you
> are simulating or rcouloumb cutoff or is dependent on Debye length.
>

As long as the minimum image convention is not violated, there is no problem. 
Using a minimum solute-box distance equal to the longest cutoff is safe, but 
overkill.  1.0 nm is typically sufficient because then there are at least 2.0 nm 
between image atoms of, e.g. a protein, which is a couple of solvent shells and 
therefore enough.

> I would be grateful if you could clear my doubts on this aspect or point me
> to the relevant literature. I am planning to use 1.25 nm for a salt conc of
> 300 mM MgCl2 with epsilon_r = 1 or should I use higher cutoff such as 1.4
> nm. Moreover is the dielectric constant of the medium 82 for TIP3P water
> model at 300K (ref http://www1.lsbu.ac.uk/water/water_models.html),
>

Do not touch epsilon_r.  Leave it at 1.  It is not the dielectric constant; it 
is the reference for dielectric screening.  If you change this value, you change 
in vacuo screening in a manner inconsistent with real physics and the force 
field parametrization.

> 3) In spite of running the NPT equilibration for 5 ns I am getting these
> values and the RMSD values are also too high when the pressure parameter is
> measured. In your tutorial you had suggested to run longer NPT
> equilibration if the density value does not converge to around 1000 kg/m^3
> and pressure value around 1 bar. I found the same issue when I ran with
> CHARMM27 force field and V-rescale thermostat. But on plotting I find the
> results are similar to your plots in the tutorial but for 5 ns in both the
> cases (CHARMM 27 and 36).
>
> for 1.25 nm protein-edge distance
>
> Pressure                    1.52207       0.52    146.108    1.62992  (bar)
> Density                     1056.64      0.029     2.7515 -0.0336808
>   (kg/m^3)
>
> for 1.4 nm protein-edge distance
>
> Pressure                   0.321022       0.94    137.028     2.4858  (bar)
> Density                     1054.45      0.053    2.67943   0.264365
>   (kg/m^3)
>

Welcome to life when trying to simulate pressure.

http://www.gromacs.org/Documentation/Terminology/Pressure

-Justin

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

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

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

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

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