Potentials based on bond stretching, bond bending, and long-ranged Coulomb interactions are widely used in molecular and organic systems. Chemists call these potentials ‘force fields.’ They cannot describe making and breaking chemical bonds, but by capturing molecular shapes, they describe the structural and dynamical properties of molecules well.

There are many commercial packages based on these force fields, for example, CHARMM35 and AMBER.36 They are primarily useful for simulating molecular liquids and solvation, but have seen little application in nuclear materials, on account of the long-range Coulomb forces, which are costly to eval­uate in large simulations.

1.10.10.2 Ionic Potentials

With no delocalized electrons, ionic materials should be suitable for modeling with pair potentials. The difficulty is that the Coulomb potential is long ranged. This can be tackled by Ewald or fast multipole meth­ods, but still scales badly with the number of atoms. The simplest model is the rigid ion potential, where charged (q) ions interact via long-range Coulomb forces and short-ranged pairwise repulsions V(r).

V (rij) +

For example, a common form of the pair potential in oxides consists of the combination of a (6-exp) Buckingham form and the Coulomb potential:

XA exp(-arj)-pi 4 + where a and b are parameters and r/ is the distance between atoms i and j.

As with other potential, various adjustments are needed in order to obtain reasonable forces at very short distance; see for example, recent reviews of UO2.37

For nuclear applications, the most commonly stud­ied material in the open literature is UO2, which is widely used as a reactor fuel. It adopts a simple fluorite structure with a large bandgap, which makes potential­fitting to get the correct crystal structure reasonably straightforward. Early work fitted the potentials to lattice parameter and compressibility, and later to elastic constants and the dispersion relation. The elas­tic constants are c11 = 395 GPa, c12 = 121 GPa, and c44 = 64 GPa.38 As previously described, the Cauchy relation generally applies to a pairwise potential,

C12 = C44, which is seldom true experimentally for oxides. However, the Cauchy relation is on the basis of the assumption that all atoms are strained equally, which is not the case for a crystal such as UO2 where some atoms do not lie at centers of symmetry. Thus, the violation of the Cauchy relation in UO2 can be fitted by attributing it to internal motions of the atoms away from their crystallographic positions. (The vio­lation of the Cauchy relation is similar in oxides with and without this effect, so it is debatable whether this is the correct physical effect.)

The earlier potentials were based on the Coulomb charge plus Buckingham described above; more recent parameterizations include a Morse potential. While this gives more degrees of freedom for fitting, having two exponential short-range repulsions with different exponents appears to be capturing the same physics twice. Comparison of the parameters39 shows that the prefactor for the U-UBuckingham repulsion varies by ten orders of magnitude when fitted. More­over, the original Catlow parameterization sets this term to zero. This difference tells us that the small U atoms seldom approach one another close enough for this force to be significant. Even the ionic charges vary between potentials by almost a factor of two, with more recent potentials taking lower values.

Despite the huge disparity in parameters, the size of cascades is similar and the recombination rate is high.

Polarizability is not incorporated in rigid ion poten­tials; they will always predict a high-frequency dielec­tric constant of 1, which is much smaller than typical experimental values. The main consequence of this for MD appears in the longitudinal optic phonon modes.

The solution is that ions themselves change in response to environment. A standard model for this is the shell model in which the valence electrons are represented by a negatively charged shell, connected to a positively charged nucleus by a spring. (Typically this represents both atomic nucleus and tightly bound electrons.) In a noncentrosymmetric environ­ment (e. g., finite temperature), the shell center lies away from the nuclear center, and the ion has a net dipole moment — it is polarized.

u=EV to-) +

у

In this case, ry may refer to the separation between nuclei i and j or the centers of the shells associated with i and i. In MD, the shells have extremely low mass, and are assumed to always relax to their equilibrium posi­tion: this is a manifestation of the Born-Oppenheimer approximation used in DFT calculation.

Shell model potentials,40 which capture the dipole polarizability of the oxygen molecules, were devel­oped by Grimes and coworkers, and have been through many extensions and reparameterizations since then. Again, there have been many successful parameterizations with wildly differing values for the parameters; even the sign of the charge on the U core and shell changes.41

A particular issue with ionic potentials is that of charge conservation. A defect involving a missing ion will lead to a finite charge. If the simulation is carried out in a supercell with periodic boundary conditions, this will introduce a formally infinite contribution to the energy. The simple way to deal with this is to ignore the long wavelength (k = 0) term in the Ewald sum. This, under the guise of a ‘neutralizing homo­geneous background charge’ is routinely done in first principles calculation. Alternately, a variable charge approach can be used42 in which the extra charge is added to adjacent atoms. The original approach then involved minimizing the total energy with respect to these additional charges, which is computationally demanding. A promising new development is to limit the range of the charge redistribution.43 While this screening approximation is difficult to justify fully in an insulator, it is very computationally efficient or a system involving dilute charged impurities, and appears to reproduce most known features of AlO.