Parameterising a Lennard-Jones interaction: Which atoms should I scan to build potential energy surface?

Question

I have a question related to parameterize Lennard-Jones interaction to get sigma $\sigma$ and epsilon $\epsilon$ to include in amber force field to do the Molecular Dynamics.

I searched a lot, they always said you should do scan between two atoms and then fit with Lennard-Jones Equation, but they don't explain which atom should I select or depend on what should I select the atoms. I am a bit confused.

In my case, I have a complex where I have two transition metals (Ru, Pt), where Ru have 6 bond with Nitrogen and Pt have two bonds with Chloro and two bonds with Nitrogen. I want to get the $\sigma$ and $\epsilon$ for both. But I don't understand which distance I should select to do the scan.

Can anyone help me to understand this point?

Here is the structure.

Answer 1

I searched a lot, they always said you should do scan between two atoms and then fit with Lennard-Jones Equation, but they don't explain which atom should I select or depend on what should I select the atoms. I am a bit confused.

Normally, pairs of atoms connected by chemical bonds are excluded from computation of non-bonded interactions because bonded energy terms replace non-bonded interactions. Also, it is common that in biomolecular force fields all pairs of connected atoms separated by up to 2 bonds (1-2 and 1-3 pairs) are excluded from non-bonded interactions. For example, in the image below, only the atoms 1 and 4 are considered for non-bonded interactions [1]:

But I don't understand which distance I should select to do the scan.

You create a new system with only the two atoms you want to calculate the interaction energy. You put them apart a fixed distance. Giving the x, y, z coordinates: you can select for example: $x_1=0$, $y_1=0$ and $z_1=0$ for the first atom, and $x_2=d$, $y_2=0$ and $z_2=0$ for the second one with $d$ being the distance between them. You start with $d=1$, do a single point energy, then increase $d$, do a single point energy, and so. At the end, you will have the data to plot the energy vs $d$ and then do the fitting.

Be aware that, if your atoms are not parametrized by the fore field you are interested in, you will need to obtain all the other force fields parameters as well (take a look here).

Answer 2

You linked to my answer about parameterizing force fields. For reference, this isn't my area of expertise, although a colleague's group does this frequently.

My understanding of your question is that you want to get the non-bonded interactions between the Ru and Pt atoms in the larger complex.

What you describe is going to be tricky, even with automated force field methods like I described on Chem.SE.

In your case, you don't have neutral metal ions - you have transition metals with a surrounding ligand sphere.

Ordinarily, you'd pick the relevant elements, maybe start from the initial distance between them, scan closer and scan farther until you see a full minima. (You want at least parts of the repulsive and attractive regions to fit the potential energy curve.)

For metals, it's more complicated. There are a few approaches.

I'd probably go with the "dummy atom" method - you surround the metal ions (Ru and Pt in your case) with partially charged dummy atoms, so the net charge is identical, but the ion is properly surrounded with the right coordination environment and has some charge delocalization.

A fairly readable article is "Force Field Independent Metal Parameters Using a Nonbonded Dummy Model".

So you'd add some dummy atoms to the Ru and Pt at the positions relevant in the complex (e.g., remove the carbon atoms and hydrogen around both, then change the nitrogen atoms around the Ru to dummy (Du) and the two nitrogen and two chlorine atoms around the Pt atom.

Then map the curve .. again starting near the distance in the complex, going shorter and longer until you've mapped out the LJ curve.

Since it looks like you're using Avogadro, you can use the "Align" tool, click on the Ru atom to make that the origin, then the Pt atom to make that along one of the axes (e.g, x-axis). Then you can select the Pt+dummy "molecule" and use the manipulate tool to translate it forward and backwards.

Merz also has a series of articles talking about the development of parameters for water models (TIP3P, etc.). These are probably less relevant to your needs, but discuss some of the challenges.

• "Rational Design of Particle Mesh Ewald Compatible Lennard-Jones Parameters for +2 Metal Cations in Explicit Solvent"
• "Parameterization of Highly Charged Metal Ions Using the 12-6-4 LJ-Type Nonbonded Model in Explicit Water"

I'm probably leaving off relevant papers.