Non-local Thermoelasticity Based on Equilibrium Statistical Thermodynamics

The purpose of this work is the formulation of energetic constitutive relations for thermoelasticity of non-simple materials based on atomistic considerations and equilibrium statistical thermodynamics (EST). In particular, both (unrestricted) canonical, and (re-stricted) quasi-harmonic, formulations are considered. In the canonical case, (spatial) non-locality results from relaxation of the assumption that atoms subject to continuum deformation change position uniformly and afﬁnely. In the quasi-harmonic case, the analogous assumption on mean atomic position (i.e., Cauchy-Born) is relaxed. Two types of spatial non-locality, i.e., strong and weak, are considered. In the former case, atomic position (or mean position) is a functional of the deformation gradient F , while in the latter, this functional is approximated by a function of F and its higher-order gradients ∇ 1 F ,..., ∇ n F . On this basis, canonical and quasi-harmonic non-local model relations are obtained for the thermoelastic free energy, entropy, internal energy, and stress. In addition, such relations are formulated for thermoelastic material properties (e.g., stiffness). In the second part of the work, basic relations from the continuum thermodynamics of (non-polar) simple materials are generalized to higher-order deformation gradient (i.e., weakly non-local) continua and applied to energetic thermoelasticity. The corresponding formulation is based in particular on (i) Euclidean frame-indifference of the energy balance and (ii) the dissipation principle. As in the standard case, necessary for (i) is linear momen-B tum balance and the symmetry of the (generalized) Kirchhoff stress (i.e., angular momentum balance). In the context of (ii), the free energy density determines in particular the ﬁrst Piola-Kirchhoff stress P , the higher-order stress measures P 1 ,..., P n conjugate to ∇ 1 F ,..., ∇ n F , as well as the generalized Kirchhoff stress. Modeling the phenomenological free energy on the corresponding (weakly non-local) canonical free energy yields EST-based constitutive forms for the entropy, all stress measures, and thermoelastic material properties. Alterna-tively, one can model the former energy as an approximation to the latter. An example of this for the second-order ( n = 2) case is discussed both theoretically and computationally in the last part of the work.


Introduction
The formulation of continuum constitutive relations based on discrete ab initio and/or atomistic considerations and methods rather than on continuum phenomenology alone has been pursued for many classes of materials in the literature. In the case of solids crystals, for example, Wallace [30, formulated zero-temperature hyperelasticity for simple materials based on an interatomic potential (see also, e.g., [23, §8.1.2]). This treatment was extended to finite temperature and energetic thermoelasticity (i.e., without heat conduction and viscosity) for simple materials by Wallace [30, and more recently [31] in the context of (quantum) phonon thermodynamics (see also, e.g., [23, §8.1.3]). Among other things, the establishment of (quantum) density functional theory for the quantitative determination of material properties (for a detailed review, see, e.g., [23,), as well as continuum modeling of nanoscopic systems and processes, has maintained interest in the formulation of continuum constitutive models based on ab initio and/or atomistic considerations. For example, a quantitative elastic stored energy model for simple materials was recently formulated in [24] based on atomistic considerations and quantum density functional theory (DFT) in terms of material-symmetry-adapted strain tensor components and determined the corresponding elastic stiffnesses.
One purpose of the current work is the generalization of such formulations for energetic thermoelasticity to non-simple materials in the context of equilibrium statistical thermodynamics (EST). As detailed in what follows, this results in general in strongly non-local (SNL) constitutive relations, i.e., relations which are functionals of the deformation gradient F . Also treated in this work is the weakly non-local (WNL) approximation of these, i.e., functions of the deformation gradient F and its higher-order gradients ∇ 1 F , . . . , ∇ n F . Both non-local formulations represent broad non-local generalizations of existing treatments (e.g., [ To employ such constitutive relations in the continuum modeling of non-simple materials, a formulation of corresponding basic field and balance relations is required. To this end, a second purpose of the current work is the phenomenological formulation of such relations for higher-order deformation gradient (i.e., WNL) continua in the context of continuum thermodynamics (e.g., [29]). Since the pioneering works of Mindlin (e.g., [14,15]) or Toupin (e.g., [26]), a number of extensions and generalizations have been pursued. For example, direct generalization of the formulation of [14] to geometrically non-linear isothermal gradient hyperelasticity has been carried out in [10]. In contrast to the variational formulation common in these works, a direct formulation is pursued in the current work in the context of the Euclidean frame-indifference of the energy balance (e.g., [9]). As such, the current work represents a generalization of the second-order case in [20] to arbitrary order. Since the focus in this work is on energetic thermoelasticity and EST, additional kinetic / dissipative constitutive relations (e.g., in the second-order case: [20]) are not considered here.
The current work begins in Sect. 2 with the constitutive formulation of energetic thermoelasticity in SNL form based on the (unrestricted) canonical ensemble and corresponding ensemble averaging. The corresponding WNL canonical formulation is given in Sect. 3. SNL and WNL constitutive formulation based on the quasi-harmonic (QH) approximation to the canonical ensemble is carried out in Sect. 4. This is followed in Sect. 5 by the phenomenological formulation of balance and field relations for higher-order deformation gradient continua based on continuum thermodynamics. This is then applied to the case of energetic thermoelasticity in the context of the dissipation principle, resulting i corresponding energetic thermoelastic constitutive relations, e.g., for stress. More detailed relations for these and related material properties (e.g., elastic stiffness) are obtained in Sect. 6 with the help of the free energy models from EST. Lastly, as a computational example, the EST-based WNL free energy is compared in Sect. 7 to second-order (i.e., n = 2) with its approximation via higher-order deformation gradient thermoelasticity in strain-gradient form. The work ends with a summary and discussion in Sect. 8. For completeness, reduced forms of the interatomic potential in the EST-based formulations are summarized in Appendix A in the context of material frame-indifference (e.g., [27,29]). Corresponding material frameindifference, reduced forms for the free energy density in higher-order deformation gradient thermoelasticity are discussed in Appendix B. Finally, the boundary-value problem for higher-order deformation gradient thermoelasticity is briefly summarized in Appendix C in variational form.
In this work, Euclidean vectors are represented by lower-case bold italic characters a, . . . , z, second-order Euclidean tensors by upper-case bold italic characters A, . . . , Z, and calligraphic characters A, . . . , Z for Euclidean tensors of arbitrary order. The notation A · B is used for the scalar product of arbitrary tensors. Given this product on vectors, (a ⊗ b)c := (b · c)a defines the tensor product a ⊗ b of a and b, and A T a · b := a · Ab the transpose A T of A. Let sym A := 1 2 (A + A T ) represent the symmetric part, and skw A := 1 2 (A − A T ) the skew-symmetric part, of A in what follows. Unless otherwise stated, upper-case subscripted slanted sans-serif characters A i , B i , . . . represent tensorial and non-tensorial quantities of order i + 2 for i 0 in this work. In particular, A 0 is then second-order. For i 1, any A i satisfying (A i b)a = (A i a)b is referred to as symmetric in what follows. Further concepts and notation will be introduced as needed along the way.

Strongly Non-local Canonical Formulation
As stated above, the current formulation is restricted to the simplest case of unary solids, primitive unit cells, and purely bulk relations (i.e., periodic system). In the corresponding canonical ensemble for a system of N of mass points at temperature θ , let r a represent the position of mass point a (a = 1, . . . , N) and p a = m aṙ a its momentum. As usual, the system Hamiltonian H (r, p) = K(p) + U(r) consists of kinetic K and potential U parts, with r := (r 1 , . . . , r N ) and p := (p 1 , . . . , p N ). Given these, the partition function and free Fig. 1 χ maps the position x r + s r of any point in a neighborhood of the reference position x r of a material element to x c + s c = χ(x r + s r ) in a neighborhood of its current position x c = χ(x r ). In the context of (3), this is the case in particular for r ra = x r + s ra and r a = x c + s a = χ(x r + s ra ) respectively, are determined, as well as ensemble averaging with respect to the canonical distribution function w. The short-hand notation dv(r) := dv(r 1 ) · · · dv(r N ) and dv(p) := dv(p 1 ) · · · dv(p N ) is employed here and in what follows, with dv(x) the volume element induced by dx.
Central to the current canonical formulation is the finite, non-affine generalization 1 of the "standard" infinitesimal, affine transformation r ra := F −1 r a (in the current notation) of phase-space position "coordinates" (e.g., [30, §7], or [23, Eq. (8.45)]). Here, χ is the continuum deformation field, F = ∇χ the deformation gradient, and x r the reference location of an arbitrary continuum material element with current location x c = χ (x r ) (Fig. 1). In the context of (3) and Fig. 1, the formulation to follow is with respect to a fixed, but otherwise arbitrary, continuum material element located at x r in the reference, and at x c in the current, configuration of the material in question. Analogous to its infinitesimal, affine counterpart (e.g., [23, §8.1.3]), (3) couples atomic and continuum kinematics. In particular, (3) induces the transformationṙ a = F (r ra )ṙ ra of velocities, and so the "canonical" transformation of K, with a := N a=1 . This results in the "discrete" functional of F . In addition, integration of dr a = F (r ra ) dr ra = F (x r + s ra ) ds ra (6) of (3) yields the functional of F , and so the corresponding one for the potential energy. Further, (6) 1 and dp a = F −T (r ra ) dp ra (at constant r ra ) from (4) 2 imply dv(r a ) = det F (r ra ) dv(r ra ), dv(p a ) = det F −T (r ra ) dv(p ra ), for the corresponding volume elements. Then dv(r a ) dv(p a ) = dv(r ra ) dv(p ra ) is invariant, Z from (1) and holds for the system free energy (1) 2 . System properties derived from (11) include as usual the entropy and internal energy respectively. In addition, the (functional) derivatives from (5) and (8) determine via ensemble averaging. Note that (D F r a ) T a · Z := a · (D F r a ) Z. In addition, (D F r a ) Z = sra 0 Z(x r + s r ) ds r from (7) for any Z(r r ). In particular, the choice Z(r r ) = AF (r r ) yields (D F r a ) AF = Ar a , and so via "push-forward" of D F K r and D F U r , respectively, from (13). In turn, then holds from (14).
In the context of material frame-indifference, the reduced form of (8) from (7). In particular, (17) results in the reduction of (13) 1 in terms of the bond force U ab := ∂ r ab U d and bond direction d ab := r ab /r ab , with and so (16) to the symmetric form Clearly, the symmetry of (D F U r )F T , and so that of (D F Ψ r )F T , here is a direct consequence of the material frame-indifference of U r [F ; x r , s r ] in (17). This is also true in the WNL formulation of energetic thermoelasticity, to which we now turn.

Weakly Non-local Canonical Formulation
This is based on the Taylor series expansions 3 where (Π i a)A := 1 i! (· · · (Aa)a · · · )a (i times) is a projection. From (22) follow analogous to r a [F ; x r , s ra ] and F (r ra ), respectively, with In turn, (23) result in analogous to (5) and (8), respectively, as well as analogous to (10), and (11), respectively. As for the SNL case (12), determine the entropy and internal energy, respectively, now with respect to (25) and (26).
On the other hand, in the context of these, "replace" D F Ψ r in the WNL case. Here, from (25) with via (23). The "push-forward" forms of (29) determine that via (30). Whereas g i ra has units of length, note that G i ra is dimensionless (e.g., like F and F ra in (23)). In contrast to ∂ F i Ψ r , which is of order (2 + i), note that (∂ F i Ψ r )F T i is of order 2 for all i. Related to this is the fact that the former has (SI) units of Jm i , while those J of the latter are independent i.
Like in the SNL case above, we have the reduced form of (25) In turn, (34) results in the reductions of (29) 1 and (31) 1 , respectively, with g i rab := g i ra − g i rb . Likewise, holds for the reduced form of (32). As shown in the continuum thermodynamic formulation of WNL energetic thermoelasticity in Sect. 5 below, (∂ F 0 Ψ r )F T 0 , . . . , (∂ Fn Ψ r )F T n determine for example the WNL form of the Kirchhoff stress. More on this later.
Combination of (40) and (41) then results in the QH approximations (11), and to Ψ r (θ, F 0 , . . . , F n ; x r , s r ) from (26) 2 , respectively. From (45) 2 and the first of these, we have ∂r ab ln λ r r ab ⊗r ab /r ab (47) as the QH approximation to the SNL relation (D F Ψ r )F T from (21). Similarly, (45) 2 and (46) 2 result in the QH approximation in particular the WNL QH form of the Kirchhoff stress in the continuum thermodynamic formulation of WNL energetic thermoelasticity.

Continuum Thermodynamic Formulation
Since the current treatment neglects heat conduction, the following is restricted to spatially uniform temperature θ . In this context, a direct generalization of the (referential or "Lagrangian") formulation of the first-order case in [20] to the current nth-order one is pursued here. For further simplicity, quasi-static conditions are assumed, and supplies are neglecting.
Given these conditions and assumptions, the total energy balance reduces tȯ in terms of the internal energy density ε r and the (mechanical) referential energy flux density h r . In the zeroth-order case of simple materials (e.g., [29, §3.6]), h r = P Tχ is determined by the material velocity 4χ and the first Piola-Kirchhoff stress P . Direct generalization of this to the current case yields 5 in terms of the higher-order hyperstresses P 1 , . . . , P n . As done here, we work with the notation P 0 := P analogous to F 0 := F in what follows for simplicity. By analogy with the zeroth-and first-order cases (e.g., [20]), the form of the energy balance resulting from (49) 4 For simplicity, the notation χ is used for both the deformation and motion fields in this work.
and (50) is Euclidean frame-indifferent (e.g., [29,Chap. 6]) iff linear and angular momentum balance respectively, hold, where 6 is the (generalized) Kirchhoff stress. As in the first-order case (e.g., [20]), the standard (i.e., zeroth-order) forms (51) of momentum balance apply in the nth-order case as well. Note that K from (52) reduces to its standard form P 0 F T 0 = P F T for n = 0. Together with energy balance (49) and linear momentum balance (51) 1 , the entropy balanceη r = π r for the current case of uniform temperature implies for the dissipation rate density in terms of the free energy density ψ r := ε r − θη r . Except perhaps for the energy flux form (50) relevant to higher-order deformation gradient continua without heat conduction, the above formulation is independent of any (further) constitutive assumptions. Restricting attention now to energetic thermoelasticity in this context, we have ψ r (θ, F 0 , . . . , F n ; x r ), η r = −∂ θ ψ r .
Substituting these into (53), we have In the current energetic thermoelastic context, exploitation of the dissipation principle (e.g., [29, §9.5 for supply-free case]) then yields the thermoelastic relations P n = ∂ Fn ψ r , and so π r = 0. Here, represents the ith-order variational derivative operator. On the basis of (54) and (56), we also have the reduced form for the Kirchhoff stress from (52).

Results Based on Free Energy Models from EST
Specific forms of (54) are determined by either the WNL canonical relations (26) relative to the reference volume V r . Since n i=0 g i rab = r ab and n i=0 G i ra = F ra in the context of (23), the right-hand side of (59) simplifies to n i=0 a<b Substituting this back into (58) yields the WNL canonical form of the Kirchhoff stress. This is the same form as that of (D F ψ r )F T resulting from (16) in SNL canonical case. Both of these can be compared with the WNL QH relation ∂r ab ln λ r r ab ⊗r ab /r ab (62) for K resulting from (48), and with the SNL-QH relation for (D F ψ qh r )F T resulting from (47). Since all of these relations are constitutively symmetric via material frame-indifference and (A.3), they satisfy angular momentum balance (51) 2 identically.
Finally, consider the second derivatives of ψ r based on the WNL canonical formulation, (26) 2 and (54), which determine equilibrium thermoelastic material properties (e.g., [29,Chap. 10]). These include the referential heat capacity the higher-order referential thermal expansion as well as higher-order referential elastic stiffness In particular, note that is obtained from (29) 1 and (36) 1 . The relations (64) and (65) represent direct generalizations of standard (i.e., n = 0) referential thermoelastic properties (e.g., [29, § §10.5-10.6]) to the current context. As done in the standard case, referential material properties such as ∂ θ ∂ F i ψ r and ∂ F j ∂ F i ψ r can also be expressed in spatial form. For example, the latter is transformed into the spatial form S ij c2 via the definition We then have analogous to S ij c2 . On this basis, one obtains via (30)

Example: Second-Order Gradient Elasticity
As an illustration and example of selected theoretical results in this work, consider lastly approximation of the EST-based WNL free energy (26) 2 by specific forms of the phenomeno-logical free energy density (54). For simplicity, attention is restricted to zero temperature and n 2 in this section.

Free Energy
In the context of material frame-indifference (see Appendix B), the phenomenological free energy density (54) takes the reduced form 7 via (B.2) for n = 2 at zero temperature with respect to the right Cauchy-Green deformation C 0 = F T F 0 = F T F = C as well as the first-order C 1 = F T F 1 and second-order C 2 = F T F 2 strain-gradient-like deformation measures from (B.3). Likewise, the EST-based WNL free energy (26) 2 simplifies to at zero-temperature for n = 2. Here, is determined by the interatomic potential (A.3), and follows from (35) with c ab (C 0 , C 1 , C 2 ) := C 0 s ra − C 0 s rb + (C 1 s ra )s ra − (C 1 s rb )s rb + (C 2 s ra )s ra s ra − C 2 s rb s rb s rb .
Assuming centrosymmetry (e.g., [1]), let (71) be given by the specific form based on second-order Taylor-series expansion of (73) about (C 0 , C 1 , C 2 ) = (I , 0, 0) with respect to a reference volume V r . Here, 0, I , 0, 0) and so on. In the standard n = 0 case, this is formally analogous for example to the approach of [24] based on material-symmetry-adapted strain tensor components. In particular, note that ∂ C 0 U i | 0 ≡ 0 for zero stress (e.g., Kirchhoff stress (85) below), while ∂ C 1 U i | 0 and ∂ C 1 ∂ C 2 U i | 0 vanish identically via centrosymmetry since C 1 is determined by an odd gradient of F . On the other hand, since C 2 is even in this sense, is generally non-zero. Besides this last relation, (73) also yields with U ab := ∂ r ab U d the atomic bond force, and U abcd := ∂ r ab ∂ r cd U d the atomic bond stiffness. Here, r ab (I , 0, 0) = s ra − s rb , from (74) with In particular, then, the isothermal elastic stiffness gradient" stiffnesses, are determined by (78).

Computational Comparisons
The following comparisons are based on the "plane-wave" deformation χ(x r ) = x r + a cos(kκ · x r + π/6) (81) of a fcc primitive unit cell with volume V r = a 3 /4 and central atom located at x r . Here, a is the lattice parameter, a is the displacement amplitude, κ ≡ 2π[001]/a is a wave vector in the first Brillouin zone, and k ∈ [0, 1]. In turn, (81) implies C 0 (x r ) = I + k 2 |a| 2 sin 2 (kκ · x r + π/6)κ ⊗ κ, As k increases, then, gradients of higher-order become increasingly important.
As an example comparison of V −1 r Ψ r (0, F 0 , F 1 , F 2 ) from (72) and ψ i (0, C 0 , C 1 , C 2 ) from (76), consider the results in Fig. 2 based on the EAM potential for Al of [32,33] for U d . The displacement amplitude a is assumed parallel to κ with magnitude a = 0.05 V 1/3 r (much smaller than the cutoff radius of the potential). As evident, the zero-and first-order straingradient energies (brown and green curves) begin to deviate both qualitatively and quantitatively from the second-order strain-gradient and EST-based cases (blue and red curves) for k above about 0.09. For k above about 0.34, even ψ i (0, C 0 , C 1 , C 2 ) begins to deviate quantitatively from V −1 r Ψ r (0, F 0 , F 1 , F 2 ).  ψ i (0, C 0 , 0, 0). Green curve (first strain gradient): ψ i (0, C 0 , C 1 , 0). Blue curve (second strain gradient): ψ i (0, C 0 , C 1 , C 2 ). See text for discussion (Color figure online)

Summary and Discussion
In this work, model relations for the energetic thermoelasticity of non-simple materials have been obtained with the help of equilibrium statistical thermodynamics (EST) as well as phenomenologically via continuum thermodynamics.
In the EST context, both the (unrestricted) canonical formulation, as well as the (restricted) quasi-harmonic (QH) approximation, have been considered. With respect to these, both strongly non-local (SNL) and weakly non-local (WNL) formulations of energetic thermoelasticity have been pursued. Principle results of the formulation include the SNL continuum form (11) for the EST canonical free energy Ψ r and its WNL approximation (26). In addition, the approximations (46) to these have been formulated in the context of the QH approximation to the unrestricted canonical formulation. As discussed in the text, the WNL form of either of these can be used to determine the phenomenological free energy density (54) for higher-order deformation gradient thermoelasticity.
A basic aspect of the continuum thermodynamic formulation of higher-order deformation gradient thermoelasticity in Sect. 5 is the generalization (50) of the (mechanical) energy flux h r in terms of the higher-order (hyper)stress measures P 1 , . . . , P n conjugate to F 1 , . . . , F n . Exploiting the Euclidean frame-indifference of the energy balance (e.g., [29,Chap. 6]), a direct consequence of (50) is the generalized Kirchhoff stress (52). In the context of the dissipation principle, (50) and (54) result in the generalized hyperelastic relations (56) for P 1 , . . . , P n and so the corresponding form (58) for the Kirchoff stress. As shown by the variational formulation of the corresponding boundary-value problem in Appendix C, the direct formulation of Sect. 5 is completely consistent with a variational one in the spirit of [14,15].
Central to the current treatment of atomistic-continuum coupling is the finite, non-affine relation r a = χ(r ra ) from (3) for the effect of continuum deformation on atomic position. For uniform (affine) local deformation χ(x r ) = F x r , note that (3) reduces to the "standard" infinitesimal, affine relation r a = χ(x r ) + F (x r ) s ra = F r ra (e.g., [23, §8.1.3]). From the point of view of Hamilton's equations, (3) represents a so-called "canonical" transformation of r 1 , . . . , r N (e.g., [8, 11, 18, 23, §8.1]). Indeed, the corresponding generator takes the form G(p, r r ) = − a p a · χ(r ra ) [cf. 23,Eq. (8.46)]. Then r a = −∂ p a G = χ (r ra ) corresponds to (3), and p ra = −∂ rra G = F T (r ra ) p a to (4) 2 . This is in contrast to the quasi-harmonic (QH) case, in which the mean atomic positions are assumed to deform with the continuum according to the generalized Cauchy-Born relation (39). As from (3) via in particular (7) and (23) in the canonical case, both SNL-and WNL-based QH formulations of energetic thermoelasticity ensue from (39) via (40).
The WNL formulation of atomistic-continuum kinematics in Sect. 3 and in particular (22) 1 can be compared for example with the recent treatment of [1]. They work with the polynomial-map-based approximatioñ to χ(x r ). Here, χ 0 is a constant vector, F 0 a constant second-order tensor, and F (i) 0 (i = 1, . . . , m) constant symmetric tensors of order i + 2. The distancer αβ := |r αβ | between atoms (or primitive unit cells) α and β withr αβ :=x α −x β is then determined bỹ in terms of the higher-order continuum deformation measures C (i) 0 := F T 0 F (i) 0 formally analogous to C i in (B.3). As discussed in [1], this can also be expressed with respect to E 0 and E (1) 0 , with C 0 = I + 2E 0 , and C (1) 0 = cyc E (1) 0 analogous to (B.9).
As investigated in previous work (e.g., [6,7,16,19,22]), other issues from material theory for gradient continua besides MFI include higher-order material symmetry restrictions. As discussed in detail for example in [22], since ∇ 1 F , . . . , ∇ n F do not transform tensorially under change of compatible reference configuration, such material symmetry restrictions invariably involve the consideration of higher-order jets. On the atomistic side, the material symmetry of U d from (A.3) is inherited by ψ i (C 0 , . . . , C n ). Yet another issue concerns a dependence of free energy on higher-order (e.g., anharmonic) strain terms such as E · (M 4 E)E in comparison to strain gradient terms like C 1 · G 4 C 1 . Although the corresponding moduli M 4 and G 4 are of the same order, they are different in character. Indeed, in contrast to G 4 , M 4 is lengthscale-independent. In any case, these and other issues represent work in progress to be reported on in the future.
Acknowledgements Open access funding provided by Max Planck Society. We thank the reviewers of the first version of this work for helpful comments which have lead to its improvement.
Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. For example, U(r 1 , . . . , r N ) = U t (r 12 , . . . , r N−1N ) satisfies translational MFI (recall r ab := r a − r b ). Using the Cauchy's representation theorem for simultaneous vector invariants (e.g., [27]), U t can be shown (e.g., [2,3]) to satisfy orthogonal MFI iff it can be expressed as a function of d 12 · d 23 , . . . , d N−2 N−1 · d N−1 N (recall d ab := r ab /r ab ). As usual, the simpler distance-based form

Appendix A: Reduced Forms of the Interatomic Potential
is (at least) sufficient to satisfy orthogonal MFI, and satisfies differential orthogonal MFI a r a × ∂ ra U = a<b r −1 ab (∂ r ab U d ) r ab × r ab = 0 identically. On the other hand, as discussed by [2,4,25], U d is not unique for potentials beyond 4-body interactions. For such cases, the "cluster" form of U d is appropriate. In this form, each U l n is a n-body potential, and M n (N/n) is the number of n-body interactions among N mass points. By requiring that each U l n vanish in the limit as the distance of any one of the interacting mass points to the others approaches infinity, this expansion can be formulated in a unique fashion.