In molecular mechanics, how are van der Waals forces modelled?

Question

In terms of energy, how are van der Waals forces modelled (are there formulas/laws that govern these)?

Answer 1

First I will try to directly answer this question:

In terms of energy, how are van der Waals forces modelled (are there formulas that govern these)?

The most common way to model the potential energy between two ground state (S-state) atoms that are far apart, is by the London dispersion formula:

$$V(r) = -\frac{C_6}{r^6}$$

where $C_6$ depends on the dipole polarizabilities ($\alpha_1$ and $\alpha_2$) of the two atoms. One decent approximation is the Slater-Kirkwood formula:

$$C_6 \approx \frac{3}{2}\frac{I_1I_2}{I_1+I_2}\alpha_1\alpha_2$$ where $I_1$ and $I_2$ are the first ionization potentials of the atoms.

However the London dispersion potential is not the only one:

• The Casimir-Polder potential is used in the relativistic regime, and it's often closer to $C_7/r^7$
• The resonance dipole-dipole potential: $C_3/r^3$ used between S-state and P-state atoms
• If one particle is charged you can get: $C_4/r^4$ as in Eq. 2 of this paper of mine.

In molecular mechanics, how are van der Waals forces modelled?

Most often the $C_6/r^6$ formula is used, which is reasonable unless dealing with ions, or excited states, or extremely long-range (relativistic) situations. However this formula is for two particles that are very far apart, and we need a force to go the other way when particles are too close to each other, and that is where the 1924 Lennard-Jones potential enters (it has already been written by AloneProgrammer, but in a different way):

$$V(r) = \frac{C_{12}}{r^{12}}-\frac{C_6}{r^6}$$

While the $r^6$ has rigorous theoretical foundations, the $r^{12}$ does not, but in molecular mechanics calculations, this function might need to be evaluated billions of times, so it is convenient that once you've calculated temp=r^6 in your computer program, you can just do temp2=temp*temp to get $r^{12}$. This might sound crazy now, but the earliest computers were so slow that being able to reuse the calculation $r^6$ in order to take a short cut to calculate $r^{12}$, actually made a big difference, and the most high-performance codes, even today, still use this short-cut.

---

However, now we have to address the comment of Charlie Crown:

LJ and EXP6 are merely two of the most cliche ways. They are both very simple and easy to use, but neither is all that great. LJ is too steep at close distances and EXP6 is unphysical at small distances

This is exactly what I told you: $C_6/r^6$ is only valid when the atoms are very far apart, and $C_{12}/r^{12}$ has no physical foundation at all (it is simply convenient since $(r^6)^2=r^{12}$.

AloneProgrammer gave the Morse potential (from 1929) which is actually really good when the atoms are closer together:

$$V(r) = D_e\left(1 - e^{\beta(r-r_e)}\right)^2$$

where $r_e$ is the inter-atomic distance at equilibrium, $D_e$ is the "Depth at equilibrium" and $\beta$ controls the shape. While this is good at short-range, it is bad at long-range, because if you work out the asymptotic behaviour as $r\rightarrow \infty$ you will see that it decays exponentially, when in fact we know it should decay with an inverse-power (proportional to $1/r^6$), and exponentials behave very differently from inverse-power functions.

The solution is the Morse/long-range function or MLR, which was introduced by Bob LeRoy and myself in 2009

It looks exactly like the Morse potential when $r$ is close to $r_e$ (when the system is close to equilibrium). But if you calculate the form of the function as $\lim\limits_{r\rightarrow \infty}$, you literally get $V(r) \rightarrow -u(r)$ where $u(r)$ can be anything you want: $C_6/r^6$, $C_3/r^3$, $C_4/r^4$, etc.

Therefore the MLR potential is Morse-like near equilibrium, and LJ-like when far from equilibrium, which is exactly what Charlie Crown said was problematic if you use pure Morse or pure Lennard-Jones.

The MLR potential isn't used in mainstream molecular mechanics calculations, because evaluating the function would be slower than simply using the $(r^6)^2=r^{12}$ trick (which makes calculations very fast when using the LJ potential). The MLR potential is more accurate though, and solves the problem of the LJ being wrong at equilibrium and the Morse being wrong when far from equilibrium, so it can give more accurate results. Often there's so many approximations going on in molecular mechanics that it doesn't hurt to use the LJ potential which both other answers mentioned already. The MLR tends to be used for high-precision spectroscopy more than for molecular mechanics, but it's an option if one wanted more accurate results.

Answer 2

Additional to Alone Programmer's answer. The Lennard-Jones Potential (LJ 12-6) is the standard, but is not unique, in some cases the factor of 6 is changed to 8 to simulate better the hydrogen bonds. Also, there is the Buckingham potential, where the repulsion part ($r^{12}$ term) is modified to the exponential term. But the attractive long-range term ($r^{6}$) is the same.

LJ 12-6 fits very well for the potential energy surface of a noble gas dimer. The epsilon term is the maximum depth of the curve, and sigma is the radius where the potential energy is zero (short distance). When the LJ is used to simulate the interaction of different atomic species, there is not a rule to determine the sigma and epsilon terms... and there are geometrical and arithmetic averages, using the values for the same-species interaction for each atom.

Answer 3

You are looking for Lennard-Jones potential. Basically, the interatomic interaction is modeled by this formula:

$$U(r) = 4 \epsilon \Bigg [ \Big ( \frac{\sigma}{r} \Big )^{12} - \Big ( \frac{\sigma}{r} \Big )^{6} \Bigg ]$$

Particularly the term $r^{-6}$ in the above formula describes long-range attraction force based on van der Waals theory.

Update:

I'll elaborate a bit more about my answer here. As Lucas said there is no universal model for capturing the behavior of van der Waals forces and you could generalize the Lennard-Jones potential as:

$$U(r) = 4 \epsilon \Bigg [ \Big ( \frac{\sigma}{r} \Big )^{m} - \Big ( \frac{\sigma}{r} \Big )^{n} \Bigg ]$$

As you see, due to the fact that always molecular dynamics boxes are finite in terms of length, no matter how big is your box and how far you molecules are $U(r) \neq 0$. This causes a problem when you put a periodic boundary condition and the image of an atom interacts with itself, which is obviously incorrect. You could modify this Lennard-Jones potential and define a cut off value to reduce $U(r)$ to zero for $r > r_{c}$ where $r_{c}$ is the cut off radius, but still induces some other trouble due to non-continuity of force (first derivative of the potential). Another common model to capture van der Waals forces is called soft-potential, defined as:

$$U(r) = D_{e} \Big( 1 - \exp{(-a(r-r_{e}))} \Big)^{2}$$

This potential comes from the solution of a quantum harmonic oscillator. $D_{e}$ is the height of the potential well, and $a$ controls its width.