[gmx-users] Drude force field

[Justin Lemkul]

On 7/20/17 8:08 AM, Shraddha Parate wrote:
> Hi gmx users!
> 
> For polarizable effect, I was setting drude force field in Gromacs
> I attempted running Lagrangian Dynamics simulation using Drude force field.
> mdp file format was used from the supporting data in below paper.
> 
> [2015] Implementation of extended Lagrangian dynamics in GROMACS for
> polarizable simulations using the classical Drude oscillator model.(
> http://onlinelibrary.wiley.com/doi/10.1002/jcc.23937/abstract;jsessionid=E354008D7BE6AA1A02D26C061B49163A.f02t01
> )
> 
> 
> 
> 
> -------------------------------------------------------
> 
>    ; RUN CONTROL PARAMETERS
>    integrator = md-vv
> 
>    ; Start time and timestep in ps
>    tinit = 0
>    dt = 0.001
>    nsteps = 10000
>    comm-mode = linear
>    nstcomm = 100
>    comm-grps = System
> 
>    ; OUTPUT CONTROL OPTIONS
>    nstxout = 1000
>    nstvout = 1000
>    nstfout = 1000
> 
>    ; Output frequency for energies to log file and energy file
>    nstlog = 100
>    nstcalcenergy = 100
>    nstenergy = 100
> 
>    ; NEIGHBORSEARCHING PARAMETERS
>    cutoff-scheme = verlet
>    nstlist = 10
>    ns-type = Grid
>    pbc = xyz
>    rlist = 1.2
> 
>    ; OPTIONS FOR ELECTROSTATICS AND VDW
>    coulombtype = PME
>    rcoulomb = 1.2
>    vdwtype = cutoff
>    vdw-modifier = potential-switch
>    rvdw-switch = 1.0
>    rvdw = 1.2
>    DispCorr = EnerPres
> 
>    ; TEMPERATURE COUPLING
>    tcoupl = Nose-Hoover
>    nsttcouple = 1
>    nh-chain-length = 1
> 
>    ; Groups to couple separately
>    tc-grps = Atoms Drues
>    tau-t = 0.1 0.005
>    ref-t = 298.15 1.0
> 
>    ; GENERATE VELOCITIES FOR STARTUP RUN, OTHERWISE DO
>    ; NOT GENERATE NEW VELOCITIES IF PREVIOUSLY EQUILIBRATED
>    gen-vel = yes
>    gen-temp = 298.15
>    gen-seed = -1
> 
>    ; OPTIONS FOR BONDS
>    constraints = none ; can also be h-bonds, not all-bonds
>    continuation = no
> 
> 
>    ; DRUDE PARAMETERS
>    drude = yes
>    drude-mode = Lagrangian
>    drude-t = 1.0
>    drude-hardwall = yes
>    drude-r = 0.02
>    drude-tsteps = 5
> -------------------------------------------------------
> 
> 
> I couldn't understand [tc-grps = Atoms Drude] option.

All the atoms in one group, all the Drudes in the next.  If you want to break 
the system apart into protein and non-protein:

tc-grps = Protein-Atoms Protein Drudes Non-Protein-Atoms Non-Protein-Drudes
ref_t = 298.15 1.0 298.15 1.0
tau_t = 1.0 0.005 1.0 0.005

> If the option is changed to [tc-grps = Protein SOL],I get a fatal error
> message.
> 
> ------------------------------------------------------
> Fatal error:
> Incorrect thermostat setup in relative_tstat, ti = -1
> ------------------------------------------------------
> 
> How to resolve this issue? Why the changed parameter gives a fatal error?
> 

Do what the paper (and the ones cited therein) tell you to do.  The Drude 
integration scheme requires a special dual thermostat that couples the Drude 
oscillators to a low-temperature thermostat to approximate the SCF surface.  You 
can't treat these kinds of systems like you would additive ones.  The atoms 
(anything non-Drude) and Drudes have to be in separate thermostats, in a given 
order, as explained in the paper (and shown above).

Note that domain decomposition remains broken, so you can only run with OpenMP. 
I recently identified the (hopefully) last remaining issue but it will take some 
time to get it fixed as I am currently moving and starting my own lab, so my 
attention is elsewhere for a while.

-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

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