Implementation of pair specific parameters in CustomNonBondedForce

[Robin Singh]

Hello OpenMM Developers,

Is there any way to implement pair specific parameters in the CustomNonBondedForce that can overwrite the parameters generated by the combination rule?
Suppose I want to use an energy expression (defined in CustomNonBondedForce) for my system which is salt in water. I defined all the combination rules and parameters and I want to simulate the non-bonding part of my system according to that energy expression. But for cation-O or anion-O interaction, I want to use different parameters (from what has been defined before) which can overwrite the parameters generated by the combination rule for these specific interactions only.
Is there any way to do so using CustomNonBondedForce in OpenMM?

Regards,
Robin

[Peter Eastman]

Yes, that's exactly what it does if you use the <LennardJonesForce> tag in a force field. See the documentation on it. And here is the source code for where it creates the CustomNonbondedForce. It assigns the atom type as a per-particle parameter, and uses the types of the two atoms to look up parameters in a pair of lookup tables.

[Robin Singh]

Hello Peter,

I was working on your suggestions and I have some queries. First thing is overwriting the parameters using the <LennardJonesForce> tag works when I have to use only Lennard Jones potential to represent every non-boned interaction of my system which is not what I want to do. I want to move away from Lennard Jones potential for ion-ion interaction but want to keep Lennard Jones for non bonding interactions containing tip3p water.

What I did was I defined a custom energy expression using <CustomNonBondedForce> tag and defined all the combination rules and gave all the parameters for the ions of the salt and used <LennardJonesForce> tag separately to represent water as shown below:

Note: I have used Lennard Jones potential in the <CustomNonBondedForce> tag just for the testing purpose. In my final work I'm gonna use different energy expression to represent ion-ion Non-bonded interaction.

<CustomNonbondedForce energy="4*epsilon*((sigma/r)^12-(sigma/r)^6)+1.38e-07*((charge1*charge2)/r); epsilon=sqrt(epsilon1*epsilon2); sigma=((sigma1+sigma2)/2)" coulomb14scale="0.8333333333333334" bondCutoff="3">
<UseAttributeFromResidue name="charge"/>
<PerParticleParameter name="charge"/>
<PerParticleParameter name="sigma"/>
<PerParticleParameter name="epsilon"/>
<Atom epsilon="0.8103692536" sigma="0.3037964628858557" type="tip3p_standard-K+"/>
<Atom epsilon="0.14891274399999999" sigma="0.4477656957373345" type="tip3p_standard-Cl-"/>
</CustomNonbondedForce>
<LennardJonesForce lj14scale="1.0" useDispersionCorrection="True">
<Atom epsilon="0.635968" sigma="0.31507524065751241" type="tip3p-O"/>
<Atom epsilon="0" sigma="1" type="tip3p-H"/>
<NBFixPair epsilon="0.635968" sigma="0.31507524065751241" type1="tip3p-O" type2="tip3p-O"/>
</LennardJonesForce>

And I got the below mentioned error:

File "/usr/local/miniconda3/envs/openforcefield/lib/python3.10/site-packages/openmm/app/forcefield.py", line 961, in getAtomParameters
raise ValueError('%s: No parameters defined for atom type %s' % (self.forceName, t))
ValueError: CustomNonbondedForce: No parameters defined for atom type tip3p-O

Then I used <LennardJonesForce> tag inside the <CustomNonBondedForce> tag as shown below:

Note: I have used Lennard Jones potential in the <CustomNonBondedForce> tag just for the testing purpose. In my final work I'm gonna use different energy expression to represent ion-ion Non-bonded interaction.

<CustomNonbondedForce energy="4*epsilon*((sigma/r)^12-(sigma/r)^6)+1.38e-07*((charge1*charge2)/r); epsilon=sqrt(epsilon1*epsilon2); sigma=((sigma1+sigma2)/2)" coulomb14scale="0.8333333333333334" bondCutoff="3">
<UseAttributeFromResidue name="charge"/>
<PerParticleParameter name="charge"/>
<PerParticleParameter name="sigma"/>
<PerParticleParameter name="epsilon"/>
<Atom epsilon="0.8103692536" sigma="0.3037964628858557" type="tip3p_standard-K+"/>
<Atom epsilon="0.14891274399999999" sigma="0.4477656957373345" type="tip3p_standard-Cl-"/>
<Atom epsilon="0" sigma="0.31507524065751241" type="tip3p-O"/>
<Atom epsilon="0" sigma="1" type="tip3p-H"/>
<LennardJonesForce lj14scale="1.0" useDispersionCorrection="True">
<Atom type="tip3p-O" sigma="0.31507524065751241" epsilon="0.635968"/>
<Atom type="tip3p-H" sigma="1" epsilon="0"/>
<NBFixPair epsilon="0.635968" sigma="0.31507524065751241" type1="tip3p-O" type2="tip3p-O"/>
</LennardJonesForce>
</CustomNonbondedForce>

And got the below mentioned error:

File "/usr/local/miniconda3/envs/openforcefield/lib/python3.10/site-packages/openmm/openmm.py", line 4388, in minimize
return _openmm.LocalEnergyMinimizer_minimize(context, tolerance, maxIterations)
openmm.OpenMMException: Particle coordinate is NaN.

Clearly, the usage of separate <LennardJonesForce> tag for water is not working.

I thought of defining a custom energy expression inside <CustomNonBondedForce> (different from Lennard Jones, may contain also to defines interactions of water) and then use Discrete2DFunctions to define the pair-specific parameters which can overwrite the parameters generated from the combination rules for certain pairs but could not understand the implementation of Discrete2DFunctions. I have read the API documentation also but that was not of much use.

Is this is a feasible way to overwrite parameters in <CustomNonBondedForce> tag? If yes then please let me know about its implementation in OpenMM.
If no then kindly let me now if there is any other way which can allow me to implement pair specific parameters in any custom energy expression defined in <CustomNonBondedForce> and not just <LennardJonesForce>.

Regards,
Robin

[Peter Eastman]

CustomNonbondedForce (and LennardJonesForce, which is implemented with a CustomNonbondedForce) implements a nonbonded interaction of every particle in the system with every other. That's just what it does. You must provide parameters for every particle. If you don't, you'll get an exception. What you can do, though, is set epsilon to 0 for certain types or certain pairs. It will still compute the interaction for them, but the energy of the interaction will always be 0.

You can't put the XML tag for one force inside the tag for a different force. It will just be ignored.

The LennardJonesForce tag uses Discrete2DFunctions as you describe. You can find the code for it here. The relevant lines are

Code: Select all

        self.force = mm.CustomNonbondedForce('acoef(type1, type2)/r^12 - bcoef(type1, type2)/r^6;')
        self.force.addTabulatedFunction('acoef', mm.Discrete2DFunction(numLjTypes, numLjTypes, acoef))
        self.force.addTabulatedFunction('bcoef', mm.Discrete2DFunction(numLjTypes, numLjTypes, bcoef))
        self.force.addPerParticleParameter('type')

[Robin Singh]

Hello Peter,

Thank you so much for your answers and time. It is clear to me how to overwrite the parameters using <LennardJonesForce> tag but I have an another query regarding other functional forms which I am gonna ask in a new thread.