[gmx-users] setting vdw cutoff with specific force-filed?

[Justin A. Lemkul]

Yun Shi wrote:
> Hi Justin and Mark,
> 
> Thank you very much for the reply.
> 
> I was using table 7 (Normal van der Waals Parameters) to calculate 
> non-bonded vdw interactions that are not between third neighbors, such 
> as CH1 carbons between different chains in a biomolecular system. 
> Anything wrong here?
> 
> I understand that it is the force that dictates the MD evolution, and I 
> calculated in this case as F = 12 * 9.85^2 / 1.5^13 = 5.98 kJ/mol/nm for 
> the repulsion term. The force from different directions on a atom in a 
> homogeneous system would cancel each other to some extend, but what 
> about the energy arises from this interaction? Would this considerably 
> affect the calculation of, say, binding energy of a ligand to a receptor 
> from thermodynamic integration or pulling simulation?
> 

Check your units and the column headings of Table 7.  Plugging in 9.85 as the 
energy will give you a wildly inflated result.  The C12 parameters listed are 
actually square roots and listed as 10^-3.  I think you will find the resulting 
energies and forces are vastly smaller for a simple interaction between two atoms.

> Besides, the GROMOS 53a6 paper used triple range scheme for calculations 
> of nonbonded interactions, and I guess it was rlist = 0.8 nm while rvdw 
> = rcoulomb = 1.4 nm. So is this considered to be accurate enough in 
> calculating free enthalpies of solvation since we know the interactions 
> between 0.8 and 1.4 nm were calculated every 5 steps?
> 

There is no need to update the neighbor list every single step.  Typically, 
water is the fastest-diffusing molecule in the system, but it will generally not 
have a dramatic displacement on the scale of 10 fs or so.

> The paper also used reaction-field instead PME to account for long-range 
> electrostatic interactions. I heard some people argue that PME would be 
> more accurate and it seemed to be utilized more often even in gromacs 
> tutorials. So does this mean certain accuracy could be achieved by using 
> triple range scheme and reaction-field together because the errors they 
> incur respectively somehow cancel out each other?
> 

PME is substantially more accurate.  Using it also requires rlist=rcoulomb, so 
the exact details of the Gromos96 derivation may be somewhat outdated.  Typical 
settings for Gromos96 would be something like:

rlist = 0.9
rcoulomb = 0.9
rvdw = 1.4
nstlist = 5
coulombtype = PME

Note that the value of rcoulomb and rlist can vary a bit as a consequence of PME.

-Justin

> Thanks a lot,
> 
> Yun Shi
> 
> 
> On 04/08/11, "Justin A. Lemkul"  <jalemkul at vt.edu 
> <mailto:jalemkul at vt.edu>> wrote:
>  >
>  >
>  > Yun Shi wrote:
>  > >Hi all,
>  > >
>  > >I am working with GROMOS 53a6 ff in GROMACS 4.5, and I assume a 
> Lennard-Jones interaction function was used for short-range vdw 
> interactions.
>  > >
>  > > From the reference paper /A Biomolecular Force Field Based on the 
> Free Enthalpy of Hydration and Solvation: The GROMOS Force-Field 
> Parameter Sets 53A5 and 53A6/, I found that for example,
>  > >
>  > >when rvdw = 1.5nm, the repulsion term of the interaction between two 
> CH1 type atoms (C12ij = 9.85^2) can be calculated as 9.850*9.850 / 
> (1.5^12) = 0.747786 kJ/mol. So I wonder if this value is considered to 
> be small enough to be ignored.
>  >
> 
> You should pay attention to the column headings in table 7 so that you 
> can compute the contribution correctly. However, the magnitude of the 
> energy of any particular interaction is not really of any concern. The 
> evolution of the system depends on the *forces*, and it is likely that 
> the sum of the forces on any atom from all its repulsion interactions 
> from atoms that are (say) 1.4nm to 1.5nm away is very close to zero, 
> except in highly non-homogeneous spatial distributions of particles. In 
> any case, the sum of that contribution will be much smaller than the 
> other contributions.
> 
> Mark
> 
>  >
>  > >
>  > >In addition, it seems not until 5 nm does the dispersion term become 
> larger than the repulsion term in this case, so would turning on 
> Dispersion Correction between, say 1.5 to 5 nm introduce more errors 
> than turning it off?
>  > >
>  >
>  > You should use the cutoff described the authors of the force field, 
> in this case rvdw=1.4.  Unless you can demonstrate that by using a 
> different value you can achieve superior results, stick with the 
> specifics of parameterization.  I have never seen ill effects of setting 
> rvdw=1.4 and using dispersion correction with this force field.
>  >
>  > -Justin
>  >
>  > --
>  > ==============================
> ==========
>  >
>  > Justin A. Lemkul
>  > Ph.D. Candidate
>  > ICTAS Doctoral Scholar
>  > MILES-IGERT Trainee
>  > Department of Biochemistry
>  > Virginia Tech
>  > Blacksburg, VA
>  > jalemkul[at]vt.edu <http://vt.edu/> | (540) 231-9080 
> <tel:%28540%29%20231-9080>
>  > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
>  >
>  > ========================================
>  > --
>  > gmx-users mailing list    gmx-users at gromacs.org 
> <mailto:gmx-users at gromacs.org>
>  > http://lists.gromacs.org/mailman/listinfo/gmx-users
>  > Please search the archive at 
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>  > Please don't post (un)subscribe requests to the list. Use the www 
> interface or send it to gmx-users-request at gromacs.org 
> <mailto:gmx-users-request at gromacs.org>.
>  > Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>  >
> 

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

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

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