High throughput thermal conductivity of high temperature solid phases: The case of oxide and fluoride perovskites

Using finite-temperature phonon calculations and machine-learning methods, we calculate the mechanical stability of about 400 semiconducting oxides and fluorides with cubic perovskite structures at 0 K, 300 K and 1000 K. We find 92 mechanically stable compounds at high temperatures -- including 36 not mentioned in the literature so far -- for which we calculate the thermal conductivity. We demonstrate that the thermal conductivity is generally smaller in fluorides than in oxides, largely due to a lower ionic charge, and describe simple structural descriptors that are correlated with its magnitude. Furthermore, we show that the thermal conductivities of most cubic perovskites decrease more slowly than the usual $T^{-1}$ behavior. Within this set, we also screen for materials exhibiting negative thermal expansion. Finally, we describe a strategy to accelerate the discovery of mechanically stable compounds at high temperatures.


I. INTRODUCTION
High throughput ab-initio screening of materials is a new and rapidly growing discipline [1]. Amongst the basic properties of materials, thermal conductivity is a particularly relevant one. Thermal management is a crucial factor to a vast range of technologies, including power electronics, CMOS interconnects, thermoelectric energy conversion, phase change memories, turbine thermal coatings and many other [2]. Thus, rapid determination of thermal conductivity for large pools of compounds is a desirable goal in itself, which may enable the identification of suitable compounds for targeted applications. A few recent works have investigated thermal conductivity in a high throughput fashion [3,4]. A drawback of these studies is that they were restricted to use the zero kelvin phonon dispersions. This is often fine when the room temperature phase is mechanically stable at 0 K. It however poses a problem for materials whose room or high temperature phase is not the 0 K structure: when dealing with structures exhibiting displacive distortions, including temperature effects in the phonon spectrum is a crucial necessity.
Such a phenomenon often happens for perovskites. Indeed, the perovskite structure can exhibit several distortions from the ideal cubic lattice, which is often responsible for rich phase diagrams. When the structure is not stable at low temperatures, a simple computation of the phonon spectrum using forces obtained from density functional theory and the finite displacement method yields imaginary eigenvalues. This prevents us from assessing the mechanical stability of those compounds at high temperatures or calculating their thermal conductivity. Moreover, taking into account finite-temperature effects in phonon calculations is currently a very demanding task, especially for a high-throughput investigation.
In this study, we are interested in the high-temperature properties of perovskites, notably for thermoelectric applications. For this reason, we focus on perovskites with the highest symmetry cubic structure, which are most likely to exist at high temperatures [5][6][7][8][9]. We include the effects of anharmonicity in our ab-initio calculations of mechanical and thermal properties.

II. FINITE-TEMPERATURE CALCULATIONS OF MECHANICAL STABILITY AND THERMAL PROPERTIES
Recently, several methods have been developed to deal with anharmonic effects at finite temperatures in solids [10][11][12][13][14][15]. In this study, we use the method presented in Ref. [15] to compute the temperature-dependent interatomic force constants, which uses a regression analysis of forces from density functional theory coupled with a harmonic model of the quantum canonical ensemble. This is done in an iterative way to achieve self-consistency of the phonon spectrum. The workflow is summarized in Fig. 1. In the following (in particular Section III), it will be referred as "SCFCS" -standing for self-consistent force constants. As a trade-off between accuracy and throughput, we choose a 3x3x3 supercell and a cutoff of 5 Å for the third order force constants. Special attention is paid to the computation of the thermal displacement matrix [15], due to the imaginary frequencies that can appear during the convergence process, as well as the size of the supercell that normally prevents us from sampling the usual soft modes at the corners of the Brillouin zone (see Supplementary Material). This allows us to assess the stability at 1000 K of the 391 hypothetical compounds mentioned in Section I. Among this set, we identify 92 mechanically stable compounds, for which we also check the stability at 300 K. The phonon spectra of the stable compounds are provided in the Supplementary Material.  Furthermore, we compute the thermal conductivity using the finite temperature force constants and the full solution of the Boltzmann transport equation as implemented in the ShengBTE code [16].
We list the stable compounds and their thermal conductivities in Table I. Remarkably, this list contains 37 perovskites that have been reported experimentally in the ideal cubic structure (see References in Table I), which lends support to our screening method. On the other hand, we also find that 11 compounds are reported only in a non-perovskite form. This is not necessarily indicative of mechanical instability, but instead suggests thermodynamical stability may be an issue for these compounds, at least near this temperature and pressure. 36 compounds remain unreported experimentally in the literature to our knowledge. Thus, by screening only for mechanical stability at high-temperatures, we reduce the number of potential new perovskites by a factor of 10. Furthermore, we find that 50 of them are mechanically stable in the cubic form close to room temperature.
Of the full list of perovskites, only a few measurements of thermal conductivity are available in the literature. They are displayed in parentheses in Table I along with their calculated values. Our method tends to slightly underestimate the value of the thermal conductivity, due to the compromises we made to limit the computational cost of the study (see Supplementary Material). This dis-crepancy could also be partially related to the electronic thermal conductivity, which was not substracted in the measurements. Still, we expect the order of magnitude of the thermal conductivity and the relative classification of different materials to be consistent. More importantly, this large dataset allows us to analyze the global trends driving thermal conductivity. These trends are discussed in Section IV.
We also investigate the (potentially) negative thermal expansion of these compounds. Indeed, the sign of the coefficient of thermal expansion α V is the same as the sign of the weighted Grüneisen parameter γ, following where K T is the isothermal bulk modulus, c V is the isochoric heat capacity and ρ is the density [80,81]. The weighted Grüneisen parameter is obtained by summing the contributions of the mode-dependent Grüneisen Finally the modedependent parameters are related to the volume variation of the mode frequency ω i via γ i = −(V /ω i )(∂ω i /∂V ). In our case, we calculate those parameters directly using the second and third order force constants at a given temperature [12,82,83]: This approach has been very successful in predicting the thermal expansion behavior in the empty perovskite ScF 3 [15], which switches from negative to positive around 1100 K [84]. In our list of filled perovskites, we have found only two candidates with negative thermal expansion around room temperature: TlOsF 3 and BeYF 3 , and none at 1000 K. This shows that filling the perovskite structure is probably detrimental to the negative thermal expansion.
We also examine the evolution of the thermal conductivity as a function of temperature, for the compounds that are mechanically stable at 300 K and 1000 K. There is substantial evidence that the thermal conductivity in cubic perovskites generally decreases more slowly than the model κ ∝ T −1 behavior [85,86] at high temperatures, in contrast to the thermal conductivity of e.g. Si or Ge that decreases faster than κ ∝ T −1 [87]. This happens for instance in SrTiO 3 [47,48] [88] and BaZrO 3 [53]. We also predicted an anomalous behavior in ScF 3 using ab-initio calculations, tracing its origin to the important anharmonicity of the soft modes [15]. Fig. 2 displays several experimentally measured thermal conductivities from the literature on a logarithmic scale, along with the results of our high-throughput calculations. As discussed above, the absolute values of the calculated thermal conductivities are generally underestimated, but their relative magnitude and the overall temperature dependence are generally consistent. Although the behavior of the thermal conductivity κ(T ) is in general more complex than a simple power-law behavior, we a AuMgF 3 was mentioned theoretically in Ref. [46]. b The thermal diffusivity of BaLiF 3 was measured at 300 K in Ref. [67] as α=0.037 cm 2 s −1 . Table I. List of cubic perovskites found to be mechanically stable at 1000 K and their corresponding computed lattice thermal conductivity (in W/m/K). We also report the computed lattice thermal conductivity at 300 K (in W/m/K) when we obtain stability at that temperature. We highlight in blue the compounds that are experimentally reported in the ideal cubic perovskite structure, and in red those that are reported only in non-perovskite structures (references provided in the table). When no reference is provided, no mention of the compound in this stoichiometry has been found in the experimental literature. Experimental measurements of the thermal conductivity are reported in parentheses, and in italics when the structure is not cubic.
model the deviation to the κ ∝ T −1 law by using a parameter α that describes approximately the temperaturedependence of κ between 300 K and 1000 K as κ ∝ T −α . For instance, in Fig. 2, KMgF 3 appears to have the fastest decreasing thermal conductivity with α = 0.9 both from experiment and calculations, while SrTiO 3 is closer to α = 0.6. At present, there are too few experimental measurements of the thermal conductivities in cubic perovskites to state that the κ ∝ T −α behavior with α < 1 is the general rule in this family. However, the large number of theoretical predictions provides a way to assess this trend. Of the 50 compounds that we found to be mechanically stable at room temperature, we find a mean α 0.85, suggesting that this behavior is likely general and correlated to structural characteristics of the perovskites.

III. ACCELERATING THE DISCOVERY OF STABLE COMPOUNDS AT HIGH TEMPERATURE
Through brute-force calculations of the initial list of 391 compounds, we extracted 92 that are mechanically stable at 1000 K. However, this type of calculation is computationally expensive. Thus, it is desirable for future high-throughput searches of other material classes to define a strategy for exploring specific parts of the full combinatorial space. In this section, we propose and test such a strategy based on an iterative machine-learning scheme using principal component analysis and regression.
faster than finite-temperature calculations. This gives us a list of 29 perovskites that are mechanically stable in the cubic phase at 0 K. Since this is the highest symmetry phase, they are likely also mechanically stable at high-temperatures [89]. We calculate their self-consistent finite-temperature force constants Φ SCF CS 1000 K as described in Section II. This initial set allows us to perform principal component analysis of the 0 K force constants so that we obtain a transformation that retains the 10 most important components. In a second step, we use regression analysis to find a relation between the principal components at 0 K and at 1000 K. This finally gives us a model that extracts the principal components of the force constants at 0 K, interpolate their values at 1000 K, and reconstruct the full force constants matrix at 1000 K Φ model 1000 K . We say that this model has been "trained" on the particular set of compounds described above. Applying it to the previously calculated Φ 0 K for all compounds, we can efficiently span the full combinatorial space to search for new perovskites with a phonon spectrum that is unstable at 0 K but stable at 1000 K. For materials determined mechanically stable with Φ model 1000 K , we calculate Φ SCF CS 1000 K . If the mechanical stability is confirmed, we add the new compound to the initial set and subsequently train the model again with the enlarged set. When no new compounds with confirmed mechanical stability at high temperatures are found, we stop the search. This process is summarized in Figure 3. Following this strat-If unstable, draw new candidate If really stable, add to the set Train PCA and regression of force constants from 0 K to 1000 K Set of stable compounds at 1000 K and force constants at 0 K and 1000 K Perform full calculation to verify stability and obtain the force constants at 1000 K Find new candidate for stability at 1000 K using the model Figure 3. Depiction of strategy for exploring the relevant combinatorial space of compounds that are mechanically stable at high temperature.
egy, we find 79 perovskites that are stable according to the model, 68 of which are confirmed to be stable by the full calculation. This means that we have reduced the total number of finite-temperature calculations by a factor of 5, and that we have retrieved mechanically stable compounds with a precision of 86% and a recall of 74% [90]. It allows us to obtain approximate phonon spectra for unstable compounds, which is not possible with our finite-temperature calculations scheme (see Supplementary Material). It also allows us to find compounds that had not been identified as mechanically stable by the first exhaustive search due to failures in the workflow. Considering the generality of the approach, we expect this method to be applicable to other families of compounds as well. Most importantly, it reduces the computational requirements, particularly if the total combinatorial space is much larger than the space of interest.

IV. SIMPLE DESCRIPTORS OF THE THERMAL CONDUCTIVITY
We now focus on the analysis of the thermal conductivity data provided in Table I. We note that this set contains about two times more fluorides than oxides. This was already the case after the first screening in which we kept only the semiconductors, and it can be explained by the strong electronegativity of fluorine, which generally forms ionic solids with the alkali and alkaline earth metals easily, as well as with elements from groups 12, 13 and 14. This is shown on Fig. 4, in which we display histograms of the columns of elements at sites A and B of the perovskite in our initial list of paramagnetic semiconductors and after screening for mechanical stability.
We can also see that the oxides tend to display a higher thermal conductivity than the fluorides, as shown on the density plot of Fig. 5. This is once again due to the charge of the fluorine ion, which is half that of the oxygen ion. In a model of a purely ionic solid, this would cause the interatomic forces created by electrostatic interactions to be divided by two in fluorides as compared to oxides. This is roughly what we observe in our calculations of the second order force constants. It translates into smaller phonon frequencies and mean group velocities in fluorides as compared to oxides. Fluorides also have smaller heat capacities, due to their larger lattice parameters (see Supplementary material). Those two factors mainly drive the important discrepancy of the thermal conductivity between fluorides and oxides. Following the same reasoning, it means that halide perovskites in general should have a very low thermal conductivity.
Finally, we analyze the correlations between the thermal conductivity and different simple structural descriptors. Fig. 6 displays the correlograms for fluorides and oxides between the following variables: the thermal conductivity κ, the thermal conductivity in the small grain limit κ sg [3,91], the mean phonon group velocity v g , the heat capacity c V , the root mean square Grüneisen parameter γ rms [92,93], the masses of atoms at sites A and B of the perovskite ABX 3 , their electronegativity, their Pettifor number [94], their ionic radius, the lattice parameter of the compound and its electronic gap. Remarkably, sites A and B play very different roles in fluorides and oxides. In particular, the thermal conductivity of fluorides is mostly influenced by substitutions of the atom inside the fluorine octahedron (site B ), while the interstitial atom at site A has a negligible impact. The opposite is true for the oxides. This means that when searching for new compounds with a low lattice thermal conductivity, substitutions at the A site of fluorides can be performed to optimize cost or other considerations without impacting thermal transport. It is also interesting to note that the gap is largely correlated with the electronegativity of atom B, suggesting the first electronic excitations likely involve electron transfer from the anion to the B atom.
Common to both fluorides and oxides, the lattice parameter is mostly correlated with the ionic radius of atom B rather than atom A. Interestingly, the lattice parameter is larger for fluorides, although the ionic radius of fluorine is smaller than for oxygen. This is presumably due to partially covalent bonding in oxides (see e.g. Ref. 95). In contrast, fluorides are more ionic: the mean degree of ionicity of the X-B bond calculated from Pauling's electronegativities [96] e X and e B as I XB = 100 1 − e (e X −e B )/4 yields a value of 56% for oxides versus 74% for fluorides. Ionicity is also reflected by the band structure, as can be seen from the weak dispersion and hybridization of the F-2p bands [97]. This may explain why the role of atoms at site A and B is so different between the two types of perovskites. We think that the more ionic character combined to the small nominal charge in fluorides makes the octahedron cage enclosing the atom B less rigid, such that the influence of the atom B on the thermal conductivity becomes more significant. Figure 6. Correlograms between the thermal conductivity κ, the thermal conductivity in the small grain limit κsg, the mean phonon group velocity vg, the heat capacity cV, the root mean square Grüneisen parameter γrms, the masses mA and mB of atoms at sites A and B of the perovskite ABX 3, their electronegativity eA, eB , their Pettifor scale χA, χB , their ionic radius rA, rB , the lattice parameter of the compound a latt and its electronic gap, for mechanically stable fluorides (left) and oxides (right) at 1000 K.

V. CONCLUSION
Employing finite-temperature ab-initio calculations of force constants in combination with machine learning techniques, we have assessed the mechanical stability and thermal conductivity of hundreds of oxides and fluorides with cubic perovskite structures at high temperatures. We have shown that the thermal conductivities of fluorides are generally much smaller than those of oxides, and we found new potentially stable perovskite compounds. We have also shown that the thermal conductivity of cubic perovskites generally decreases more slowly than the inverse of temperature. Finally, we provide simple ways of tuning the thermal properties of oxides and fluorides by contrasting the effects of substitutions at the A and B sites. We hope that this work will trigger further interest in halide perovskites for applications that require a low thermal conductivity.
This work is partially supported by the French "Carnot" project SIEVE. C. Oses acknowledges support from the National Science Foundation Graduate Research Fellowship under Grant No. DGF1106401. We also acknowledge the CRAY corporation for computational support.