First-principles study of phonon anharmonicity and negative thermal expansion in ScF3

The microscopic origin of the large negative thermal expansion of cubic scandium trifluorides (ScF3) is investigated by performing a set of anharmonic free-energy calculations based on density functional theory. We demonstrate that the conventional quasiharmonic approximation (QHA) completely breaks down for ScF3 and the quartic anharmonicity, treated nonperturbatively by the self-consistent phonon theory, is essential to reproduce the observed transition from negative to positive thermal expansivity and the hardening of the R4+ soft mode with heating. In addition, we show that the contribution from the cubic anharmonicity to the vibrational free energy, evaluated by the improved self-consistent phonon theory, is significant and as important as that from the quartic anharmonicity for robust understandings of the temperature dependence of the thermal expansion coefficient. The first-principles approach of this study enables us to compute various thermodynamic properties of solids in the thermodynamic limit with the effects of cubic and quartic anharmonicities. Therefore, it is expected to solve many known issues of the QHA-based predictions particularly noticeable at high temperature and in strongly anharmonic materials.


I. INTRODUCTION
Negative thermal expansion (NTE) materials are useful for technological applications and have been studied actively in the recent two decades [1].So far, various kinds of NTE materials have been discovered including metal oxides [2][3][4][5], metal fluorides [6][7][8] and metal-organic frameworks [9,10].Also, NTE materials provide a route for near-zero thermal expansion composite materials, which are desired in the field of high-precision measurements and semiconductor devices.The mechanism of NTE can be roughly categorized into two types: the vibrational effect and the non-vibrational effect [11].The former, which applies to many NTE materials, attributes NTE to existence of phonon modes having a large negative Grüneisen parameter [12].Indeed, many NTE materials comprise rigid octahedral or tetrahedral structure units, mainly formed by oxygen atoms, and their rotational modes often show the pressure-induced softening.
Scandium trifluorides (ScF 3 ) is an NTE material belonging to the family of metal trifluorides [13,14], whose structure is the ReO 3 -type "open perovskite" (Fig. 1).Greve et al. reported that the coefficient of the thermal expansion (CTE) of ScF 3 is strongly negative at low temperature [α l −10 ppm • K −1 at 200 K], and the NTE is persistent up to ∼ 1100 K [6].Although the modest NTE can also be observed in other ReO 3type structures [15,16], the magnitude of NTE is far lower than in ScF 3 .Moreover, ScF 3 is unique in that it does not show a structural phase transition even at 0.38 K [17], while the other 3d metal trifluorides transform into a rhombohedral structure with cooling [18].In ScF 3 , the R 4+ soft mode, relevant to the cubic-to-rhombohedral structural phase transition, softens with cooling but does not condensate in the low-temperature limit, as evidenced by recent inelastic X-ray scattering (IXS) experiments [18,19].
To elucidate the unique thermophysical properties of ScF 3 , accurate treatment of lattice dynamics including anharmonic effects is essential.Since the discovery of the large NTE in ScF 3 , several first-principles studies based on density functional theory have reported the CTE calculated within the quasi-harmonic approximation (QHA) [20][21][22].Since the QHA only accounts for the volume dependence of phonon frequencies and neglects higher-order anharmonicities, however, none of these QHA results correctly reproduced the experimental CTE values.Li et al. [20] have demonstrated that the inclusion of quartic anharmonicity of the R 4+ soft mode is significant and essential for explaining the experimental result.More recently, an improved treatment of anharmonic effects either by ab initio molecular dynamics (AIMD) [23] or by a stochastic implementation [24] of the self-consistent phonon (SCP) theory [25] has been shown to explain the observed CTE more quantitatively and highlighted the importance of higherorder anharmonicities.Despite the success of these numerical methods for reproducing the temperature dependence of CTE of ScF 3 somewhat quantitatively, the microscopic mechanism of the large NTE is not fully understood.In particular, the role of the quartic anharmonicity on the CTE has not been quantified in comparison with that of the cubic anharmonicity, which can also affect the vibrational free energy.Therefore, an alternative approach is desired.
In this paper, we present another first-principles approach for incorporating anharmonic effects in the vibrational free energy beyond the QHA level.Our approach is based on the recent implementation of the SCP theory that uses fourth-order interatomic force constant (IFC) [26].Since the SCP approach can compute anharmonic phonon quasiparticle non-perturbatively on a dense momentum grid, thermodynamic properties in the thermodynamic limit can be evaluated accurately, which is essential for accurate determination of phase boundaries [27].By applying the method to ScF 3 , we show that the inclusion of the quartic anharmonicity solves the issue of the QHA and reproduces the experimentally observed change from negative to positive CTE with heating, in accord with Ref. [20].In FIG. 1. Crystal structure of cubic ScF 3 with space group Pm 3m (created with [29]).Each scandium atom is surrounded by six fluorine atoms, which form corner-sharing octahedra.The structure is similar to the ABX 3 cubic perovskite, but the A site is vacant.addition, we find that the effect of the cubic anharmonicity, evaluated by the improved self-consistent (ISC) theory [28], is equally important and compensates the overcorrection made by the SCP theory with only the quartic anharmonicity.
This paper is organized as follows.First, we introduce the theoretical methods to incorporate anharmonic effects in vibrational free energy in Sec.II.We then describe computational and technical details of the DFT calculation in Sec.III.In Sec.IV B, we demonstrate the limitation of the QHA theory in ScF 3 , which can be resolved by more accurate treatment of anharmonic effects as shown in Sec.IV C.Moreover, we discuss the temperature dependence of the R 4+ soft mode in Sec.IV D. Finally, a concluding remark is made in Sec.V.

II. VIBRATIONAL FREE ENERGY
The Helmholtz free energy of a non-metallic system is expressed as a function of the volume of the unit cell V and temperature T, where E el (V) is the static internal energy of electrons obtained by a first-principles calculation and F vib (V, T) is the vibrational free energy.Assuming that the lattice vibration is well described by thermal excitation of phonons, we examine the behavior of F vib (V, T) with the three levels of approximation in this study: the QHA, SCP theory, and ISC theory.Once the free energy is obtained at various volumes, the F(V, T) curve is fitted by the equation of state (EOS) to estimate an equilibrium volume at each temperature.By repeating the procedure at different temperatures, we can obtain the volume thermal expansion α v (T) = 1 V ( ∂V ∂T ) P .For isotropic systems, the CTE is equal to 1  3 α v .For the sake of brevity, we use the shorthand notation of q = (q, j) and −q = (−q, j) where q is the phonon momentum and j is the phonon branch index.

A. Quasiharmonic theory
The quasiharmonic (QH) theory is the standard approximation of F vib (V, T), in which the vibrational free energy is given as Here, β = 1/kT with the Boltzmann constant k, and ω q (V) is the phonon frequency at volume V obtained within the harmonic approximation (HA).While the anharmonic effect is partially incorporated via the volume dependence of ω q (V), the QHA completely neglects the intrinsic anharmonic effects which are responsible for making the temperature dependence of phonon frequencies.Nevertheless, the QH theory turns out to be a good approximation at temperatures far below the melting point and has been employed to predict the thermal expansivity and phase boundary of various materials based on DFT.When the temperature reaches the melting point or the structure is strongly anharmonic, the QHA is less reliable.Moreover, the QHA is not valid for cases where phonon modes become unstable within the HA as will be discussed in Sec.IV B.

B. Self-consistent phonon theory
The SCP theory is one of the most successful approaches for calculating the temperature dependent phonon frequencies nonperturbatively [25,26].In this method, we assume the existence of effective harmonic phonon frequency Ω q (Ω 2 q ≥ 0) and polarization vector ε q , with which we can define an effective harmonic Hamiltonian H 0 as Here, A q = a q + a † −q is the displacement operator with a † q and a q being creation and annihilation operators of SCP, respectively.In the first-order SCP theory, the renormalized phonon frequencies and eigenvectors {Ω q , ε q } are determined so that the vibrational free energy within the first-order cumulant approximation is minimized.Let F vib denotes the vibrational free energy of an anharmonic system described by Hamiltonian

see Eq. (A2)] and F
vib denotes the free energy with the first-order cumulant approximation, the following inequality then holds: Here, vib is given as follows: where , . . ., ω 2 qs ) and α q = 2Ω q [1 + 2n(Ω q )] with n(ω) being the Bose-Einstein distribution function.Φ (SCP) (q; −q; q ; −q ) is the fourth-order IFC in the normal coordinate basis of the SCP.For the detailed expression, see appendix A. From the condition of ∂F vib ∂Ω q = 0 and ∂F vib ∂C † q = 0, we obtain the first-order SCP equation as Since the anharmonic interaction is relatively short-range compared with the harmonic one, it is sufficient to solve the above equation on an irreducible set of q points that are commensurate with an employed supercell.Once the change of the dynamical matrix ∆D(q) = D (SCP) (q) − D (HA) (q) is obtained on these grids, ∆D(q) can be interpolated to arbitrary q points by the Fourier interpolation [26].
In Eqs. ( 5) and ( 6), we neglected the sixth and higher-order anharmonic terms assuming that their effects are far smaller than the dominant fourth-order term.With this assumption, the third term in Eq. ( 5) can be removed by using Eq. ( 6) as Hence, evaluating the anharmonic free energy F (SCP) vib (V, T) in the thermodynamic limit is feasible owing to the interpolation technique, which is hardly achievable by the thermodynamic integration based on AIMD.The first term in Eq. ( 7) corresponds to the QH contribution [Eq.( 2)] but the harmonic frequency is replaced with the SCP frequency.The second term is a correction necessary to satisfy the correct thermodynamic relationship S = −dF/dT [30], where the entropy is given as S = k q [(n q +1) ln(n q +1)−n q ln n q ] with n q = n(Ω q ).When the anharmonic effect is small, we can assume Ω q ≈ ω q + ∆ q and C q,i j ≈ δ i j .Then, we obtain the small perturbation limit of Eq. ( 7) as This result is the same as the one reported by Allen [30].In this study, we always use Eq. ( 7) because it is more generalized and even applicable to systems where unstable (imaginary) phonon modes exist within the HA.

C. Improved self-consistent phonon theory
The SCP theory accounts for the quartic anharmonicity nonperturbatively but neglects the effect of the cubic anharmonicity.The ISC theory [28] accounts for the additional three-phonon term perturbatively as where F (B) vib (V, T) is the Helmholtz free energy from the bubble diagram [26], associated with the cubic anharmonicity, given as Here, we simply denote n(ω q i ) as n i and ∆(q) = N q δ q,mG with any integer m.In the calculation of F (B) vib (V, T), we use the SCP lattice dynamics wavefunction instead of those within the HA.Previous numerical studies using empirical model potentials showed that the ISC theory could describe various thermodynamic properties of copper [31] and noble-gas solids [28,32] accurately in a wide temperature range.In this study, we evaluate all the components entering the free energy with fully nonempirical first-principles calculation.The bubble free-energy is theoretically related to the bubble self-energy Σ (B)  q (iω m ) as where ω m = 2πm/β is the Matsubara frequency, is the SCP Green's function, and The ISC Green's function provides information of lattice dynamics in the same approximation level as Eq.(9).Equation ( 14) can be used to calculate the phonon spectral function with the effects of cubic and quartic anharmonicities, as demonstrated in cubic SrTiO 3 [33], and the imaginary part of Σ (B) q (Ω q ) is essential for thermal conductivity calculations [26,33,34].
Although the ISC theory is computationally more costly than the SCP theory because of the triplet (q 1 , q 2 , q 3 ) loop in Eq. (10), its computational complexity is the same as that of thermal conductivity calculations with three-phonon interactions, which has been applied to various materials so far including relatively complex ones.Therefore, we expect the ISC is feasible even for complex systems.Moreover, for high-symmetry structures such as ScF 3 , we can drastically reduce the computational cost by utilizing the symmetry of Φ(q; q 1 ; q 2 ) [35].

A. DFT calculation
We employed Vienna ab initio simulation package ( ) [36] for calculating the electron states of ScF 3 , which implements the projector augmented wave (PAW) [37,38] method.The adopted PAW potentials treat Sc 3s 2 3p 6 3d 1 4s 2 and F 2s 2 2p 4 as valence states [39].A kinetic energy cutoff of 700 eV and the 12×12×12 Monkhorst-Pack kpoint mesh were employed.The self-consistent field loop was continued until the total electronic energy change between two steps became smaller than 10 −8 eV.To investigate the influence of the exchange-correlation functional, we adopted the local-density approximation (LDA), Perdew-Burke-Ernzerhof (PBE) [40] functional, and a variant of PBE optimized for solids (PBEsol) [41].

B. Force constant calculation
To compare the vibrational free energy within the different levels of approximation, either the QH, SCP or ISC theory, it is necessary to extract second-, third-, and fourth-order IFCs from first-principles calculations.To this end, we employed the real-space supercell method with the 2×2×2 supercell (32 atoms), which is sufficiently large to allow the out-of-phase tilting motion of fluorine octahedra.To extract second-order IFCs, an atom in the supercell was displaced from its equilibrium site by 0.01 Å and the Hellmann-Feynman forces were calculated for the displaced configuration.From the data sets comprising the displacements and forces, we estimated the second-order terms by the least-squares fitting [42], as implemented in the [43,44] package.To check the convergence of the calculation, we changed the supercell size to 4×4×4 (256 atoms) and calculated the second-order terms.However, the calculated F (QH)  vib (V, T) was almost unaltered by this change, thus validating the use of the 2×2×2 supercell.
The anharmonic IFCs were extracted by using the compressive sensing lattice dynamics method [45].Following the prescription given in Refs.26 and 45, we generated 60 random displacement patterns from the trajectory of ab initio molecular dynamics at 500 K and calculated atomic forces for these patterns.We then fitted the Taylor expansion potential (TEP) [26] by using the least absolute shrinkage and selection operator (LASSO) with a regularization parameter optimized through the cross-validation.The adopted TEP includes anharmonic terms up to the sixth-order; the fourth-order IFCs are restricted to on-site, two-body, and three-body terms and the higher-order IFCs are restricted to on-site and two-body terms.We also verified the accuracy of the obtained TEP for independent test data sets.In the LASSO regression step, we fixed the second-order IFCs to the precalculated values and optimized the anharmonic terms only.

C. Calculation of phonons and CTE
The SCP equation [Eq.( 6)] was solved using a numerical algorithm of Ref. [26], as implemented in the code.The q mesh of the SCP was set to 2×2×2, which is commensurate with the supercell size, and the inner q mesh was increased up to 8×8×8 to achieve convergence of anharmonic phonon frequencies.After the solution to the SCP equation was found, we converted the effective dynamical matrices at the commensurate q points into the real-space ef- fective second-order IFCs, which was then used to calculate anharmonic phonon frequencies at a denser q-point grid for Eqs. ( 7) and (10).For the calculation of the vibrational free energy, we employed 20×20×20 q point, which was sufficient to reach convergence.In all of the phonon calculations, the non-analytic correction to the dynamical matrix was included by using the Ewald's method [46] with the dielectric tensor and Born effective charges obtained from density functional perturbation theory (DFPT) [47] calculations.
To calculate the CTE of ScF 3 , we changed the lattice constant from the optimized value from −0.05 Å to +0.02 Å in steps of 0.005 Å and repeated the procedure described above and in Sec.III B. The lattice constants and Bulk modulus at finite temperature were estimated by fitting the obtained F(V, T) curves with the Birch-Murnaghan EOS [48], as implemented in [49].

A. Ground state structural property
Table I shows the equilibrium lattice constant a 0 and bulk modulus of ScF 3 calculated from E el (V), which are compared with previous experimental and computational results.All of the adopted functionals give a 0 values close to the experimental results of Greve et al. [6] within an error of ∼1.0%.Among the three functionals, PBEsol best reproduces the experimental result with deviation as small as 0.1%.In contrast, the bulk modulus estimated from E el (V) overestimates the experimental value at 300-500 K [50], in accord with the previous computational studies [21,23].This indicates the essential role of lattice vibration in determining B 0 , which will be discussed further in Sec.IV C.

B. Volume dependent phonon frequency and breakdown of QHA
In this section, we demonstrate the breakdown of the QHA in ScF 3 .Figure 2  high-symmetry lines calculated with LDA, PBE, and PBEsol at volumes close to their equilibrium values.In ScF 3 , two soft modes exist: the triply-degenerate R 4+ mode and the M 3+ mode.Along the high-symmetry line M-R, the soft mode shows almost no dispersion.Also, the soft mode is unstable at the equilibrium volume within LDA.

compares phonon dispersion curves along
The frequencies of the soft modes sensitively changes with the volume as shown in Fig. 3.The frequency of the R 4+ soft mode increase sharply with increasing the cell volume, which represents a large negative mode Grüneisen parameter γ q = − V ω q ∂ω q ∂V as pointed out in the previous theoretical work [20].Assuming that all phonon modes are stable, the CTE is given as where C v , B T , γ, and c q are the isothermal heat capacity, isothermal bulk modulus, the average Grüneisen parameter, and the mode-specific heat, respectively.Therefore, the NTE originates from a negative γ, which in turn can be attributed to the large negative γ q of the low-frequency soft mode whose thermal weight c q is relatively large in a low-temperature region.
To estimate the CTE quantitatively, we need to calculate Eq. ( 15) or (2) from phonon frequencies.When an unstable mode exists, however, these equations must not be used.One may naively expect that using F (QH)  vib (V, T) only with stable phonon modes gives a reasonable estimate.However, such a simple treatment produces an irregular oscillation of F (QH)  vib (V, T) curve which affects the accuracy of the EOS fitting.The same discussion has also been made by Lan et al. [51].Therefore, within the QH theory, we argue that the EOS fitting should be performed only with the volumes where all phonon modes are stable.By this principle, we calculated the CTE with the three functionals.As shown in Fig. 4, the FIG. 4. CTE calculated within the QHA compared with the experimental results [6] and the previous QHA results [20,21].
CTE as well as the lattice constant keep decreasing with increasing the temperature.Above ∼ 300 K, the EOS fitting failed because of the lack of valid free-energy data points in the small volume region.In Fig. 4, we also compare our results with the previous QHA results [20,21], showing noticeable disagreement.While we cannot identify the origin of the disagreement, it can likely be attributed to the difference in the treatment of unstable phonon modes and the details of DFT calculation.Indeed, when we used F (QH) vib (V, T) calculated from stable phonon modes only, we obtained the temperature dependence of CTE similar to the result of Ref. 21.

C. Inclusion of anharmonic free energy
The SCP theory can circumvent the limitation of the QHA mentioned above because it can stabilize the unstable soft mode by the Debye-Waller type renormalization.Since the SCP assumes the existence of well-defined phonon modes (Ω 2 q > 0), the SCP equation ( 6), if converged, always gives stable phonons.Therefore, F (SCP)  vib (V, T) can be evaluated even when an unstable mode exists within the harmonic approximation.Besides, it is straightforward to include the additional correction term from the cubic anharmonicity by evaluating Eq. ( 10) on top of the SCP solution.
To quantify the effect of quartic and cubic anharmonicities, we first investigate the volume and temperature dependence of the anharmonic free energies as shown Fig. 5.The anharmonic free energy correction based on the SCP is defined as which can be defined only when all phonon modes are stable within the harmonic approximation.As shown in Fig.
(5)(a), the ∆F (SCP) vib (V, T) value of ScF 3 is positive and decreases gradually as increasing the cell volume.Therefore, a larger unitcell volume system is relatively more stabilized by the quartic anharmonicity.On the other hand, the correction from the bubble diagram F (B)  vib (V, T) is negative and shows an opposite volume dependence.The sign of the total correction ∆F (SCP)  vib (V, T) + F (B) vib (V, T) is positive in most of the studied volumes, but it becomes negative in the large volume systems as noticeable in the PBE results.Our result thus shows the importance of the anharmonic correction not only from the quartic anharmonicity but also from the cubic one, which competes with each other.
Figure 6 shows the temperature-dependent CTE and bulk modulus calculated within the SCP and ISC theories.The SCP theory correctly reproduced the positive sign of the CTE observed in the high-temperature region.However, the temperature at which the sign change occurs was around 400-500 K, which underestimates the experimental value of ∼ 1100 K.This underestimation was partially cured by the correction from the bubble diagram, and the sign change temperature became ∼ 600-800 K. Considering the fact that the CTE is sensitive to the exchange-correlation functional and pseudopotentials, the agreement between our ISC results with the experimental CTE of Greve et al. [6] is reasonable.In addition, we successfully reproduced the isothermal bulk modulus of Morelock et al. [50] semiquantitatively including its weak temperature dependence as shown in Fig. 6(b).

D. Temperature dependence of soft mode frequency
Finally, we discuss the influence of quartic and cubic anharmonicities on the frequency of the R 4+ soft mode.Figure 7 shows the temperature dependence of the squared phonon frequency of the R 4+ soft mode calculated with the PBEsol functional.The results shown in Fig. 7(a) were calculated with the equilibrium lattice constant of the ISC calculation, and the ISC phonon frequency was calculated by Ω (ISC) q = (Ω 2 q + 2Ω q ∆ (B) q (Ω q )) 1 2 where ∆ (B) q (ω) = −ReΣ (B)  q (ω).The harmonic and anharmonic IFCs at these new volumes were obtained by linearly interpolating the values calculated at the nearest two volumes.As shown in the figure, the QHA frequency softens as a function of temperature because of the NTE and the negative Grüneisen parameter.This temperature dependence is opposite to the IXS results [18,19].The hardening of the soft mode frequency can be reproduced only by the SCP and ISC theories.Unfortunately, however, the SCP and ISC phonon frequencies systematically overestimate the experimental values.Since the frequency is sensitive to the adopted volumes, the deviation may be mitigated by improving the agreement of the PBEsol and experimental lattice constants.To examine this point, we adjusted the external pressure P in such a way that the PBEsol lattice constant agrees with the experimental result at 0 K, resulting in P ≈ 0.7 GPa.Then, the agreement was improved considerably as shown in Fig. 7(b).In the low-temperature limit, the ISC result best agrees with the IXS results, while the SCP theory well explains the slope of the hardening.The IXS data lies in between the SCP and ISC results in Fig. 7(b).This may indicate the necessity of fully self-consistent treatment of the cubic and quartic anharmonicities, which is left for a future study.
It is interesting to observe in Fig. 7 (a) and (b) that the SCP frequency is smaller than the QHA result near ∼ 0 K, meaning that the quartic anharmonicity decreases the soft mode frequency.This behavior seems to contradict with the large positive quartic potential of the R 4+ mode [20], which must harden the frequency.To understand this unusual behavior, we investigated the mode-dependent contribution to the quartic renormalization Φ (SCP) (q; −q; q ; −q )α q [see Eq. ( 6)] for q = R 4+ and q on the 8×8×8 uniform grid.(q; −q; q ; −q )α q (meV 2 shows the mode dependence of Φ (SCP) (q; −q; q ; −q )α q at 0 K.The quartic coupling is strongly positive in the lowfrequency phonon modes, which are the rigid unit motion of the fluorine octahedra, in accord with the previous numerical evaluation [20].In contrast, the coupling coefficients are negative for high-frequency optical modes in 40-80 meV, which are mostly dominated by the vibration of fluorine atoms as shown in the inset of the figure.Since the total negative contribution is slightly larger than the total positive one, the SCP frequency becomes smaller than that of the QHA at 0 K.With increasing the temperature, the factor α q increase more rapidly for lowfrequency phonon mode, and the total contribution becomes positive around 60 K [see Fig. 7(a)].We note that this unique phenomenon results from the zero-point motion of atoms.If we replace n q with n q − 1 2 in Eqs. ( 6) and ( 12), the SCP frequency is larger than the QHA frequency of the temperature as shown in Fig. 7(c).

V. CONCLUSION
We investigated the role of the cubic and quartic anharmonicities in the large NTE of cubic ScF 3 by using many-body theoretical approaches based on first-principles anharmonic force constants.We showed that the quartic anharmonicity, which was included by the SCP theory, is essential to reproduce the observed transition from negative to positive CTE as heating, which is consistent with the previous work [20].In addition, we newly found that the inclusion of the cubic anharmonicity, which was achieved by the ISC theory, improves the quantitative agreement between theoretical and experimental thermal expansivities.Therefore, both cubic and quartic anharmonicities are equally important for a quantitative understanding of the unusually large NTE of ScF3 3 .Moreover, we calculated the temperature-dependence of the R 4+ soft mode within the SCP and ISC theories and obtained results that agree semiquantitatively with the IXS data, particularly when the experimental lattice constant is employed.We showed that the R 4+ mode softens by the quartic anharmonicity near 0 K, and revealed that it results from the zero-point vibration of atoms and the negative quartic coupling between the R 4+ mode and high-frequency optical modes.
The present first-principles approach enables us to obtain vibrational free-energy of solids with cubic and quartic anharmonicities and even applicable to systems where unstable phonon modes exist within the harmonic approximation (e.g.high-temperature phase of solids).Since the method is based on the reciprocal space formalism, thermodynamic properties in the thermodynamic limit can be calculated efficiently by using the interpolation technique, which is particularly important for a computational estimation of CTE and phase boundaries.

FIG. 3 .
FIG. 3. Volume dependence of the phonon frequency of the R 4+ soft mode.The imaginary frequencies are shown as negative values.The filled symbols indicate the volumes corresponding to the results shown in Fig. 2.

FIG. 6 .
FIG.6.(a) CTE and (b) bulk modulus calculated within the SCP and ISC theories compared with the experimental results of Greve et al.[6] and Morelock et al.[50].

FIG. 7 .
FIG. 7. Temperature dependence of squared phonon frequency of the R 4+ soft mode calculated within PBEsol in comparison with the experimental result of Handunkanda et al. [18].(a) The calculation at each temperature is conducted with the corresponding ISC volume.(b) The calculation is conducted with the lattice constant adjusted to match the experimental lattice constant at 0 K. (c) Same as (b) but the classical approximation is made.

)FIG. 8 .
FIG.8.Mode-dependent contribution to the quartic renormalization of the R 4+ soft mode at 0 K.The open circle shows the contribution from phonon mode q and the dotted line shows its cumulative value.The inset shows the displacement patterns of the phonon modes having large positive and negative quartic interaction with the R 4+ soft mode.

TABLE I .
Comparison of calculated and experimental lattice constant and bulk modulus of ScF 3 .