LAMMPS Implementation of Rapid Artificial Neural Network Derived Interatomic Potentials
11]D.Dickel , , M. Nitol , C. D. Barrett , Abstract
While machine learning approaches have been successfully used to rep-resent interatomic potentials, their speed has typically lagged behind con-ventional formalisms. This is often due to the complexity of the structuralfingerprints used to describe the local atomic environment and the largecutoff radii and neighbor lists used in the calculation of these fingerprints.Even recent machine learned methods are at least 10 times slower thantraditional formalisms. An implementation of a rapid artificial neuralnetwork (RANN) style potential in the LAMMPS molecular dynamicspackage is presented here which utilizes angular screening to reduce com-putational complexity without reducing accuracy. For the smallest neu-ral network architectures, this formalism rivals the modified embeddedatom method (MEAM) for speed and accuracy, while the networks ap-proximately one third as fast as MEAM were capable of reproducing thetraining database with chemical accuracy. The numerical accuracy of theLAMMPS implementation is assessed by verifying conservation of energyand agreement between calculated forces and pressures and the observedderivatives of the energy as well as by assessing the stability of the poten-tial in dynamic simulation. The potential style is tested using a force fieldfor magnesium and the computational efficiency for a variety of architec-tures is compared to a traditional potential models as well as alternativeANN formalisms. The predictive accuracy is found to rival that of slowermethods.
Keywords.
Machine Learning; Artificial Neural Networks; MolecularDynamics; a r X i v : . [ c ond - m a t . m t r l - s c i ] F e b AMMPS Implementation of Rapid ArtificialNeural Network Derived Interatomic Potentials [February 23, 2021
Artificial neural networks (ANNs) and other machine learning techniques aregaining prominence as a method for developing classical interatomic poten-tials for use in molecular statics (MS) and molecular dynamics (MD) simu-lation [1–10]. Because of their ability to accurately fit large and complex datasets and because they can be rapidly optimized, ANNs have shown the poten-tial to predict atomic energies and forces more accurately than most existingclassical formalisms [6, 7]. By using density functional theory (DFT), whichis relatively expensive computationally, to produce large training databases forthe ANN, potentials which reproduce this data to high precision at much fastercomputational speeds can be created. While only a limited number of large scaledynamic simulations have been performed using these potentials, static simula-tions for a variety of materials have been able to reproduce ab initio forces andenergies with chemical accuracy [2, 5, 11, 12].ANNs used to model interatomic potentials behave similarly to traditional for-malisms with an energy being defined for each individual atom based on thelocal atomic arrangement and the atomic species of it and its neighbors, map-ping the local environment directly to numerical values for the energy and forces.This encoding of the local environment is referred to as the structural finger-print. Thompson et al. developed the Spectral Neighbor Analysis Method(SNAP) [13] which uses machine-learning techniques on components of the lo-cal neighbor density. Behler and Parrinello [14, 15] defined a set of functionsbased on Gaussians which could be used to generate the fingerprint. This fin-gerprint formalism has since been used by a number of authors as the basisfor an ANN. The celebrated Gaussian Approximation Potentials (GAP) [16]have been used to develop potentials with chemical accuracy for a variety ofsystems including tungsten [17], lithium-carbon [18], iron [19] and water [20].Artrith and Urban used the Behler-Parrinello basis functions to create a poten-tial for titanium dioxide using an ANN with two hidden layers with 10 neuronseach which accurately described the various cystal phases of TiO [2]. Artrith,Urban, and Ceder have also recently proposed an efficient scheme for the gener-ation of machine-learned potentials for composite materials with many atomic2pecies [21]. Takahahsi et al [11] used a combination of the Gaussian formalismalong with insights from the modified embedded atom method (MEAM) [22] togenerate a fingerprint of several thousand values which was then used to fit apotential for titanium using linear ridge regression. Recently, a method usingthe embedding functions from MEAM has been used to generate potentials fortitanium and zinc using a single hidden layer [23, 24]. Singraber, Behler, andDellago [25] recently presented the n2p2 neural network package which uses theBehler symmetry functions in the construction of the fingerprint. This packagehas been used to develop a number of potentials including for Al-Cu [26] andMg [27].While the optimal structural fingerprint and network architecture for an ANN isstill unclear and indeed is probably strongly dependent on the material systembeing examined, there is a need to be able to perform large scale dynamic cal-culations with these potentials. Additionally, while machine learned potentialsgreatly outperform ab initio calculations in terms of efficiency, they still lag be-hind conventional formalisms such as the embedded atom method (EAM) [28]and MEAM, which have been the standard for analysis of metals since theirintroduction over 25 years ago. Additionally, MEAM takes advantage of angu-lar screening to limit the size of the neighbor list and incorporate the knownshielding effect to improve performance. To this end, we have created an im-plementation of a rapid artificial neural network (RANN) potential formalismfor use in LAMMPS [29]. This implementation accepts arbitrary Multi-layerPerceptron (MLP) architectures and structural fingerprints of arbitrary size.Sample potentials have been generated from the MEAM based structural fin-gerprints with angular screening, though novel functions or fingerprint typescan be added independently. Section 2 of this paper reviews the formalismfor ANN interatomic potentials and provides the specific structural fingerprintsimplemented including a discussion of the use of angular screening. Section 3describes the framework which has been implemented in LAMMPS. Section 4demonstrates a number of simple static and dynamic calculations using an ANNpotential. In particular, conservation of energy is shown along with accuratecalculation of the forces and pressure, and a comparison is done between thecomputational efficiency of this potential style as compared to other commonstyles for metal systems. It is demonstrated that, for the network architecturescapable of reproducing a magnesium training database to within 1meV/atom,performance approximately one third as MEAM can be achieved. Additionally,the performance of the potential as a function of neural network architectureand neighbor list size is considered. Most machine-learned interatomic potentials make use of MLPs to predict theenergy of individual atoms in a system given their environment. Details ofthe training and optimization of these networks in reference to ab initio datacan be found elsewhere [2, 23, 30], but the general structure of a network as an3nteratomic potential is as follows. The energy of a particular atom i , determinedby its environment, is the last of N layers of the neural network. The values forany particular layer, A n , after the first is determined by the previous layer andthe weight and bias matrices W n and B n : Z nl n = (cid:88) l n − W nl n l n − A n − l n − + B nl n (1) A nl n = g n ( Z nl n ) (2)Where l n is the number of neurons in layer n and g n ( x ) is a nonlinear activationfunction. The input layer, A is given by the structural fingerprint of the localatomic environment. The ouput layer, Z N , will always contain a single node,representing the energy, and so W N will always be a row vector and B N alwaysa single number.A number of different functions can be used to describe this structural finger-print. Behler and Parrinello [15] consider two types of fingerprint functions:radial symmetry functions which are a pairwise sum of Gaussians, and angu-lar terms over all triplets of atoms. For the RANN style, however, we use thefingerprint style introduced by Dickel, Francis, and Barrett, motivated by theModified Embedded Atom Method (MEAM) formalism [22, 31, 32], with theaddition of angular screening. This style seems particularly relevant for metalpotentials for the same physically motivated reasons that the MEAM formalismis effective. In this style, two different kind of input fingerprints are considered.First, simple pair interactions are considered and summed over all the neighborsof a given atom. For a given atom labeled i , we define a set of pair potentialsinteractions with the form F n = (cid:88) j (cid:54) = i ( r ij r e ) n e − α n rijre f c ( r c − r ij ∆ r ) S ij (3)Where j labels all the neighbors atoms of i within a cutoff radius r c , n is an in-teger, different for each member of the pairwise contributions to the fingerprint, r e is the equilibrium nearest neighbor distance, S ij is an angular screening termdescribed below, and α n are metaparameters, which can be tuned to better op-timize the potential. In principle, α n are related to the dimensionless α used inthe MEAM formalism which is related to the bulk modulus of the material [32]. f c ( x ) is a cutoff function which smoothly varies the weight of the interactionfrom 1 for atoms inside the cutoff radius to 0 once r ij > r c where ∆ r defines thewidth of the transition. The cutoff function used here is the same as employedin MEAM and is given by r c ( x ) = , x > − (1 − x ) ) , ≤ x ≤ , x < (4)The second kind of fingerprint function considers three body terms, with a formsimilar to the partial electron densities used in MEAM. There are two parame-4ers, m and k which can be varied to generate more inputs for the fingerprint: G m,k = (cid:88) j,k cos m θ jik e − β k rij + rikre f c ( r c − r ij ∆ r ) f c ( r c − r ik ∆ r ) S ij S ik (5)where θ jik is the angle between r ij and r ik , m is a non-negative interger and β k is a set of metaparameters which determine the length scale of the differentterms. For large numbers of neighbors, it can be more efficient to calculate thisas a sum over a single list of neighbors as follows. G ,k = (cid:88) j e − β k rijre f c ( r ij ) S ij (6) G ,k = (cid:88) α (cid:88) j r α i,j r i,j e − β k rijre f c ( r c − r ij ∆ r ) S ij (7) G ,k = (cid:88) α ,α (cid:88) j r α i,j r α i,j r i,j e − β k rijre f c ( r c − r ij ∆ r ) S ij (8) G ,k = (cid:88) α ,α ,α (cid:88) j r α i,j r α i,j r α i,j r i,j e − β k rijre f c ( r c − r ij ∆ r ) S ij (9) G m,k = (cid:88) α (cid:88) α ... (cid:88) α m (cid:88) j r α i,j r α i,j ...r α m i,j r mi,j e − β k rijre f c ( r c − r ij ∆ r )) S ij (10)where r α m ij is the x , y , or z component of r ij for α m = 1 , G ,k for example, there will be a calculation forxy as well as for yx. By taking advantage of this redundancy, we can reducethe number of terms that need to be calculated for G m,k from 3 n to ( n +1)( n +2)2 .The simplified expression for G m,k now reads: G m,k = (cid:88) α ≥ α ≥ ... ≥ α m n ! n x ! n y ! n z ! (cid:88) j r α i,j r α i,j ...r α m i,j r mi,j e − β k rijre f c ( r c − r ij ∆ r ) (11)where n x , n y , and n z are the total number of x-, y-, and z-components inthe set ( α , α , ..., α m ), respectively. Note that which of these forms (Eq. 5or Eq.11 will be more efficiently calculated will depend on the length of theneighbor list and the magnitude of m . Structural fingerprints of arbitrary sizecan be constructed by varying the number of different n , m , and k used givingmore or less information to the ANN. For the potentials tested here, we haveused consecutive values for n and m and totaling n tot and m tot respectively( n ∈ ( − , , , ..., n tot − m ∈ (0 , , ..., m tot − g N ( x ) = x (12) g n ( x ) = x
10 + 910 log( e x + 1) for n < N (13)The size of the weight and bias matrices will depend on the length of the fin-gerprint and the length of each hidden layer.As discussed above, one of the major computational barriers in the implemen-tation of ANN methods is the large neighbor list included in the calculationof the structural fingerprint. In order to obtain good convergence with DFTresults, cutoff radii used can be larger than 12 ˚A [27], which for the groundstate of Mg includes over 300 atoms. By comparison, recent MEAM potentialsfor Mg [33] have used a cutoff distance of 5.875 ˚A, which includes only 38 atomsin the ground state. Since the computation time scales at least linearly withthe number of neighbors, minimizing the length of the neighbor list withoutsacrificing performance should be a key goal for optimal efficiency. The radialscreening function f c can accomplish this, but leads to nonphysical results atlarge separations (see, for example [33]). Angular screening, whereby the ef-fective interaction between atoms is reduced or eliminated by the presence ofan atom located between them can effectively limit the neighbor list in a sim-ilar way without the unphysical results. Such a method has been utilized byMEAM [34], and the same screening method has been employed here. Briefly,the screening between two atoms is determined from the product of all thescreening interactions by other atoms in the neighborhood: S ij = (cid:89) k (cid:54) = i,j S ikj (14)where S ikj is calculated from a geometric construction considering the ellipseformed by the 3 atoms with r i,j one of the axes. The screening parameter, C ikj is then given by: C ikj = 1 + 2 r ij r ik + r ij r jk − r ij r ij − ( r ik − r jk ) (15)and then screening value is S ikj = f c (cid:18) C ikj − C min C max − C min (cid:19) (16)where f c is the same cutoff function used for the radial cutoff and C max and C min are metaparameters which can be tuned to determine which neighbors canbe excluded from calculations.The effect of including angular screening can be demonstrating by considering6igure 1: A sample value for a fingerprint in the presence (orange) and absence(blue)of angular screening as the volume of a cell is changed isotropically. Be-cause the neighbor list frequently changes without angular screening, the curveis seen to be less smooth, requiring more data to provide predictive informationto the ANN.the change of an individual fingerprint as the length scale is changed contin-uously. For this example, we consider the value of a single fingerprint for aMg atom in a perfect bulk hcp lattice as the lattice parameter is changed con-tinuously. In the absence of angular screening, the value of the fingerprintwill change rapidly at particular values of the lattice constant as new neigh-bors enter the radial screening distance. If angular screening is included withmetaparameters such that, for example, only 3rd nearest neighbors are ever in-cluded regardless of lattice parameter, the change in the value is considerablysmoother. Figure 1 demonstrates the difference between the two cases. Notonly is the unscreened value more expensive to calculate as it requires moreneighbors for most cases, it can also be seen that the predictive power of themodel is hampered as the relation between two states which are physically sim-ilar – ideal 3D lattices with different lattice constants, given markedly differenttrends and values as the number of neighbors changes. This can be overcomethrough the action of the neural network and precise cancellation among dif-ferent fingerprint terms, but the quality of fit and the extrapolative power ofthe screened functions appears more natural. Indeed unphysical deviations atrelatively moderate isotropic compressions, as seen for example in [27] can beavoided entirely through angular screening as the fingerprint values change con-tinuously in the screened case even at high compression.For the full network, forces are calculated by taking a gradient of the energyover the entire system. This includes the individual ANN contributions fromevery atom in the system, so the force on atom i is determined not only by theANN output of that atom, but of all of its neighbors as well.7 Numerical Implementation
The ANN potential file written for LAMMPS includes all of the details of theneural network. It begins with a header which lists the types and numbers offingerprints to be used. While we have only included the two varieties, the bondpower and radial power terms, F n and G m,k , in this work, the LAMMPS poten-tial style can be expanded in a straightforward way to include other descriptors.For each of the types used, the potential file then specifies the number used andany metaparameters. For F n this includes the range of values for n as well as themetaparameters α n , r c , C max , C min and r e . For G m,k , m tot and k tot are givenalong with the metaparameters β k , r e , r c , C min , and C max . Note that while weuse the same cutoff radius and angular screening parameters for both varieties,they can in principle be different as every fingerprint type has its own set ofmetaparameters. Next, the network architecture is specified by giving the totalnumber of layers and the length of each layer. Figure 2 shows a sample headerfor the potential file. Following the header is the weight and bias matrices foreach layer and finally the activation functions used for each layer. Using thepotential file, the ANN subroutine then calculates the feature list for every atomand uses these to calculate the energy and force for every atom in the usual way.The subroutine and potential file also have the capability of assigning multipleANNs to different atoms types to be used for multi-element simulations. In thisway, the framework can be used for arbitrary structural fingerprints as inputs toarbitrary MLP architectures for systems with an arbitrary number of elementtypes.The subroutine is included as a normal pair-style in LAMMPS called ‘rann’ andcan be called in the usual way:pair style rannpair coeff * * potential file.nn Element1 Element2 ... With the pair style implemented in LAMMPS, several potentials were fit for areference database for magnesium (Mg). The database, containing over 300,000atomic environments was constructed using the open-source density functionaltheory (DFT) software, Quantum Espresso [35]. The configurations used totrain the networks are shown in Table 1. These configurations included sim-ple triaxial strain perturbations of the low energy SC, FCC, BCC, and HCPphases of Mg along with larger supercells with atoms displaced from equilib-rium to mimic the effect of thermal excitation. Free surface data and vacancydata, with atoms also perturbed from equilibrium positions were also includedto improve stability. Effective temperatures up to 1000K were included in the8igure 2: Sample header for the neural network potential style. The headerspecifies the different element types as well as the style and number of fin-gerprints used for each type. Metaparameters needed for each style are alsoincluded as well as the overall network architecture.9hermally excited database. Training of the network was carried out by min-imizing the RMSE of the system energies compared to the DFT results usingthe Levenerg-Marquadt method [36, 37]. 10% of the training data from eachconfiguration type was left out of training and used for validation. The mostsignificant variation however, came in the size and construction of the finger-print used as the input layer. Along with the number and type of fingerprintsincluded, the most important factor in determining accuracy and efficiency wasthe number of atoms included in the neighbor list. This is a function of thecutoff radius r c and the angular screening parameter C min . Table 2 summa-rizes the different fingerprints considered and the RMSE for both the trainingand validation datasets for each potential. As can be seen, even the smallestarchitectures had RMSE of less than 3 meV/atom and the best were less than1 meV/atom. Table 2 also includes computational performance for a samplecalculation which will be described in detail below.Table 1: DFT simulations used to generate training and validation data for theAl potential. Training database attempts to consider all relevant low energystructures which may be encountered during MD simulations. Sample Description Atoms per simulations Number of Simulations Total atomic environmentsSC cubic cell w/ strains up to ±
15% 1 500 500FCC primitive cell w/ strains up to ±
15% 1 500 500BCC primitive cell w/ strains up to ±
15% 1 500 500HCP unit cell w/ strains up to ±
10% 2 2700 5400FCC 2x2x2 orthogonal supercell 32 500 16000BCC 3x3x3 orthogonal supercell 54 500 27000HCP 3x3x3 primitive supercell 54 2000 108000HCP 4x3x3 primitive supercell with vacancy 67 100 6700HCP 0001 free surface 54 500 27000HCP 1010 free surface 54 100 5400Totals 8100 307000
Tests of the validity of the potential and its implementation in LAMMPS wereall performed using the 13 member input fingerprint with a single hidden layercontaining 20 neurons. The metaparameters for this potential are shown in Ta-ble 3. C min was chosen exclude atoms beyond the third nearest neighbors in theequilibrium hcp Mg structure. First, to confirm agreement with first principlesresults beyond the RMSE in energy, the bulk properties of hcp magnesium weretested for this potential. Table 4 shows a comparison of these bulk propertiesamong the sample ANN, the first principles database, and experimental results.Agreement can be seen on the simple bulk properties targeted.As an initial test of the fidelity of the implementation, a number of simplestructures were generated and tested to confirm that the forces and energieswere correct and in agreement with one another. As a simple test, we displacea single atom from its equilibrium location in a minimized HCP Mg lattice.Figure 3 shows the force on the displaced atom as a function of displacement10able 2: Input Fingerprints considered in this work. k tot and m tot are the totalnumber of values uses for k and m in Eq. 11. For the pair-style input, n tot was always 5. All ANNs contained a single hidden layer of 20 neurons. r c and C m in are the radial cutoff and angular screening parameter that determine thenumber of neighbors considered in calculation of the fingerprints. k tot m tot r c C min n tot k tot m tot α β k r c AC min C max a (nm) 3.209 [38] 3.194 3.196 c (nm) 5.211 [38] 5.184 5.189 c/a (nm) 1.623 1.623 1.623 E c (eV) 1.51 [38] 1.45 1.46 C (GPa) 63.5 [39] 61.6 63.4 C (GPa) 25.9 [39] 23.8 24.3 C (GPa) 21.7 [39] 21.2 18.5 C (GPa) 66.5 [39] 64.3 62.9 C (GPa) 18.4 [39] 17.3 18.4As a second test, a periodic fcc cubic unit cell was generated and allowed torelax following a conjugate gradient minimization scheme. Symmetry is main-tained by the structures throughout minimization. The simulated cell is thenuniaxially compressed 10% in the x direction and slowly pulled in tension toa total positive strain of 10% and the changes in energy and stress tensor areexamined. The results are shown in Figure 4. We see, again, that the numericalderivative of the energy as a function of strain agrees exactly with the stresscalculated by the algorithm.Having demonstrated that the pair style correctly reproduces the energy, forces,and stresses in static configurations, we now confirm that dynamic behavior canalso be correctly produced. A dynamic simulation was carried out with a 256atom periodic hcp system (4x4x4 orthogonal unit cells). The system was ini-tialized from the equilibrium ground state with atoms given a random velocitydistribution producing an average temperature of 300K. The system was thenallowed to evolve at constant energy and volume. Figure 5 shows the evolutionof the kinetic, potential and total energy of the system. We see that total energyis conserved. While ANN potentials have the capability to reproduce experimental or firstprinciples results with greater accuracy than traditional formalisms, if the com-putational performance suffers too greatly, the benefits become negligible. Thespeed at which the ANN pair style will perform depends on the size and com-plexity of the structural fingerprint, the number of neighbors within the cutoffradius and the network architecture. 12igure 3: The force as determined by the numerical derivative of the systemenergy (solid red) and the calculated force (dashed black) as a single atom isdisplaced in the (1120) direction in a 4x4x4 hcp Mg lattice. (inset) Differencebetween derivative and direct force calculation.13igure 4: The xx component of the stress tensor (solid red) and the derivativeof the system energy (dashed black) as uniaxial tension/compression is appliedin the x direction. (inset) Difference between derivative and stress calculations.14
500 1000 1500 2000Time (fs)05101520 E n e r g y ( e V ) Figure 5: The potential, (red) kinetic, (black) and total (blue) energy relativeto the 0K minimum for a 4x4x4 hcp lattice of magnesium atoms with initialvelocity distribution producing a temperature of 300K. The system was allowedto evolve under NVE conditions. Conservation of energy can be observed forthis system. 15s a benchmark of the performance, we compare the performance of the MgANN potential to an magnesium Modified Embedded Atom potential [33] and arecent n2p2 style neural network potential [27]. To test the speed of the MEAMpotential and the various NN potentials, a sample system of hcp magnesiumcontaining 4000 atoms was minimized and each atom given a random velocitymatching a distribution with a temperature of 300K. The systems was thenallowed to evolve at constant energy and volume for 1000 timesteps. The num-ber of timesteps calculated per second is shown in Table 2. Calculations wereperformed on one core of a Intel Xeon E5-2690 2.0GHz processor. We see thatall of the potentials are slower than MEAM, with a strong dependence on thenumber of neighbors considered and on the cutoff radius r c , but the fastest po-tential, used in the validation study, runs approximately 1/3 as fast as MEAMand four times faster than competitive n2p2 potentials. As can be seen the firstcolumns of Table 2, there is a strong dependence on the cutoff radius even whenthe same number of atoms are included due to the angular screening parameter.This is due to the computational time required to calculate the angular screen-ing needed to remove atoms from the neighbor list. This is still seen to be moreefficient than considering all the atoms within the larger cutoff radius.For the size of networks considered here, the computation time is dominatedby the calculation of the structural fingerprint with the propagation throughthe network only a limiting factor for the smallest fingerprints and largest ar-chitectures. This suggests that, assuming the structural fingerprint containssufficient information to describe the relevant configurations, it is more efficientto increase the size of the network to increase the accuracy of the potential, us-ing this formalism, rather than expanding the fingerprint. This is particularlytrue for extending the length of the cutoff radius. While a large cutoff radius isoften desirable and is used in existing ANN potentials as it should improve thesmoothness of the potential as atoms move away from the target atom, the highcomputational cost of increasing the neighbor list suggests that ideally onlyfirst, second, and possibly third neighbors should normally be considered formetals. Angular screening combines both of these aspects, limiting the numberof neighbors in bulk settings while still allowing relatively long-range interac-tions for studies involving voids, free surfaces, or fracture.The RANN formalism presented here has the advantage of requiring similarcomputer time as MEAM to run simulations using the developed potentialsand, thanks to the additional of angular screening, faster than other existingANN formalisms, but the development time should be considerably smaller, asfitting to an arbitrary number of targets is automated with many free parame-ters available to tune to all available data. We have presented an implementation of the Rapid Artificial Neural Network(RANN) interatomic potential pair style for use in the LAMMPS molecular16ynamics code. The potential is suitable for static and dynamic calculations,conserving total energy and correctly calculating the pressure and forces of bulksolids. Such potential styles have been shown to be capable of reproducingtraining data from first principles calculations at a level which exceeds previ-ous semi-empirical formalisms [6, 7, 10]. The developed pair style is capableof accommodating network architectures of arbitrary dimension and activationfunction, and the calculation of the structural fingerprint requires only singlylooped summations over the neighbor list of a target atom, improving compu-tational efficiency. Angular screening has also been introduced both to improveefficiency by limiting the number of included neighbors and to improve predic-tive power by introducing the physically motivated phenomenon of shielding andsmoothing fingerprints as atoms move through the radial cutoff. Computationalefficiency is of particular importance for this formalism as the improvements inaccuracy over existing formalisms such as MEAM are greatly diminished if theruntime is closer to that of first principles calculations. As such, we demonstratethat for a network with 13 input fingerprints and a cutoff radius of 6.0 ˚ A forthe neighbor list, and angular screening restricting interactions to third nearestneighbors in bulk HCP the implementation performs at one third the speed ofMEAM. Larger networks, fingerprints, or neighbor lists will hamper this perfor-mance, with the most significant cost associated with increasing the neighborlist through the length of the radial cutoff. Ultimately, a flexible, scalable for-malism for ANNs, which can utilize the native parallel processing available inLAMMPS has been demonstrated with a formalism which produces high accu-racy potentials which operate on speeds comparable to those of MEAM. AcknowledgementsReferencesReferences [1] Hongxiang Zong, Ghanshyam Pilania, Xiangdong Ding, Graeme J Ackland,and Turab Lookman. Developing an interatomic potential for martensiticphase transformations in zirconium by machine learning. npj Computa-tional Materials , 4(1):48, 2018.[2] N. Artrith and A. Urban. An implementation of artificial neural-networkpotentials for atomistic materials simulations: Performance for tio2.
Com-putational Materials Science , 114:135–150, 2016.[3] Gabriele C Sosso, Giacomo Miceli, Sebastiano Caravati, J¨org Behler, andMarco Bernasconi. Neural network interatomic potential for the phasechange material gete.
Physical Review B , 85(17):174103, 2012.174] J¨org Behler. Representing potential energy surfaces by high-dimensionalneural network potentials.
Journal of Physics: Condensed Matter ,26(18):183001, 2014.[5] Daniele Dragoni, Thomas D Daff, G´abor Cs´anyi, and Nicola Marzari.Achieving dft accuracy with a machine-learning interatomic potential:Thermomechanics and defects in bcc ferromagnetic iron.
Physical ReviewMaterials , 2(1):013808, 2018.[6] Atsuto Seko, Akira Takahashi, and Isao Tanaka. First-principles inter-atomic potentials for ten elemental metals via compressed sensing.
PhysicalReview B , 92(5):054113, 2015.[7] Tran Doan Huan, Rohit Batra, James Chapman, Sridevi Krishnan, LihuaChen, and Rampi Ramprasad. A universal strategy for the creation of ma-chine learning-based atomistic force fields.
NPJ Computational Materials ,3(1):37, 2017.[8] Volker L Deringer and G´abor Cs´anyi. Machine learning based interatomicpotential for amorphous carbon.
Physical Review B , 95(9):094203, 2017.[9] Albert P Bart´ok, James Kermode, Noam Bernstein, and G´abor Cs´anyi.Machine learning a general-purpose interatomic potential for silicon.
Phys-ical Review X , 8(4):041048, 2018.[10] Ryo Kobayashi, Daniele Giofr´e, Till Junge, Michele Ceriotti, and William ACurtin. Neural network potential for al-mg-si alloys.
Physical Review Ma-terials , 1(5):053604, 2017.[11] A. Takahashi, A. Seko, and I. Tanaka. Conceptual and practical bases forthe high accuracy of machine learning interatomic potentials: Applicationto elemental titanium.
Physical Review Materials , 1:063801, 2017.[12] G. P. P. Pun and Y. Mishin. Optimized interatomic potential for siliconand its application to thermal stability of silicene.
Physical Review B ,95(22):224103, 2017.[13] Aidan P Thompson, Laura P Swiler, Christian R Trott, Stephen M Foiles,and Garritt J Tucker. Spectral neighbor analysis method for automatedgeneration of quantum-accurate interatomic potentials.
Journal of Com-putational Physics , 285:316–330, 2015.[14] J¨org Behler. Atom-centered symmetry functions for constructing high-dimensional neural network potentials.
The Journal of Chemical Physics ,134:074106, 2011.[15] J¨org Behler and Michele Parrinello. Generalized neural-network representa-tion of high-dimensional potential-energy surfaces.
Physical review letters ,98(14):146401, 2007. 1816] Albert P Bart´ok, Mike C Payne, Risi Kondor, and G´abor Cs´anyi. Gaussianapproximation potentials: The accuracy of quantum mechanics, withoutthe electrons.
Physical review letters , 104(13):136403, 2010.[17] Wojciech J Szlachta, Albert P Bart´ok, and G´abor Cs´anyi. Accuracy andtransferability of gaussian approximation potential models for tungsten.
Physical Review B , 90(10):104108, 2014.[18] So Fujikake, Volker L Deringer, Tae Hoon Lee, Marcin Krynski, Stephen RElliott, and G´abor Cs´anyi. Gaussian approximation potential modeling oflithium intercalation in carbon nanostructures.
The Journal of chemicalphysics , 148(24):241714, 2018.[19] D. Dragoni, T. D. Daff, G. Cs´anyi, and N. Marzari. Achieving dft accu-racy with a machine-learning interatomic potential: Thermomechanics anddefects in bcc ferromagnetic iron.
Physical Review Materials , 2(1):013808,1991.[20] Thuong T Nguyen, Eszter Sz´ekely, Giulio Imbalzano, J¨org Behler, G´aborCs´anyi, Michele Ceriotti, Andreas W G¨otz, and Francesco Paesani.Comparison of permutationally invariant polynomials, neural networks,and gaussian approximation potentials in representing water interac-tions through many-body expansions.
The Journal of chemical physics ,148(24):241725, 2018.[21] Nongnuch Artrith, Alexander Urban, and Gerbrand Ceder. Efficient andaccurate machine-learning interpolation of atomic energies in compositionswith many species.
Physical Review B , 96(1):014112, 2017.[22] M. I. Baskes. Modified embedded-atom potentials for cubic materials andimpurities.
Physical Review B , 46(5):2727, 1992.[23] D Dickel, DK Francis, and CD Barrett. Neural network aided developmentof a semi-empirical interatomic potential for titanium.
Computational Ma-terials Science , 171:109157, 2020.[24] Mashroor S Nitol, Doyl E Dickel, and Christopher D Barrett. Artificialneural network potential for pure zinc.
Computational Materials Science ,188:110207, 2021.[25] Andreas Singraber, J¨org Behler, and Christoph Dellago. Library-basedlammps implementation of high-dimensional neural network potentials.
Journal of chemical theory and computation , 15(3):1827–1840, 2019.[26] Daniel Marchand, Abhinav Jain, Albert Glensk, and WA Curtin. Machinelearning for metallurgy i. a neural-network potential for al-cu.
PhysicalReview Materials , 4(10):103601, 2020.1927] Markus Stricker, Binglun Yin, Eleanor Mak, and WA Curtin. Machinelearning for metallurgy ii. a neural-network potential for magnesium.
Phys-ical Review Materials , 4(10):103602, 2020.[28] M. S. Daw and Baskes M. I. Embedded-atom method: Derivation andapplication to impurities, surfaces, and other defects in metals.
PhysicalReview B , 29:6443, 1984.[29] S. J. Plimpton. Fast parallel algorithms for short-range molecular dynam-ics.
Journal of Computational Physics , 117(1):1–19, 1995.[30] J¨org Behler. Perspective: Machine learning potentials for atomistic simu-lations.
The Journal of chemical physics , 145(17):170901, 2016.[31] M. I. Baskes. Determination of modified embedded atom method parame-ters for nickel.
Mater. Chem. Phys. , 50:152, 1997.[32] B.-J. Lee and M. I. Baskes. Second nearest-neighbor modified embeddedatom method potential.
Physical Review B , 62:8564, 2000.[33] Z Wu, MF Francis, and WA Curtin. Magnesium interatomic potential forsimulating plasticity and fracture phenomena.
Modelling and Simulationin Materials Science and Engineering , 23(1):015004, 2015.[34] Michael I Baskes. Determination of modified embedded atom method pa-rameters for nickel.
Materials Chemistry and Physics , 50(2):152–158, 1997.[35] Paolo Giannozzi, Stefano Baroni, Nicola Bonini, Matteo Calandra, RobertoCar, Carlo Cavazzoni, Davide Ceresoli, Guido L Chiarotti, Matteo Cococ-cioni, Ismaila Dabo, et al. Quantum espresso: a modular and open-sourcesoftware project for quantum simulations of materials.
Journal of physics:Condensed matter , 21(39):395502, 2009.[36] Kenneth Levenberg. A method for the solution of certain non-linear prob-lems in least squares.
Quarterly of Applied Mathematics , 2(2):164–168, jul1944.[37] Donald W. Marquardt. An algorithm for least-squares estimation of non-linear parameters.
Journal of the Society for Industrial and Applied Math-ematics , 11(2):431–441, jun 1963.[38] N Saunders, AP Miodownik, and AT Dinsdale. Metastable lattice stabilitiesfor the elements.