[gmx-users] Charmm 36 forcefield with verlet cut-off scheme

[Justin Lemkul]

On 9/19/13 2:21 PM, akk5r wrote:
> Hey All,
>
> I am running a simulation of a POPC bilayer in water using a Charmm 36
> forcefield. I am having trouble using cutoff-scheme=Verlet in my nvt
> equilibration.
>
> My .mdp file is as follows:
>
> title	 = OPLS Lysozyme NVT equilibration
> define	 = -DPOSRES	; position restrain the protein
> ; Run parameters
> integrator	= md	 ; leap-frog integrator
> nsteps	 = 50000	 ; 2 * 50000 = 100 ps
> dt	 = 0.002	 ; 2 fs
> ; Output control
> nstxout	 = 0	 ; save coordinates every 0.2 ps
> nstvout	 = 0	 ; save velocities every 0.2 ps
> nstenergy	= 0	 ; save energies every 0.2 ps
> nstlog	 = 0	 ; update log file every 0.2 ps
> ; Bond parameters
> continuation	= no	 ; first dynamics run
> constraint_algorithm = lincs	; holonomic constraints
> constraints	= all-bonds	; all bonds (even heavy atom-H bonds) constrained
> lincs_iter	= 1	 ; accuracy of LINCS
> lincs_order	= 4	 ; also related to accuracy
> ; Neighborsearching
> ns_type	 = grid	 ; search neighboring grid cells
> nstlist	 = 5	 ; 10 fs
> rlist	 = 1.2	 ; short-range neighborlist cutoff (in nm)
> rcoulomb	= 1.2	 ; short-range electrostatic cutoff (in nm)
> rvdw	 = 1.2	 ; short-range van der Waals cutoff (in nm)
> rcoulomb-switch          = 0
> rvdw-switch              = 1.0
> ; long-range cut-off for switched potentials
> rlistlong                = 1.4
> ; Electrostatics
> coulombtype	= PME	 ; Particle Mesh Ewald for long-range electrostatics
> pme_order	= 4	 ; cubic interpolation
> fourierspacing	= 0.16	 ; grid spacing for FFT
> ; Relative dielectric constant for the medium and the reaction field
> epsilon_r                = 1
> epsilon_rf               = 1
> ; Temperature coupling is on
> tcoupl	 = V-rescale	; modified Berendsen thermostat
> tc-grps	 = POPC SOL	; two coupling groups - more accurate
> tau_t	 = 0.1	0.1	; time constant, in ps
> ref_t	 = 300 300	; reference temperature, one for each group, in K
> ; Pressure coupling is off
> pcoupl	 = no ; no pressure coupling in NVT
> ; Periodic boundary conditions
> pbc	 = xyz	 ; 3-D PBC
> ; Dispersion correction
> DispCorr	= No	; account for cut-off vdW scheme
> ; Method for doing Van der Waals
> vdw-type                 = switch
> ; Velocity generation
> gen_vel	 = yes	 ; assign velocities from Maxwell distribution
> gen_temp	= 300	 ; temperature for Maxwell distribution
> gen_seed	= -1	 ; generate a random seed
> cutoff-scheme   = Verlet
>
> *Here is the output when I run grompp to make the .tpr file:*
>
> grompp -f nvt.mdp -c em_POPC.gro -p topol.top -o nvt.tpr
>
> ERROR 1 [file nvt.mdp]:
>    With Verlet lists only cut-off LJ interactions are supported
>
>
> NOTE 1 [file nvt.mdp]:
>    With Verlet lists the optimal nstlist is >= 10, with GPUs >= 20. Note
>    that with the Verlet scheme, nstlist has no effect on the accuracy of
>    your simulation.
>
>
> NOTE 2 [file nvt.mdp]:
>    The switch/shift interaction settings are just for compatibility; you
>    will get better performance from applying potential modifiers to your
>    interactions!
>
>
> Generated 21528 of the 21528 non-bonded parameter combinations
> Generating 1-4 interactions: fudge = 1
> Generated 18355 of the 21528 1-4 parameter combinations
> Excluding 3 bonded neighbours molecule type 'POPC'
> turning all bonds into constraints...
> Excluding 2 bonded neighbours molecule type 'SOL'
> turning all bonds into constraints...
> Setting gen_seed to 7358
> Velocities were taken from a Maxwell distribution at 300 K
> Removing all charge groups because cutoff-scheme=Verlet
>
> There were 2 notes
>
> -------------------------------------------------------
> Program grompp, VERSION 4.6.1
> Source code file: /home/ali/Downloads/gromacs-4.6.1/src/kernel/grompp.c,
> line: 1593
>
> Fatal error:
> There was 1 error in input file(s)
> For more information and tips for troubleshooting, please check the GROMACS
> website at http://www.gromacs.org/Documentation/Errors
> -------------------------------------------------------
>

As the error message says, with "cutoff-scheme = Verlet" you must use "vdwtype = 
cutoff."  The notes above suggest reasonable methods for implementing switched 
interactions.

-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

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