# How to create a lookup table of ϵ and σ values for Lennard-Jones Potentials?

So in my last question How may I calculate the bond length between two atoms? I mentioned that I want to use Lennard-Jones Potential to simulate atom/molecule interaction for my sandbox-game-like application. However to use this I need to know ϵ and σ which are unique for each atom pair. I think the best advised idea (performance wise) is to generate a lookup table that lists ϵ and σ for every atom pair. It would be 13,924 records, which is okay. Now the question is what software could I use to fetch this data:

H+C | ϵ = x | σ = y
...

• You probably don't need every element and infact you cannot reliably even do all of the elements. If this is for a sandbox game, you might be overcomplicating things for yourself. Apr 28 at 16:00
• If there is a software, where I can input atoms and get these params it doesn't matter. It°s then up to me how much I need. Apr 28 at 16:22
• The problem is these parameters will be incorrect, so you should either spend the real work to get real values and a more complicated model, or just use values that get close to a real thing and ignore the real scientific model behind it. You will not get proper values for all elements on the periodic table that make any sense. Apr 28 at 16:56
• Yes ok. Do you have a proposal (s)? Apr 28 at 17:04
• John: I think at the top of one of your questions, if you give us a thorough explanation of what you're trying to do (explain this "sanbox-like-game") people would be able to help you much, much better. We have no idea what you're trying to do, so we're unable to advise your appropriately, which leads to you asking more and more questions and still not getting the answer you want. Apr 28 at 19:09

# Universal Lennard-Jones

The OpenKIM project is an Open Knowledgebase of Interatomic Models. It's integrated into a variety of simulation codes and allows you to pick out parameters for force field and interatomic potential simulations of materials.

The "Universal" Lennard-Jones model includes $$\sigma$$ and $$\epsilon$$ parameters for every element:

[It] was automatically fit using Lorentz-Berthelot mixing rules. It reproduces the dimer equilibrium separation (covalent radii) and the bond dissociation energies. It has not been fitted to other physical properties and its ability to model structures other than dimers is unknown.

While more specific atom types will certainly produce more accurate parameters (e.g., specialized kinds of carbon atoms), these may be a useful starting point.

They also have a tutorial on automated fitting a Lennard-Jones model if you wish to tweak the parameters.

• Do you mean this openkim.org/files/MO_959249795837_003/… ? Apr 29 at 15:24
• This combined with the answer by WaterMolecule is likely the best possible answer you will get and this should be accepted. Apr 29 at 18:53
• @JohnT - yes, that's the file, but I'm using the DOI as the correct citation link. Apr 29 at 19:30
• So H and H, means sigma etc is for the combination? And I can use combine rules to get individual? May 2 at 17:13
• @JohnT - right. If I have C-H, then I'd use the combination rules with the C and H parameters. May 2 at 18:40

## Lennard-Jones combining rules

Most force fields define a single σ and ε for each atom and calculate them for combinations using Lorentz-Berthelot combining rules. These rules specify a arithmetic mean for σ and a geometric mean for ε:

$$\sigma_{ij} = \frac{1}{2}\left( \sigma_i + \sigma_j \right)$$

$$\epsilon_{ij} = \sqrt{\epsilon_i \epsilon_j}$$

This reduces the number of parameters you need. Amber and CHARMM force fields use these rules, while OPLS uses the geometric mean for σ as well as for ε.

However, for force fields for organic molecules, there are usually more than one type (with different LJ σ and ε values) for each element. Below is an example that I got from "ffnonbonded.itp" of the Gromacs format parameter files for the CHARMM force field. σ and ε values are in nm and kJ/mol, respectively. This force field is mostly intended for biomolecules and drug-like molecules and will produce meaningless results when used outside this context.

  type  atnum sigma     epsilon  (description of the atom type)
AG  47    0.26326   19.07904 (silver metal)
AL  13    0.26059   16.81968 (aluminum metal)
ALG1  13    0.35636    2.71960 (aluminum tetrafluoride)
ALSI  13    0.39200    2.71960 (aluminum in aluminosilicates)
AU  79    0.26290   22.13336 (gold metal)
BAR  56    0.33676    0.62760 (Ba2+ ion)
BG201   5    0.45436    0.15899 (boron in boronic acid)
BRGA1  35    0.35101    2.00832 (organic bromine bromoethane-like)
BRGA2  35    0.36527    2.21752 (organic bromine 1,1-dibromoethane)
C   6    0.35636    0.46024 (peptide carbonyl carbon)
CA   6    0.35501    0.29288 (phenyl carbon)
CAI   6    0.35458    0.30543 (aromatic C next to CPT in tryptophan)
CAL  20    0.24357    0.50208 (Ca+2 ion)
CC   6    0.35636    0.29288 (carbonyl C, side chains of Asn, Asp, Gln, Glu amino acids)
CS   6    0.39200    0.46024 (thiolate carbon)


As you can see, there are different Lennard-Jones parameters based on the chemical environment of the atom. You will not be able to do anything useful with these parameters without putting the atoms in their expected context and implementing the other terms in the force field. The closest thing to a contextless force field that exists is ReaxFF and it needs to be quite complicated to handle the complexity of atomic interactions.

My recommendation is to you look at the force-fields data from software that implement Molecular Mechanics.

One example: at the Github site from the PYMATGEN project, you can get some library potentials files (in this case, used by GULP package).
For the DREIDING force-field, you have a file with data with bond, torsions and Lennard-Jones parameters: (note: there are other force-fields available like UFF, for example).

lennard epsilon 12 6 kcal x13 geometric all
0.0 12.5
hydrogen-bond kcal 12 10 4 taper dreiding
H_     N_3   N_3  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H_     N_2   N_3  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H_     N_2   N_2  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H_     N_R   N_3  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H_     N_R   N_2  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H_     N_R   N_R  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H_     N_1   N_3  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
...
...
...
H___b  I_    S_3  3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H___b  I_    F_   3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H___b  I_    Cl   3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H___b  I_    Br   3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0
H___b  I_    I_   3741298.1709492207 593660.5362167358 105.0 115.0 &
0.0 3.0 0.0 3.0 0.0 5.0

• Can you point out "Lennard-Jones parameters" ϵ and σ in your log? I can't read it. Apr 28 at 12:21
• I think it is needed a little more reading (software manuals, force-field papers) to understand how they are coded. Right now, I am running out of time. Sorry.
– Camps
Apr 28 at 12:28
• Looks like this is no good. Apr 28 at 16:09
• By the way, since you mentioned PyMatgen, our highest voted unanswered question is about PyMatgen and some other software: mattermodeling.stackexchange.com/q/1921/5. Do you have any idea how to answer that? Apr 29 at 2:56

Keep the following things in mind:

1. The Lennard-Jones potential will not give you chemical bonding
2. Morse parameters can only really be found for pairs of atoms that form diatomic molecules, as it does not describe polyatomic structures.

I am not exactly sure what you want, as you have received almost a dozen really good ideas in your previous questions, but none of them really fits your needs. If you have the background you can probably appreciate the work of this group, who probably do exactly what you want: take a larger amount of atoms, and just let them do their thing. However, such a simulation is not really possible in general, and it needs explicit potential energy surfaces for all the species involved.

You could also use the reactive force fields as many people (including me) have already suggested. However, if you can't use them, then I don't think what you are asking is possible. I don't know about any formula that could parametrize chemical bonding between arbitrary atoms. I just don't think such thing exists.

• While I tried to answer the question somewhat literally, this is the right answer. Reactive force fields are probably the only way to do what the OP wants and currently none of these can simultaneously cover large portions of the periodic table. May 4 at 16:22