Thermal expansion of FeNi Invar and zinc-blende CdTe from the view point of local structure

Thermal expansion of FeNi Invar and zinc-blende CdTe was investigated from the view point of local structure using the extended x-ray absorption fine structure (EXAFS) spectroscopic data and the path-integral effective classical potential (PIECP) Monte Carlo computational simulations. In this Review article, first the quantum statistical perturbation theory is intuitively described to see different character concerning thermal expansion between the quantum and classical theories. The diatomic Br2 molecule is employed as a simple example. Second, the PIECP theory is briefly described to note advantages and disadvantages of this simulation technique. Historical background is also discussed for the EXAFS investigation of thermal expansion based on the quantum statistical theories. The results of the FeNi Invar alloy are then summarized. The origin of zero thermal expansion in the FeNi alloy is ascribed to the so-called Invar effect that implies the variation of the electronic structure of Fe atoms depending on temperature. Zero thermal expansion at low temperature is however found to originate from the vibrational quantum effect. It is also noted that the interatomic distances of Fe-Fe, Fe-Ni, and Ni-Ni pairs are slightly but meaningfully different from each other, although the alloy exhibit a simple fcc crystal. Such a pairdependent difference is also true for thermal expansion and we will discuss thermal expansion from the local point of view, which is interestingly different from the lattice thermal expansion significantly. Finally, the results of the zinc blende (or diamond) structure are presented. Although the origin of negative thermal expansion in these tetrahedral crystals is known as a result of classical vibrational anomaly within the Newton dynamics theory, the quantum statistical simulation is found to be essential to reproduce the negative thermal expansion of CdTe. It is emphasized that the vibrational quantum effect and classical anharmonicity are of great importance for the understanding of low-temperature thermal expansion as well as the elastic constants. Page 2 Yokoyama et al . Microstructures 2021;1: I http://dx.doi.org/10.20517/microstructures.2021.001


INTRODUCTION
Since the discovery of anomalously large negative thermal expansion (NTE) in ZrW 2 O 8 over an extremely wide temperature range in 1996 [1] , NTE has revived as a hot topic in structural solid-state chemistry and physics [2] , and novel NTE materials such as Zn(CN) 2 have extensively been reported [3] . It is essentially important to understand detailed mechanisms of the material function of interest in order to further develop new functional materials. The origin of zero thermal expansion (ZTE) and NTE is known to be categorized into two mechanisms. The first one is a so-called Invar effect, in which the electronic structure of some atoms in the system is temperature dependent. The atomic radius becomes smaller with the temperature rise due to the transformation of the electronic structure, and the Fe 65 Ni 35 Invar alloy [4] is one of the best examples in this category. The second one is derived from vibrational anomaly. The materials with zinc blende or diamond structure mostly show NTE at low temperature due to the presence of vibrational anomaly [5] . Well known giant NTE materials such as ZrW 2 O 8 and Zn(CN) 2 also belong to this category.
It is well known that thermal expansion is caused by anharmonicity in the interatomic potential, especially the third order [6] . Although the lattice vibrational dynamics is often approximated within the assemble of harmonic oscillators, anharmonicity must be taken into consideration to account for thermal expansion. Moreover, in order to evaluate thermal expansion quantitatively especially at low temperature, it is essentially important to treat the anharmonic lattice system quantum mechanically. The quantum mechanical anharmonic vibration is usually treated as a perturbation from a harmonic oscillator, and it is well known that the presence of the zero-point vibration drastically makes the quantum effect important at low temperature. In order to describe the finite temperature effect in the lattice dynamics, the quantum statistical perturbation theory should be employed [7] . Although the quantum statistical perturbation theory can be applied to simple systems such as diatomic molecules, it is too complicated to construct analytical formalism and to perform numerical calculations of the solid NTE and ZTE materials. In real NTE and ZTE systems, computational simulation techniques should be much more suitable rather than analytical formulation methods. A computational simulation technique based on quantum mechanics, which treats quantum mechanical nuclear motion is known as the path integral theory [7,8] . The path-integral Monte-Carlo simulation method could be ideally the most appropriate to describe anharmonic lattice dynamics at finite temperature. It is however practically still too complicated to evaluate lattice dynamics of complicated solid NTE and ZTE systems. More practical approximation method should be helpful. In 1990s, the pathintegral effective classical potential method (PIECP) was developed exactly for the purpose of quantum mechanical anharmonic lattice dynamics and is practically quite useful [9,10] . In this theory, the path-integral quantum effect is included as a correction in the classical Monte-Carlo simulations, and we do not have to perform complicated quantum mechanical path integral calculations, being suitable to mostly describe thermal expansion of normal solids at low temperature.
The other concept other than the quantum effect in the present article is a local point of view in thermal expansion. Although thermal expansion itself is a macroscopic character, there can be seen several essential differences in thermal expansion between microscopic and macroscopic geometrical structures. In disordered alloy systems such as the Invar alloy, although the x-ray diffraction often exhibits simple crystal structure like fcc, different atomic pairs give different interatomic distances. In case of the Fe 65 Ni 35 Invar alloy, the crystal structure looks simple fcc, while the interatomic distances are not equivalent among Fe-Fe, Fe-Ni, and Ni-Ni pairs. In zinc blende materials as well, the first nearest neighbor (NN) interatomic distance does not correspond to the equilibrium distance calculated from the lattice constant, due to the presence of vibrational motion perpendicular to the bond direction at the equilibrium geometry. In order to reveal the local structure, especially to clarify the difference between the local structure and the equilibrium crystal structure, extended x-ray-absorption fine structure (EXAFS) spectroscopy is quite useful [11,12] . The EXAFS method provides detailed information on anharmonicity or asymmetric radial distribution quantitatively.
Let us here briefly review the EXAFS spectroscopy for anharmonic lattice dynamics from the historical point of view. In the conventional EXAFS formula, a narrow Gaussian function is assumed for the radial distribution function around x-ray absorbing atom. Although this approximation is often appropriate and is assumed even in the present standard EXAFS analysis, erroneous results such as too short interatomic distances would be obtained in the analysis of a large anharmonic system. This caution was firstly pointed out by Eisenberger and Brown [13] , and subsequently a general formulation using cumulant expansion was proposed by Bunker [14] . In the cumulant expansion, the radial distribution function is parametrized as cumulant moments up to the third or fourth order (up to the second-order expansion corresponds to the conventional Gaussian approximation), which will be given in detail in the next section. It is now recognized that appropriate consideration of the EXAFS cumulants (usually up to the third-or fourthorder cumulants) will advantageously provide additional important information on anharmonicity of the lattice vibration, rather than disadvantage of the EXAFS technique.
Theoretical formulation of the EXAFS cumulants has also been developed. The classic statistical dynamics theory easily provides the EXAFS cumulants assuming the Boltzmann distribution function, although the overall potential in the solid is rather complicated even in the case of summation of two-body interatomic potentials [15] . Frenkel and Rehr [16] showed the third-order EXAFS cumulant expression for a diatomic molecule using the quantum statistical first-order perturbation theory. More recently, Yokoyama [17,18] similarly derived the fourth-order EXAFS cumulant for a diatomic molecule using the quantum statistical second-order perturbation theory. Although the EXAFS cumulants obtained in this theory can be applied to rather simple polyatomic systems as triatomic [19] , octahedral [20] , and tetrahedral [21] molecules and a onedimensional infinite chain [22] , it is too complicated to adopt anharmonic vibration of crystalline solids, as mentioned above. Fujikawa et al. [23] firstly applied the PIECP method for the theoretical description of the EXAFS cumulants, where the one-dimensional chain system was again examined. Subsequently, Yokoyama [24,25] utilized the PIECP method to describe practical three-dimensional crystalline system. This extension allows one to apply the PIECP method to the ZTE and NTE materials rather straightforwardly.
In the present article, we will focus our attention on the quantum effect and the local structural point of view for thermal expansion, especially for the understanding of ZTE and NTE. We will recognize that these effects are quite important for thermal expansion phenomena of ZTE and NTE materials. In the subsequent section, after the definition of the EXAFS cumulants, the quantum statistical perturbation theory is intuitively recalled to understand thermal expansion in detail. Differences in thermal expansion between the quantum mechanical and classical theories is emphasized. The analytical formulas for the EXAFS cumulants will be a great help to recognize the relationship between the local and macroscopic average thermal expansion. The other anharmonicity of temperature dependent elasticity (softening with a temperature rise) is also discussed briefly. As an example, the experimental Br K-edge EXAFS result of diatomic Br 2 is presented and comparatively discussed in the quantum statistical second-order perturbation theory. In a subsequent section, the path-integral effective classical theory is described. It is here important for readers to understand advantages and disadvantages in this theoretical calculation, because the method contains many approximations to realize numerical simulations with low computational cost. In the next section, the local thermal expansion of the FeNi Invar alloy is discussed by using the experimental Fe and Ni K-edge EXAFS spectra, together with the PIECP and classical Monte Carlo (MC) computer simulations. The contributions of the quantum effect and the different atomic pairs to thermal expansion are elucidated. Subsequently, the local thermal expansion of CdTe is discussed by using the experimental Cd and Te K-edge EXAFS spectra, together with the PIECP and classical MC computer simulations. The origin of negative thermal expansion in these tetrahedral crystals is known as a result of classical vibrational anomaly that can be described within the Newton dynamics theory. In spite of that, the quantum statistical simulation is found to be essentially indispensable to reproduce the negative thermal expansion of CdTe appropriately. Finally, in the last section, the present review will be conclusively summarized and some perspectives on future directions will be discussed including the role of microstructures of the materials.

EXAFS cumulants
In this section, we will evaluate theoretical description of the structural parameters of the EXAFS cumulants. First, let us recall the EXAFS formalism, the structure parameters of which will be described by the following quantum statistical perturbation or path-integral theories. EXAFS is an oscillatory structure in the x-ray absorption spectra that appears above the threshold edge photon energy. The experimental EXAFS oscillation is interpreted by the quantum mechanical scattering theory [11,12] . In case of sufficiently large photoelectron energies, the EXAFS can appropriately be described within the single-scattering theory, except a collinear case where strong focusing forward scattering occurs. The EXAFS oscillation function χ(k) is defined by χ(k) ≡ [μ(k)μ 0 (k)]/μ 0 (k), where k is the photoelectron wave number given by  (E 0 the threshold photon energy corresponding to the absorption edge energy, E the photon energy, m e the electron mass, and h the Planck constant divided by 2p, and μ(k) and μ 0 (k) are the observed and hypothetical atomic background x-ray absorption coefficients as a function of k, respectively.
Within the single-scattering theory [11,12] , the EXAFS function χ(k) for a single shell is given as (1) where N is the coordination number of the equivalent surrounding atoms, R the corresponding interatomic distance, 2 0 S the intrinsic loss factor due to the many-electron effect on the core-hole creation, F(k) the backscattering amplitude of the scattering atoms, Φ(k) the total phase shift of the absorber and scattering atoms, and λ(k) the extrinsic loss factor due to the inelastic photoelectron scattering. C 2 , C 3 , and C 4 are the cumulant expansion coefficients (cumulants), defined as (2) where < > av denotes the thermal average [12,14] . Higher order cumulants than the fourth order are omitted here, and the interatomic distance R is regarded as the first-order cumulant, also given by R = <r> av . When the radial distribution function of the coordination shell around the x-ray-absorbing atom is sufficiently narrow in width and is expressed as a Gaussian function, all the cumulants higher than the second order are known to become zero (C n = 0 in case n ≥ 3). Practically, by evaluating F(k), Φ(k), λ(k), and 2 0 S theoretically and/or empirically, the structural parameters of N, R, C 2 , C 3 , and C 4 are obtained experimentally. Note here that R, C 2 , C 3 , and C 4 correspond to relative quantities between the absorber and scattering atoms and it is sometimes important to distinguish the difference between the local structure determined by EXAFS and the periodic crystal structure given by x-ray diffraction (XRD).

Quantum statistical perturbation theory
In this subsection, the EXAFS cumulants defined in eq. (2) is analytically derived, using the quantum statistical perturbation theory [7] for the description of thermal expansion and other thermodynamical properties of a simple two-body system such as a diatomic Br 2 molecule. The interatomic potential V(r) of a diatomic molecule is usually close to a harmonic oscillator and is given as a Taylor-expanded polynomial form: ( where the coefficients κ 0 , κ 3 , and κ 4 are the harmonic (second-), third-, and fourth-order force constants, respectively, and r 0 the distance at the potential minimum. Alternatively, we can quite often employ the Morse function as (4) where D e is the dissociation energy and α shows the width of the potential. It is known that the Schrödinger equation for the steady state can exactly be solved in the Morse potential function, and this is sometimes quite useful to discuss the physical properties of the molecules.
According to the quantum statistical perturbation theory, the perturbation factor f(β) is usually introduced, which is defined using a density function as . Using the expression, an expected value of any physical quantity M is calculated by (7) When the second and third anharmonic terms in eq. (3) are regarded as perturbed Hamiltonian and H' with the first harmonic term as nonperturbed Hamiltonian H 0 , the second-order quantum mechanical partition function Z qu is obtained by calculating eq. (6), and the result is given as (harmonic zero-point vibrational amplitude) and 0 / m ω κ = (harmonic eigenfrequency). The quantum mechanical interatomic distance R qu and third-order cumulants C 3,qu were are calculated within the first-order perturbation as [16] (9) Since the third-order force constant κ 3 is usually positive (the interatomic potential is steeper at a shorter distance side and is looser at a longer distance), the interatomic distance R qu exhibits thermal expansion, and correspondingly the third-order cumulant C 3,qu is also positive. The positive third-order force constant is thus the dominant origin of the thermal expansion. On the other hand, the second-and fourth-order cumulants C 2,qu and C 4,qu are given in the second-order perturbation by (10) and (11) respectively [17,18] . In C 2,qu , the first term corresponds to the nonperturbed term due to the harmonic oscillator, while in both C 2,qu and C 4,qu , the first-order perturbation terms originating from κ 3 are zero. Also note here that the first-order perturbation terms due to κ 4 [second and third terms in eq. (10) and first and second terms in eq. (11)] are negative because of a positive κ 4 , while the second-order perturbation terms due to 2 3 κ [fourth and fifth terms in eq. (10) and third and fourth terms in eq. (11)] are always positive. This implies that the sign of experimentally obtained C 4 allows us to recognize which term in eq. (11) contributes more dominantly to experimental C 4 , the first-order perturbation term due to κ 4 or the secondorder perturbation term due to 2 3 κ .
The classical interatomic distance R cl and the cumulants C 2,cl , C 3,cl , and C 4,cl , can immediately be obtained by calculating the convergent value at T→∞. This yields: The classical thermal expansion R cl -r 0 and harmonic C 2,cl are proportional to T, and the perturbed term of C 2,cl and C 3,cl are proportional to T 2 , while C 4,cl is to T 3 . Figure 1 shows the Morse potential of diatomic Br 2 , where D e = 2.113 (eV), α = 0.1661 (nm), and r e = 0.2284 (nm) in eq. (4) [26] . The light-blue lines correspond to the eigenstates exactly obtained by solving the Schrödinger equation. The black dot denotes the average interatomic distance at each eigenstate. It is easily understandable that the interatomic distance becomes longer as the eigenenergy increases, implying the occurrence of the thermal expansion in diatomic Br 2 . The force constants in eq. (3) are given as κ 0 = 2.459 × 10 2 (N/m), κ 3 =1.756 × 10 12 (N/m 2 ) and κ 4 =1.058 × 10 22 (N/m 3 ). Figure 2 gives the experimentally obtained EXAFS cumulants C 2 , C 3 , and C 4 of diatomic Br 2 by measuring temperature dependent Br K-edge EXAFS spectra of gaseous Br 2 , together with the calculated theoretically using eqs. (9)(10)(11)(12)(13) [17,18] . Agreements between the experiment and the calculation are fairly good. It is noted that C 2,qu shows nonzero value at T = 0 due to harmonic zero-point vibration. C 3,qu also exhibits nonzero value at T = 0, although the difference between the quantum and classical results is much smaller than in C 2 . Experimental C 4 is found to be positive, implying that the second-order perturbation term originating from the third-order force constant dominantly contributes to C 4 . Figure 3 shows the theoretically calculated thermal expansion and temperature dependent effective force constants. Thermal expansion exhibits a similar behaviour to C 2,qu , as expected in eqs. (9) and (12). At 2 3 κ low temperature, the deviation between the two lines is enlarged, and even at T = 0, the quantum result exhibits elongation with respect to r 0 due to the finite zero-point vibrational amplitude, while the classical result yields a linear function exactly proportional to T. The decrease in the effective force constant with a temperature rise implies the vibrational softening due to the third-and fourth-order anharmonic force constants. The classical result shows a linear decrease, while the quantum result exhibits a concave feature as in thermal expansion. From positive C 4 in Figure 2, the second-order perturbation due to the third-order force constant is concluded to be a dominant contribution to the temperature dependence of the effective force constant (vibrational softening).

General formulation
In the previous section, we found that the quantum statistical perturbation theory successfully gives analytical expressions for the EXAFS cumulants of a diatomic system. It becomes however more complicated rapidly as the number of atoms in a system increases. The quantum statistical first-order perturbation theory has been applied to triatomic [19] , octahedral (seven atoms) [20] , and tetrahedral (five atoms) [21] systems. The second-order perturbation theory however may not be practically applicable even in small molecular systems. To simulate anharmonic behaviours theoretically in many-atom systems such as solids, a computational simulation method should be more suitable. In this section, we will recall Feynman's path integral theory [7,8] and see a practical simulation method of the PIECP theory [9,10,23-25] for a possible application to ZTE and NTE solid materials.
According to the Feynman's path-integral theory [7,8] , the density matrix ρ(R) (R denotes the 3N-dimensional Cartesian coordinates of atomic positions) is given by using the real-space representation of the wave function as (14) where A[R(u)] is the Euclidean action such that (15) where M is the diagonal matrix that represents atomic masses, t X denotes the transpose matrix of X, and V[R(u)] is the total potential energy of the system. The PIECP method assumes a trial Euclidean action A 0 [R(u)] with some variational parameters. Since the present purpose is a description of vibrational properties of solids, the harmonic approximation such that (16) should be a suitable trial function. Here, R is the average path defined by  and F is the matrix that gives the second-order force constants under the Cartesian coordinate. The scalar quantity ( ) w R and the matrix F are the variational parameters to be optimized. Prior to the integral evaluation, the Cartesian coordinate R should be transformed into the normal coordinate Q using a linear transformation where U is the unitary matrix that is given as eigenvectors of the matrix M -1/2 FM -1/2 . This immediately yields (18) where ω 2 is the diagonal eigenfrequency matrix given as The path-integral evaluation can exactly be performed for the multidimensional harmonic oscillators, and the resultant density opt ( ) ρ R within the framework of the PIECP theory is given as (19) where α kμ denotes the pure quantum fluctuation (difference of fluctuation between the quantum and classical estimations) of the normal vibrational mode k and the phonon branch μ defined by (20) Here, although we implicitly assume the periodicity by using the wave number vector k and the phonon branch μ for convenience, the periodicity does not have to be imposed. The thermal average of a certain physical quantity M within the PIECP approximation is optimized as (21) where ( ) eff V R is the effective classical potential defined as (22) and << >> denotes the 3N-dimensional integral average concerning the pure quantum fluctuations.
For the optimization of ( ) w R and ω 2 , one can employ the Jensen-Feynman inequality [7] as (23) where F and F 0 are, respectively, the true and trial free energies. The eventual variational conditions can simply be given by (24)

Low coupling approximation
Although considerable approximations have already been employed and the resultant density opt ( ) ρ R has been reduced to rather a simple expression, eqs. (19) and (21) still look too complicated to employ in numerical calculations, because they contain the 3N-dimensional integration originating from the quantum fluctuations. For the one-dimensional system, no further approximation is necessary, while in the multidimensional system like solids (N~10 23 ) the quantum fluctuation integrals should be separated out. Here, we should come closer to the harmonic approximation, which is the so-called low coupling approximation. This approximation consequently assumes that the variational parameters of ( ) w R and ω 2 are independent of R .
For the evaluation of the effective potential V eff in eq. (22), let us assume the pairwise potential and the Bravais lattice. In such a case, a final expression for V eff is derived explicitly as (27) and (28) where u ij (R ij ) is the classical temperature-independent pairwise potential between atom pair I and j with the instantaneous distance R ij , 0 ij R is the equilibrium distance, and is well known as the EXAFS Debye-Waller factor C 2 when one replaces α kμ with the usual mean-square fluctuations of (h/2ω kμ )coth f kμ [27] and thus (2)L ij σ can be regarded as a contribution of the quantum fluctuation to C 2 .

Embedded-atom method
For the description of the thermal properties of metals, the pair-potential approximation is not sufficient since the metallic bond should inherently contain many-body interaction due to the presence of delocalized electrons. The embedded-atom method (EAM) idea was firstly developed [28][29][30][31] . The effective medium theory [32] and the glue model [33] belong to the same category, these being based on the density-functional theory. In the EAM, it is argued that the total adiabatic potential energy E of the system can be written as (31) where ρ h,I is the total electron density at atom I due to the host (the rest of the atoms in the system), F i is the embedding energy for placing atom I into the electron density, and Φ ij is the short-range pairwise corecore repulsion that is defined as a function of the interatomic distance R ij . Although the embedding energy F i should be described as a complicated density functional, a so-called local-density approximation was introduced in eq. (31) and F i was simply given as a normal function of the host density ρ h,I at the position of atom i. In case of isotropic metallic crystals, ρ h,I can be given as a sum of the spherically averaged atomic densities , , and angular dependence can also be included by modifying the summation process, which is useful for instance to describe bcc metals appropriately. In most EAM simulations, the numerical values of F i and Φ ij are determined empirically, instead of the density functional calculations. The EAM does not require the three-dimensional periodicity and has thus been applied to alloys, defects, and surfaces.
When one tries to apply the EAM scheme to the PIECP method, the dynamical matrix should at first be calculated within the harmonic approximation. For this purpose, the EAM energy E i of atom I is simply Taylor expanded around the equilibrium host electron density , h i ρ . The resultant formula is given as (32) where u ij is the pairwise potential between atoms I and j, which is given as (33) Here, we have reached a very important conclusion. Although the EAM scheme inherently treats manybody forces of metallic bonds, the pair potential scheme is guaranteed in the Taylor expansion. This comes from the assumption that the embedding energy F i is described only as a function of the sum of the spherically averaged atomic density of electrons, which depends only on the interatomic distance. No angular dependence such as a bending force needs to be taken into account within this approximation. This allows us to apply the EAM to the PIECP theory in a straightforward manner. Namely, eq. (33) can be employed by replacing the first term of the total adiabatic potential in eq. (27) with eq. (32).

Summary of PIECP
At the end of this section, let us summarize the advantages and disadvantages of the PIECP method when we apply the method to describe thermal expansion or other vibrational thermal properties of solids like anharmonic elastic constants. There can be found several advantages as mentioned above. The vibrational quantum effect can be taken into account, especially the zero-point vibration within the harmonic approximation. Anharmonicity can also be included because of the applicability of many kinds of empirical potential schemes. We can employ for instance the EAM to describe metallic systems and the bondorder potential for covalent systems, as exemplified below. Lennard-Jones and Morse potentials can also be utilized, of course. Moreover, since we do not have to perform complicated quantum mechanical pathintegral evaluation, the computer cost should be almost equivalent to the classical simulations.
We have however imposed many approximations to reach the final PIECP formulism, especially in solid materials. The most serious approximation may be a low-coupling approximation introduced to do without 3N-dimensional integral calculations due to the quantum fluctuation. The low-coupling approximation neglects the interaction between phonons with different modes or branches, discarding most quantum mechanical anharmonicity. Consequently, the anharmonicity after the low-coupling approximation is almost classical. Fortunately, quantum mechanical anharmonicity is usually quite small except for instance helium or hydrogen, and the low coupling approximation can be applied to most metallic systems without losing high reliability. It is noted that classical anharmonicity is guaranteed in the low coupling approximation. The interatomic distance can correctly be taken into consideration, as in the quasi- harmonic approximation of periodic systems. Superiority to the quasi-harmonic approximation is found in the applicability to the nonperiodic systems, although the dynamical matrix calculation should be performed with the assumption of some averaged periodicity to estimate quantum mechanical correction in the effective classical potentials. Moreover, the PIECP method is applicable to the double well potential system that cannot be investigated by the quasi-harmonic approximation. Although there can be found several limits in the PIECP method, this approach will be the best candidate in many cases to evaluate thermal expansion of solids.

INVAR ALLOY
Anomalously small thermal expansion over a wide temperature range in an iron-nickel alloy with a nickel concentration of around 35% was discovered by Guillaume [4] in 1897. Figure 4 shows temperature dependence of the interatomic distance and Young moduli of typical metals. The Fe 65 Ni 35 Invar alloy actually exhibits almost ZTE, while stainless steel (SUS304, Fe 72 Ni 9 Cr 18 Mn 1 ) gives much larger normal thermal expansion. Interestingly, in Figure 4B, the Invar alloy shows an increase in the Young modulus with the temperature rise, which is essentially different from other normal metals as SUS304, Ni, and Cu. In Figure 4, the behaviours of Ni Span C (Fe 50 Ni 42 Cr 6 Ti 4 ) are also given, which is sometimes called the Elinvar alloy exhibiting almost no temperature dependence of the elastic constant, accompanied by smaller thermal expansion. The Invar and Elinvar alloys have been utilized in various kinds of industrial products.
The ZTE behaviour of Fe 65 Ni 35 is well known as a result of the so-called Invar effect. Other physical properties such as elastic constants and magnetization show also anomalous behaviours [35][36][37][38] . A basic concept of the Invar effect is that there exist at least two types of electronic states in Fe, typically highspin (HS) and low-spin (LS) states. In this two-state model, the equilibrium potential energy is lower in the HS state than in the LS one, while the equilibrium atomic radius is larger in the former. This results in the compensation of thermal expansion due to increasing density of the LS state at higher temperature. A recent ab initio electronic structure calculation at 0 K, however, suggests much more complicated electronic  [34] .
configurations in a smaller volume region [39] . Computational simulations at finite temperatures have also been carried out for the understanding of magnetization and thermal expansion [40][41][42][43][44] .
In this section, let us see the experimental and computational results of the local thermal expansion and the anharmonic behaviour [45] . Temperature dependent Fe and Ni K-edge EXAFS spectra of the Fe 65 Ni 35 Invar alloy were recorded using synchrotron radiation (Photon Factory, Tsukuba, Japan) and the MC simulation results based on the PIECP and classical theories within the simple two-state (HS + LS) model. Since we discuss thermal expansion, the MC simulations were performed under a constant pressure condition (NPT condition). 500 atoms (5 3 fcc cubic unit cells) were employed and the distributions of Fe and Ni were chosen randomly. Eleven types of the superlattices were simulated and the results were averaged to provide consequent physical quantities. Figure 5 shows the first NN interatomic distances around Fe and Ni obtained by the experimental EXAFS analysis [45] and the lattice constant given by the X-ray diffraction [46] , together with those simulated by the PIECP and classical MC methods. Figure 5A gives the potential energies of Invar Fe 65 Ni 35 and fcc Fe at T = 0 as a function of the first NN distance. In the present atomic potentials, the fcc Fe system shows that the LS state is slightly more stable by 8.0 meV than the HS state with R HS = 0.2530 (nm) and R LS = 0.2492 (nm), while the Invar case exhibits a more stable HS state by 25.0 meV with R HS = 0.2530 (nm) and R LS = 0.2490 (nm). In Figure 5B-D, the agreement between the PIECP and experiments is fairly good: almost no thermal expansion around Fe and small but meaningful thermal expansion around Ni. On the contrary, the classical method is found to give fatal discrepancy at low temperature below ~100 K; the interatomic and lattice distances significantly increase with the temperature rise. These findings imply the essential importance of the quantum effect, which is recognized as a result of zero-point vibration.  [44] . To get further insights into local thermal expansion, the first NN interatomic distance of each component (Fe-Fe, Ni-Ni and Ni-Fe) pair is shown in Figure 6. In this plot, the PIECP MC results, by using only the HS Fe state, are also depicted to recognize hypothetical normal thermal expansion in this system. As expected, the Fe-Fe pair shows the largest discrepancy between the two-state (HS + LS) and the HS-only models. This is caused by the increasing population of the Fe LS state with the temperature rise, yielding compensation of thermal expansion with the one originating from anharmonic vibration. The most important finding in Figure 6 is that even Ni-Ni and Ni-Fe pairs exhibit significant suppression of thermal expansion. This is consistent with the above experimental finding that thermal expansion around Ni is noticeably suppressed compared to that of fcc Ni. Although Ni does not change its electronic configuration depending on temperature and may tend to expand because of anharmonic vibration, the Ni-Ni or Ni-Fe expansion is significantly suppressed due to the anomalously small expansion of the lattice. Interestingly, the suppression of thermal expansion seems more significant in the Ni-Fe pair than in the Ni-Ni one. This can be explained by the fact that the Fe atom surrounded by many Ni atoms tends to maintain the HS state. In the present EAM potentials, the energy difference between the HS and LS states is determined by the sum of the d electron density of surrounding atoms at the Fe site. Although the number of the d electrons of Ni is larger than that of Fe, the d electron density at the first NN distance is smaller in Ni than in Fe because of larger nuclear charge in Ni and thus less delocalized electronic cloud in Ni. This leads to stabilization of the ferromagnetic exchange interaction of Fe and thus more probability of the HS state surrounded by more Ni atoms. Furthermore, the Ni-Ni interaction is noticeably softer than the Ni-Fe one and is more likely to match the lattice parameter. These effects consequently yield smaller thermal expansion in the Ni-Ni pair than the Ni-Fe one.
Let us also discuss the third-order anharmonicity in the Invar alloy. Figure 7 shows C 3 for the average first NN shells around Fe and Ni, respectively. The agreement between the experiments and the PIECP simulations is not perfect, but C 3 increases gradually with the increase in temperature for both the Fe and Ni data. Clear anharmonicity in the Invar alloy is confirmed by comparing the PIECP simulations between the two-state (HS + LS) and HS-only models. Both the simulated data exhibit essentially the same C 3 values, implying no suppression of C 3 due to the contribution of the LS state, as observed in the thermal expansion. Since the asymmetric radial distribution for the first NN shell almost exclusively originates from the anharmonic interatomic potential, the present result implies that the third-order anharmonicity clearly exists even in the case of no thermal expansion. Note that in the present PIECP simulations within the low  [44] .
coupling approximation, the quantum effect in C 3 is not properly taken into consideration because of the neglect of the phonon-phonon coupling. The important qualitative finding of the presence of C 3 without thermal expansion is, however, clearly exemplified.
At the end of this section, the interatomic potentials around Fe are summarized in Figure 8 to compare between the three different Fe-based alloys of the Invar (Fe 65 Ni 35 ), SUS316L (Fe 67 Cr 18 Ni 12 Mo 2 Mn 1 ) [47] , and Elinvar (Fe 50 Ni 42 Cr 6 Ti 4 ) [34] alloys. The stainless steel SUS316L alloy exhibits normal thermal expansion that is quantitatively almost the same as that of SUS304 in Figure 4A. The Elinvar (Ni Span C) alloy shows intermediate thermal expansion as shown in Figure 4A due to the partial Invar effect. In the Invar alloy, the HS state of Fe is located below the LS state, ΔE ≡ E LS -E HS = 0.034 (eV), and the LS state gives a shorter interatomic distance of R LS = 0.2512 (nm) and a stiffer bulk modulus of B LS = 182 (GPa) than those of the LS state with R HS = 0.2553 (nm) and B HS = 147 (GPa) (Note here that the potential parameters are slightly different from those in Figure 4). Since the energy difference between the HS and LS states is not so large compared to the thermal energy, the population of the LS state increases gradually with the temperature rise, yielding the compression of the interatomic distance that compensates for the thermal expansion due to anharmonicity and the enlargement of the bulk modulus. The enlargement of the elastic constant is an interesting peculiar character compared to normal metals that show the softening of the elastic constant.
On the contrary, SUS316L shows that the HS state of Fe is located DE = -0.095 (eV) above the LS state, and the LS state gives a shorter interatomic distance of R LS = 0.2544 (nm) and a stiffer bulk modulus of B LS = 177 (GPa) than those of the HS state with R HS = 0.2549 (nm) and B HS = 143 (GPa). A much larger energy difference between the LS and HS states does not induce noticeable thermal excitation, leading to normal thermal expansion and elasticity softening of SUS316L at high temperature. In case of the Elinvar alloy, the situation is more complicated. There exist two types of Fe atoms, the one of which gives more unstable ] as in the case of the SUS316 alloy. These findings imply that some Fe shows the Invar effect, while the other does the anti-Invar effect, resulting in an incomplete Invar effect in the Elinvar alloy. Such a different behavior in the Fe atoms is attributed to different coordination around Fe. Fe is surrounded by 35% Ni in the Invar alloy, and by 18% Cr and 12% Ni in the SUS316L. In the EAM scheme, surrounding Cr is likely to stabilize the LS state of Fe, while Ni is to stabilize the HS state of Fe, because the Cr 3d electrons are much more delocalized and interact more strongly with Fe 3d than the Ni one due to the difference in the nuclear charge. In the Elinvar alloy, Fe is surrounded by 42% Ni, 6% Cr and 4% Ti, yielding coexistence of more stable HS and LS states in Fe, depending on the coordination of Ni or Cr and Ti.  [44] .
Other related works concerning the local thermal expansion and the elastic constant using experimental EXAFS analysis combined by the PIECP computational simulation can be found in literatures: solid Kr [48] , martensitic Mn 88 Ni 12 alloy [49] and (Zr,Nb)Fe 2 [50] . At the end of this section, a sophisticated hot article published very recently is cited [51] , which gives the experimental EXAFS results of the FeNi Invar alloy at extremely high pressure including the magnetic phase transition.

ZINC BLENDE
Most zinc blende and diamond structure materials show NTE at low temperature, except diamond. The mechanism of the NTE in these tetrahedral materials is rather simple and seems well established already. Let us consider the transverse acoustic phonon mode around the X point in the first Brillouin zone in diamond or zinc-blend structure, as shown in Figure 9. Assuming that the phonon propagating direction is along the [001] axis and the vibrational direction of atoms is [110], the mode yields a vibrational motion between atoms 1 and 2 perpendicular to the bond direction (see Figure 9B). In this situation, the equilibrium distance between atoms 1 and 2, R 0 in the left of Figure 9C, is forced to be the shortest bond length. This lattice configuration is unstable because the force between atoms 1 and 2 is always attractive as long as the equilibrium distance corresponds to the natural length of the spring, consequently leading to lattice contraction to provide the most stable lattice constant depending on the vibrational amplitude of the phonon mode as in the right panel of Figure 9C. This NTE model has extensively been discussed since the first discovery of NTE in Si [5] . More recently, the first-principles DFT calculations concerning the NTE of Si were performed [52] . Although in the electronic structure calculations the Si atoms move classically, they applied the quantum mechanical quasi-harmonic approximation to the lattice dynamics and elucidated a negative Grüneisen constant and a resultant NTE at low temperature. It should be noted here that the mechanism of the NTE is purely classic, although the term of phonon is used for convenience, because the lattice instability due to the vibration perpendicular to the chemical bond, as discussed above, can be described completely within the classical Newton equation of motion in the lattice dynamics theory.
More recently, Dalba and Fornasini group performed systematic investigations concerning NTE in tetrahedral and other systems using EXAFS spectroscopy quite in detail [53][54][55][56][57] . They observed excellently significant discrepancies between the interatomic distance given by EXAFS and the equilibrium distance (lattice constant) given by XRD, where the interatomic distances show quite normal thermal expansion as in ordinary materials with positive thermal expansion. This discrepancy is however straightforwardly understandable from the above discussion, because in the presence of sufficiently large vibrational amplitude perpendicular to the chemical bond the interatomic distance is always larger than the equilibrium distance, as seen in Figure 9C. They successfully derived the vibrational amplitude perpendicular to the bond from the interatomic distances, which was found to be several times larger than the one parallel to the bond. The origin of the NTE in tetrahedral systems is well understood more quantitatively through the experiments. alloys. FM and NM denote the ferromagnetically coupled high-spin and nonmagnetic low-spin Fe, respectively [47] .
In this study, to realize the NTE of zinc-blend structure computationally, classical and quantum mechanical simulations have been performed for CdTe. There have been no reports concerning classical or quantum mechanical lattice dynamics simulations for the NTE in diamond or zinc-blend structure to the best of the author's knowledge. As discussed above, the mechanism of the NTE can thoroughly be described within the classical dynamics theory, and the NTE in CdTe should be reproduced by the classical MC simulations. Classical and PIECP MC simulations have been performed under the NPT condition with the total number of atoms N = 512 (4 3 cubic unit cells). The bond-order potential (BOP) [58,59] of CdTe reported recently in the literature [60] is sufficiently reliable for the present investigation and was used without any modification. Figure 10 shows temperature dependence of the calculated interatomic distances in CdTe, together with the experimental results of EXAFS [55,56] and XRD [46] . The distances are scaled to match the first NN Cd-Te one; for instance, the lattice constant is multiplied by √3/4. The XRD result shows the NTE of the lattice constant between 20 and 80 K, while the EXAFS result exhibits rather normal thermal expansion with larger temperature dependence even at higher temperature. The calculated result using the classical dynamics provides almost linearly increasing thermal expansion. Although the classical dynamics may expect the NTE in CdTe as mentioned above, it is clarified that the NTE cannot be realized within the classical theory. On the contrary, the PIECP result exhibits the NTE at low temperature, implying semiquantitative agreement with the XRD experiment. Moreover, the calculated first-NN interatomic distance shows rather normal thermal expansion, which is also in semiquantitative agreement with the EXAFS experiment. As the shell is going to higher neighbors, the temperature dependence is gradually suppressed and approaches the lattice behavior. It is strikingly concluded that vibrational quantum fluctuations are essentially required to realize the NTE in CdTe, although the NTE mechanism is purely classic. It is here noted that similar simulations have been performed for diamond and CdSe, which are found to exhibit positive and negative thermal expansions, respectively, also in qualitative agreement with the experimental findings. Figure 11 shows temperature dependence of the mean square relative displacements for the first to third NN shells in CdTe obtained by the lattice dynamics calculations within the harmonic approximation, together with the EXAFS results. Both the components parallel and perpendicular to the corresponding atom pairs are depicted (parallel component corresponds to C 2 ). It is found that in the first NN Cd-Te shell the perpendicular MSRD is much larger than the parallel one due to the presence of strong chemical bond between Cd and Te, semiquantitatively consistent with the experimental findings. The large perpendicular MSRD may result in normal local thermal expansion for the first NN Cd-Te distance in spite of the NTE for the lattice constant. The second NN shell exhibits smaller difference between the parallel and perpendicular components due to the absence of the direct chemical bonds in the second NN shell. Note here that although there exist two kinds of the second NN shells of Cd-Cd and Te-Te, the difference is found to be negligibly small and only the Cd-Cd calculated results are depicted in Figure 11. The third Figure 10. Calculated temperature dependence of the distances in CdTe, together with the experimental XRD (light blue solid circle) [46] and EXAFS (orange open circle) [55,56] results. The lattice constant (black), the first (red), second (green), third (blue), and fourth (magenta) NN distances are depicted, obtained by the PIECP (solid) and classical (dashed) MC simulations. All the distances are scaled to match the first NN distances. Figure 11. (A) Temperature dependence of the first-NN MSRD parallel (red) and perpendicular (blue) to the Cd-Te bond, obtained by the lattice dynamics calculations within the quantum mechanical (solid lines) and classical (dashed lines) harmonic approximations, together with the experimental EXAFS (symbols) results [55,56] . (B) Temperature dependence of the second (green) and third (blue) NN MSRD parallel (solid lines) and perpendicular (dashed lines) to the atomic pair direction, together with the experimental EXAFS (symbols) results for the parallel components. The perpendicular component is divided by 2 to compare with the parallel one directly.
NN shell gives almost identical MSRD between the parallel and perpendicular components, implying no anisotropic correlated motions in higher NN shells because of almost no interactions in the third NN shell.
The most striking consequence in this work is that vibrational quantum fluctuation is essentially important to realize the NTE, although the mechanism of the NTE is well established and is regarded as a purely classical origin. In the classical model, normal vibrational modes such as stretching modes with higher frequencies that contribute to positive thermal expansion dominate the lattice potential, because the vibrational excitations take place irrespective of the vibrational frequencies in the classical theory. The vibrational amplitude perpendicular to the chemical bond, which is a key role of the NTE and was experimentally given by the previous sophisticated EXAFS works, is successfully reproduced in the present computational study. Other thermal properties such as bulk moduli and higher-order anharmonicity (EXAFS cumulants) show quite normal behaviors, which are also reasonably and reliably understood in the present simulations.

CONCLUSIONS
In this article, fundamental quantum mechanical theories of thermal expansion were at first described in the quantum statistical perturbation and path-integral frameworks. The EXAFS cumulants that characterize thermal vibration and thermal expansion are formulated using the Br K-edge EXAFS of diatomic Br 2 as a practical experimental example. The importance of the quantum effect is exemplified for the description of the low-temperature thermal properties such as thermal expansion. Subsequently, zero thermal expansion of FeNi Invar and negative thermal expansion of zinc-blende CdTe were discussed from the view point of local structure. These two examples exhibit anomalous thermal expansion due to different origins of the so-called Invar effect (temperature variation of electronic structure) and vibrational anomaly. The local thermal expansion is discussed based on experimental EXAFS spectroscopic data and the path-integral effective classical potential (PIECP) Monte Carlo computational simulations.
In the FeNi Invar alloy, in which zero thermal expansion is seen over a wide temperature range including very low temperature, the vibrational quantum effect is found to be essentially important; quantum mechanical thermal expansion tends to be zero at the absolute temperature, while classic mechanically thermal expansion is constant, irrespective of temperature. It is also noted that the interatomic distances and thermal expansion of different atomic pairs (Fe-Fe, Fe-Ni, and Ni-Ni in the FeNi Invar alloy) are slightly but meaningfully different from each other, although the alloy exhibit a simple fcc crystal. Moreover, the third-order anharmonicity is found to be present as in normal metals, although usually the third-order anharmonic force constant is believed to be proportional to the thermal expansion coefficient. In the zinc blende and diamond structures, the origin of negative thermal expansion in these tetrahedral crystals is known as a result of classical vibrational anomaly within the Newton dynamics theory. The quantum statistical simulation is however found to be essential to reproduce the negative thermal expansion of CdTe.
Finally, let us briefly discuss some perspectives on future directions of EXAFS spectroscopy for the investigation of thermal properties of the materials such as ZTE and NTE. Since the synchrotron radiation techniques have substantially been developed in recent years and are still being exploited, especially for the three-dimensional spectromicroscopic computer-tomographic imaging method. It is essentially important to investigate structures within single domains and also their grain boundaries not only in the magnetic ZTE and NTE materials but also in other functional materials including multi-phase materials that exhibit interesting thermal and other physical/chemical properties. The spatial resolution of the computertomographic EXAFS imaging may reach less than 10 nm, leading to visualization of three-dimensional crystal and local atomic-scale structure of functional materials in near future.