It is clear that numerical methods essential when trying to understand the detailed atomic structure and dynamics. The most popular quantum mechanical approach is density functional theory [1] (DFT), which exploits an isomorphism between the ground state wave function and ground state electron density to approach the many-body Schrödinger equation variationally, often allowing efficient evaluation of complex structures and energy landscapes. However, even with significant simplifications to the underlying Hamiltonian such techniques are currently limited to less than 500 metallic atoms, meaning that they only be used to inform coarse grained models of atomic interaction. At present, short segments of screw dislocation cores can be simulated in DFT allowing accurate parametrisation of the Peierls barrier and dislocation formation energies, along with similar properties for point defect formation energies, as well as bulk properties such as the lattice parameter and bulk modulus [2]. It is the aim of any coarse grained model to quantitatively reproduce these features as accurately as possible. Born and Oppenheimer [3] were the first to realise that much insight into atomic rather than electronic motion can be gained if one works on a coarse grained timescale of femtoseconds, where the electronic co-ordinates are assumed to instantaneously assume the minimum energy configuration consistent with the given atomic positions. In this way one can argue that only the atomic positions need be retained along with some phenomenological interaction potential which attempts to recreate the energy landscape generated by the full quasistatic quantum mechanical system.In practice, molecular dynamics (MD) simulation aims to use an interatomic potential that best reproduces the features from DFT. The simplest choice is some isotropic pair potential φ(r ) such thatwhere the factor of 1/2 accounts for the exchange symmetry between i, j. Despite the appeal of this simple and efficient form, it can be shown that this can only produce