[gmx-users] Implementing the NERD Force Field in GROMACS 2019.3

Alessandra Villa alessandra.villa.biosim at gmail.com
Wed Feb 19 11:27:10 CET 2020


Hi,

On Wed, Feb 19, 2020 at 11:06 AM Robert Cordina <robert.cordina at strath.ac.uk>
wrote:

> Dear Alessandra,
>
> Many thanks for your reply.
> Re comment on nrexcl, I'm using nrexcl = 3 as the Sum paper states that
> "In the NERD force field, atoms/sites within a molecule that do not
> interact by any other intramolecular potential are also allowed to interact
> though the Lennard-Jones potential".  In the paper all bond stretching,
> bond bending and torsion angle parameters are given, so is nrexcl = 3
> correct?  LJ potentials are 1-4 interactions, so this looks the right way
> to do it to me, but I stand to be corrected.
> Re compressibility, that was a typo from my part in the email, it's
> correctly set at 0.00001.
> Re vdwtype, vdw-modifier, rvdw-switch, rvdw... I have read this section
> multiple times but I'm still not sure how to translate the description in
> the Sum paper to these settings.  Any suggestions as to what "Interactions
> were truncated and shifted at rc = 11 Å with energies shifted at this
> distance" would mean for these settings?
>
> cutoff-scheme   = Verlet
> rcoulomb                = 1.1
> coulombtype     = Reaction-Field (I have to use this because I'm using the
> epsilon-rf setting too right?)
> epsilon-rf              = a.bcd
>
> vdwtype = Cut-Off (?)
> vdw-modifier    = Potential-shift-Verlet or Potential-switch (?)
>

The sentence, you have reported, is not clear. It could be that the authors
used
vdwtype = cutoff
vdw-modifier    = Potential-shift   (from mdp parameter of gromacs2020).
rvdw = 1.1
But you should check with them, which treatment they use for the non-bonded
interactions

Best regards
Alessandra



rvdw-switch     = 0 (?)
> rvdw            = 1.1 (will this just be the same as rcoulomb?)
>
>
> Thanks again for your help.
>
> Best regards,
>
> Robert
>
> -----Original Message-----
> From: gromacs.org_gmx-users-bounces at maillist.sys.kth.se <
> gromacs.org_gmx-users-bounces at maillist.sys.kth.se> On Behalf Of
> Alessandra Villa
> Sent: 19 February 2020 08:47
> To: gmx-users at gromacs.org
> Cc: gromacs.org_gmx-users at maillist.sys.kth.se
> Subject: Re: [gmx-users] Implementing the NERD Force Field in GROMACS
> 2019.3
>
> HI,
>
> Below some suggestions assuming that the energy terms are correctly
> implemented (force field)
>
> On Tue, Feb 18, 2020 at 3:13 PM Robert Cordina <
> robert.cordina at strath.ac.uk>
> wrote:
>
> > Hi, I’m trying to use the NERD forcefield (Sum et al. J. Phys. Chem. B
> > 2003, 107, 14443-14451) in GROMACS 2019.3 for a pure triacylglyceride
> > system, however I’m not entirely sure that I’m setting up the
> > equilibration parameters in the mdp file correctly as I’m getting
> > different simulation results from others who have published papers
> > using this force field.  All the bonding and non-bonding interactions
> > have been typed up in the respective force-field files, and assuming
> > those are correct (I’ve checked and re-checked them multiple times),
> > the only thing I can think that is not correct are the mdp parameters.
> >
> > The Sum et al. paper, and another paper (Pizzirusso et al.
> > J.Am.Chem.Soc.2018, 140, 12405−12414) state the following (text from
> > Sum paper in italics, text from Pizzirusso paper in bold, what I’ve
> > put in the mdp file in regular font)
> >
> > Both NVT and NPT calculations were performed using 40 molecules in a
> > cubic box with periodic boundary conditions
> > pbc                         = xyz
> >
> >
> If you do not have a pre-built structure, your equilibration time may be
> longer that in the original paper, the best is to ask the authors the
> coordinate of pre-equilibrated system.
>
>
>
> > Interactions were truncated and shifted at rc = 11 Å with energies
> > shifted at this distance
> > rcoulomb             = 1.1 (am I missing something here?)
> >
> >
> here you miss the vdwtype; vdw-modifier; rvdw-switch¶; rvdw
>
> (see for definition:
>
> http://manual.gromacs.org/documentation/current/user-guide/mdp-options.html?highlight=vdw#Van
> )
>
>
>
> > The system evolved with a leapfrog algorithm using a 2 fs timestep
> > integrator            = md
> > dt                            = 0.002
> >
> > Constant temperature and constant pressure simulations were maintained
> > with the Berendsen’s thermostat and barostat, respectively, by uniform
> > scaling of the atomic velocities (temperature) and by uniform scaling
> > of the atomic positions and box length (pressure) Atomistic MD was
> > performed using GROMACS 4.5 under an NPT ensemble. In this ensemble,
> > temperature and pressure were kept constant using a Berendsen
> > thermostat/barostat.  Temperature and pressure couplings of 1 and
> > 10 ps were used, respectively. Compressibility was fixed at 1 × 10-5
> > bar-1, and anisotropic pressure coupling was used in the MD simulations.
> > tcoupl                   = Berendsen
> > tau-t                      = 1.0
> > tc-grps                  = XXX
> >
>
> system
>
> > ref-t                       = 350
> >
> > pcoupl                  = Berendsen
> > pcoupltype         = anisotropic
> > compressibility  = 0.0001 0.0001 0.0001 0.0001 0.0001 0.0001
> >
>
> compressibility is   10-5 bar-1
>
> tau-p                     = 10
> > ref-p                      = 1.01325 1.01325 1.01325 1.01325 1.01325
> > 1.01325
> >
> > For all the calculations, a reaction-field correction was applied with
> > continuum dielectric εRF to correct for long-range interactions due to
> > electrostatics cutoff-scheme  = Verlet
> > coulombtype     = Reaction-Field
> > epsilon-rf            = a.bcd (I used actual numbers here of course)
> >
> >
> > In addition I am using
> > gen-vel                 = yes
> > gen-temp            = 350
> > gen-seed             = -1
> > continuation      = no
> >
> >
> > Should I be adding anything with respect to vdw?
> >
> >
> > In the forcefield.itp file I’m setting the following
> > nbfunc                  = 1
> > comb-rule           = 2 (Sum paper states that they use Lorentz-Berthelot
> > combining rule for the LJ parameters)
> > gen-pairs             = yes
> > fudgeLJ                = 1
> > fudgeQQ             = 1
> >
> >
> check that  nrexcl in topol.top is in line with the force field.
>
>
> > In ffbonded.itp
> > bondtypes func = 1
> > angletypes func = 1
> > dihedraltypes func = 5
> >
> >
> >
> Best regards
> Alessandra
>
>
>
> >
> > Thanks in advance for your help.
> >
> > Robert
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
>
>
> On Tue, Feb 18, 2020 at 3:13 PM Robert Cordina <
> robert.cordina at strath.ac.uk>
> wrote:
>
> > Hi, I’m trying to use the NERD forcefield (Sum et al. J. Phys. Chem. B
> > 2003, 107, 14443-14451) in GROMACS 2019.3 for a pure triacylglyceride
> > system, however I’m not entirely sure that I’m setting up the
> > equilibration parameters in the mdp file correctly as I’m getting
> > different simulation results from others who have published papers
> > using this force field.  All the bonding and non-bonding interactions
> > have been typed up in the respective force-field files, and assuming
> > those are correct (I’ve checked and re-checked them multiple times),
> > the only thing I can think that is not correct are the mdp parameters.
> >
> > The Sum et al. paper, and another paper (Pizzirusso et al.
> > J.Am.Chem.Soc.2018, 140, 12405−12414) state the following (text from
> > Sum paper in italics, text from Pizzirusso paper in bold, what I’ve
> > put in the mdp file in regular font)
> >
> > Both NVT and NPT calculations were performed using 40 molecules in a
> > cubic box with periodic boundary conditions
> > pbc                         = xyz
> >
> > Interactions were truncated and shifted at rc = 11 Å with energies
> > shifted at this distance
> > rcoulomb             = 1.1 (am I missing something here?)
> >
> > The system evolved with a leapfrog algorithm using a 2 fs timestep
> > integrator            = md
> > dt                            = 0.002
> >
> > Constant temperature and constant pressure simulations were maintained
> > with the Berendsen’s thermostat and barostat, respectively, by uniform
> > scaling of the atomic velocities (temperature) and by uniform scaling
> > of the atomic positions and box length (pressure) Atomistic MD was
> > performed using GROMACS 4.5 under an NPT ensemble. In this ensemble,
> > temperature and pressure were kept constant using a Berendsen
> > thermostat/barostat.  Temperature and pressure couplings of 1 and
> > 10 ps were used, respectively. Compressibility was fixed at 1 × 10-5
> > bar-1, and anisotropic pressure coupling was used in the MD simulations.
> > tcoupl                   = Berendsen
> > tau-t                      = 1.0
> > tc-grps                  = XXX
> > ref-t                       = 350
> >
> > pcoupl                  = Berendsen
> > pcoupltype         = anisotropic
> > compressibility  = 0.0001 0.0001 0.0001 0.0001 0.0001 0.0001
> > tau-p                     = 10
> > ref-p                      = 1.01325 1.01325 1.01325 1.01325 1.01325
> > 1.01325
> >
> > For all the calculations, a reaction-field correction was applied with
> > continuum dielectric εRF to correct for long-range interactions due to
> > electrostatics cutoff-scheme  = Verlet
> > coulombtype     = Reaction-Field
> > epsilon-rf            = a.bcd (I used actual numbers here of course)
> >
> >
> > In addition I am using
> > gen-vel                 = yes
> > gen-temp            = 350
> > gen-seed             = -1
> > continuation      = no
> >
> >
> > Should I be adding anything with respect to vdw?
> >
> >
> > In the forcefield.itp file I’m setting the following
> > nbfunc                  = 1
> > comb-rule           = 2 (Sum paper states that they use Lorentz-Berthelot
> > combining rule for the LJ parameters)
> > gen-pairs             = yes
> > fudgeLJ                = 1
> > fudgeQQ             = 1
> >
> > In ffbonded.itp
> > bondtypes func = 1
> > angletypes func = 1
> > dihedraltypes func = 5
> >
> >
> >
> > Thanks in advance for your help.
> >
> > Robert
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.


More information about the gromacs.org_gmx-users mailing list