[gmx-users] Box size problem

[Justin Lemkul]

On 3/15/17 3:09 PM, Davit Hakobyan wrote:
> Dear All,
>
> I try to run a simulation with protein-membrane system in charmm36 ff with
> Gromacs5.0.4 and the simulation runs without complain.The problem is in the unit
> cell size which seems to grow during a short 500 ps simulation from ~15 to ~28
> nm in Z direction although the system atoms remain confined in the~15 nm in Z
> direction.
>
>
> The GRO files before and after the 500 ps simulation can be downloaded from the
> following link: https://uni-muenster.sciebo.de/index.php/s/GBHjED2t1KzZZ2o
>

Your initial coordinates are badly fragmented, with regions of vacuum and 
isolated water molecules.  You're going to get bizarre behavior as the system 
collapses to try to close those voids.  This outcome likely isn't a function of 
.mdp settings, but unphysical initial coordinates.

<snip>

Bonus comments!

> MDP Parameters (OLD)
> -------------------------------------------------
> integrator          = md
> dt                  = 0.002
> nsteps              = 250000000
>
> nstxout             = 500000
> nstvout             = 500000
> nstlog              = 100000
> nstenergy           = 500000
> nstxtcout           = 500000
> xtc_grps            =
>
> energygrps          = Protein CHL1 POPC DOPC POPS Water_and_ions
> ;energygrp_table    =
> nstcalcenergy       = 500
> nstlist             = 10
> nstcomm             = 50
> comm_mode           = Linear
> comm-grps           =
> ns_type             = grid
> rlist               = 1.0
> rlistlong           = 1.4
> rvdw_switch         = 0.8
> vdwtype             = Switch

This is a potential switch, not a force switch, so lipid simulations done with 
these settings will have artificial behavior.  Force switching is standard for 
all CHARMM additive force field simulations, and is absolutely required for lipids.

> coulombtype         = pme
> rcoulomb            = 1.0
> rcoulomb_switch     = 0.0
> rvdw                = 1.2
> fourierspacing      = 0.15
> pme_order           = 6
> ;ewald_rtol         = 1e-6
> tcoupl              = V-rescale ;nose-hoover
> nhchainlength       = 1
> tc-grps             = Protein CHL1 POPC DOPC POPS Water_and_ions

Coupling all different types of molecules in the system to their own thermostats 
is inappropriate.

-Justin

> tau_t               = 0.1 0.1 0.1 0.1 0.1 0.1
> ref_t               = 293 293 293 293 293 293
> Pcoupl              = parrinello-rahman ;berendsen ;parrinello-rahman
> Pcoupltype          = semiisotropic
> tau_p               = 5.0
> compressibility     = 4.5e-5 4.5e-5
> ref_p               = 1.0 1.0
> pbc                 = xyz
> gen_vel             = no ;yes
> gen_temp            = 293
> optimize_fft        = yes
> constraints         = hbonds
> continuation        = no
> constraint_algorithm = Lincs
> lincs-order         = 4 ; 8 is needed for BD with large time-steps.
> lincs-iter          = 1 ; 1 is fine for normal simulations, but use 2 to
> conserve energy in NVE runs.
> -------------------------------------------------
>
> In the new parameters I also tried to change the tcoupl from Nose-Hoover to
> V-rescale, but the outcome is the same, the box size becomes large.
>
> ProbablyI am missing something trivial.
> Does anybody have a clue what might be wrong?
>
> Thanks very much for any help.
> Davit

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

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

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