Complementary evaluation of structure stability of perovskite oxides using bond-valence and density-functional-theory calculations

Abstract Estimation of structure stability is an essential issue in materials design and synthesis. Global instability index (GII) based on bond-valence method is applied as a simple indication, while density functional theory calculation is adopted for accurate evaluation of formation energy. We compare the GII and total energy of typical ABO3-type perovskite oxides and rationalize their relationship, proposing that the criteria for empirically unstable structures (GII > 0.2 valence unit) correspond to the difference in total energy of 50–200 meV per formula unit.


Introduction
Structure stability is an essential factor for evaluating and predicting materials including hypothetical compounds. For ionic compounds, ionic radii [1] and the radius ratios of cation to anion are primary indices to categorize and predict stable coordination polyhedra in crystal. Goldschmidt tolerance factor t [2] is extensively adopted to evaluate structure stability of perovskite structure: t = r A + r X √ 2(r B + r X ) , where r A , r B , and r X are ionic radii of A, B, and X ions in ABX 3 perovskite, respectively. When t is close to the unity, cubic perovskite (space group Pm3m) is stable at ambient conditions. When t is larger or smaller than the unity, distorted perovskite structures, typically GdFeO 3 -type orthorhombic perovskite, or phases other than perovskite (e.g. ilmenite) are expected instead of undistorted cubic perovskite. This index also serves to interpret systematic trends in structures and various properties for ABO 3 oxides [3][4][5][6]. Bond valence method is widely applied to verify adequacy of structure refinement and also to estimate oxidation states of constituent ions [7,8]. Bond valences (s ij ) are calculated by a formula: s ij = exp . l ij is a bond length, l 0 is a bond valence parameter empirically determined using experimental room-temperature structure data, and b is a parameter concerning bond softness [9]. The bond valence sum (BVS) is calculated as the sum of bond valences around the relevant ion: BVS is approximately consistent with oxidation number of ion, thus enabling estimation of valences of constituent elements and cation/anion orderings in complex metal oxides [10][11][12]. BVS-based indices are utilized for evaluation of structural distortions; the bond discrepancy factor (d) is defined as the deviation of BVS from formal valence: where V M is the formal valence of M ion. The d M serves to estimate bond strains. When d M > 0, the bonds around M ion are under compression stress (overbonding), whereas d M < 0, under tensile stress (underbonding). The global instability index (GII) [13] is defined as root mean square of bond discrepancy index in the unit cell: , where N is the number of ions. GII is adopted for comparing structure instability in various perovskite forms including hypothetical phases [14,15]; GII = 0 means perfectly stable structures without steric distortions although structures with GII > 0.2 valence unit (v.u.) are uncommon [8]. The GII-based structure prediction program, SPuDS, developed by Lufaso et al.
can calculate the degree and order of structure instability for various perovskite structures [15,16], as well as the tolerance factor. The close relationship between GII and tolerance factor was demonstrated by Zhang et al. [17]. The indices d and GII are also utilized for interpretation of structure-property relationship. Quadruple perovskite oxides RCu 3 Fe 4 O 12 (R: rare-earth metals) undergo two distinct electronic phase transitions, charge disproportionation (8Fe 3.75+ → 5Fe 3+ + 3Fe 5+ ) and charge transfer (3Cu 2+ + 4Fe 3.75+ → 3Cu 3+ + 4Fe 3+ ) at low temperatures [18,19], the type of which is predominantly determined by strains on R-O and Fe-O bonds [20]. The above-mentioned indices derived from crystal structures are helpful because of their concise procedures, but precise evaluation of structure stability needs first-principle calculations. The most stable structure at the ground can be estimated by comparing total energies of polymorphs and even structural phase transition pressures are also available by using common tangent in total energy versus volume curves [21,22]. The stability of a crystal structure is evaluated by performing phonon calculations with harmonic approximation [23]. The presence of imaginary phonon modes tells the structure is dynamically instable, and their eigenvectors in the reciprocal space gives us possible transition pathway [24]. However, these procedures are rather complicated and require enormous computational time.
To our knowledge, most of structure stability evaluation studies are performed using either BVS-based indices or first-principle calculations. Complementary evaluation and integration of multiple indices improve structure prediction and resulting materials design. The aim of this paper is to evaluate structural stability of typical perovskite oxides by comparing GII and total energy. The empirical criteria (GII < 0.2 v.u. in stable structure) correspond to energy scale of < 50-200 meV per formula unit (f.u.) for perovskite structure using first-principles calculations.
Single-point energy calculations for SrTiO 3 and BaZrO 3 in cubic perovskite structure were performed using cubic unit cell of length sa exp with fixed internal coordinates, where a exp is the experimental lattice constant, and s is the scaling factor changed by 0.001 step around unity. The calculations were performed using the plane-wave-based augmented projector wave (PAW) method as implemented in the VASP code [25][26][27][28][29]. The exchange-correlation interaction was treated as a generalized gradient approximation (GGA) using the function formulated by Perdew, Burke, and Ernzerhof [30]. The PAW potential data-set with radial cutoffs of 1.45, 2.3, 2.8, 2.5, 2.8, 1.9, 2.8, 2.5, 2.9 and 1.52 Bohr for Na, Ca, Ba, Sr, La, Al, Ti, Zr, Ta and O were employed, respectively, where Na-3s, 3p, Ca-3s, 3p, 4s, Ba-5s, 5p, 6s, Sr-4s, 4p, 5s, La-5s, 5p, 6s, 5d, Al-3s, 3p, Ti-3d, 4s, Zr-4s, 4p, 5s, 4d, Ta-6s, 5d, and O-2s, 2p were described as valence electrons. The plane wave cut-off energy was set to 550 eV, and 12 × 12 × 12 k-points mesh with the Monkhorst-Pack scheme was used for integration in the reciprocal space. The minimum value was set to zero in the relative total energy. Figure 1 shows calculated GII for the selected ABO 3 -type oxides. Only two oxides, SrTiO 3 and BaZrO 3 , had the minimum GII values close to zero, whereas others did not reach near zero (< 0.01 v.u.). SrTiO 3 (t = 1.001) and BaZrO 3 (t = 1.004), which crystallize in cubic perovskite at ambient conditions, displayed linear behavior in GII versus a. LaAlO 3 (t = 1.009) exhibited an almost linear curve down to the minimum value (GII = 0.027 v.u.), whereas it crystallizes in the rhombohedral perovskite structure at ambient conditions. This implies that the cubic perovskite phase is nearly stable for this oxide, as expected from the t value close to the unity. NaTaO 3 (t = 0.967, orthorhombic perovskite at ambient conditions), which had a smaller GII value (0.102 v.u.) next to LaAlO 3 , demonstrated slightly nonlinear behavior. Other oxides with larger GII values (>0.1 v.u.) at the minima showed quadratic-like dependences and the deviations from a linear curve were apparently enhanced as the minimum GII increased. These arguments imply that the GII values hardly apply to these oxides because the energy is proportional to square of GII in principle [8]. This result indicates that the optimized crystal structure could be obtained for this oxide, which is in agreement with the literature [31]. Also, BaZrO 3 displayed similar behavior with SrTiO 3 ; GII = 0.003 v.u. at a 0 = 4.175 Å when the structure optimization was achieved. The a 0 values obtained at GII minima were almost identical with the experimental ones for room-temperature cubic perovskite phases (Table 1).

Results and discussion
Since tolerance factors for these two oxides are almost the unity, the GII is simple index for stable structures. When t is smaller than the unity [e.g. t = 0.966 in CaTiO 3 , see Table 1 and Figure 2 . As a result, structure transformations including rotations of the TiO 6 octahedra occur, as evidenced by the fact that CaTiO 3 crystallizes in GdFeO 3 -type orthorhombic perovskite structure at ambient conditions (room temperature and atmosphere pressure). BaTiO 3 , whose t is larger than the unity, also did not converge into stable cubic perovskite structure. Unlike CaTiO 3 , overbonding of Ba-O bonds (d Ba = 0.570 v.u.) and underbonding of Ti-O bonds (d Ti = −0.576 v.u.) were predominant in the GII-optimized structure. This oxide actually crystallizes in tetragonal perovskite structure at ambient conditions. The above results demonstrate that GII and t are consistent, confirming their usefulness as structure instability indices. Figure 3 shows comparison between GII 2 and relative total energy (E total in eV/f.u.) for SrTiO 3 and BaZrO 3 . The optimized lattice constants of SrTiO 3 and BaZrO 3 by DFT calculations with minimum total energy were overestimated from experimental value by about 1.5%. This is mainly ascribed to the GGA used in the calculation to express the exchange-correlation energy. However, the change of total energy induced by the lattice strain can be calculated in better accuracy because significant part of the errors in exchange-correlation energy introduced by GGA is canceled out. In both oxides, GII 2 and E total could be fitted by using the quadratic function fitting: GII 2 = k 1 (a -a 0 ) 2 + G 0 2 and E total = k 2 (a -a 0 ) 2 , where a 0 is the unit cell edge when GII 2 or E total is the minimum, and G 0 is the minimum GII value. Fitting range were limited within the GII < 0.2 v.u. (GII 2 < 0.04 v.u. 2 ). The obtained coefficients k 1 and k 2 , and resulting k 2 /k 1 are listed in Table 2. The k 1 values for SrTiO 3 and BaZrO 3 are almost identical, and k 2 also similar. The k 2 /k 1 values are approximately 1.3 eV/(v.u. 2 f.u.). Note that this is the first example of the conversion coefficient between GII experimentally observed structures (GII 2 exp ) [32][33][34][35][36][37][38][39]. The instability in total energy was estimated by the difference between the energy at the minimum for cubic structure (E cubic ) and the energy at relaxed non-cubic perovskite structure by DFT. Figure 4 demonstrates the correlation between the DFT-and GII-based instability. Rough estimation by linear fitting represents a and E total in cubic perovskite structure, enabling total energy estimation as E total = ~1.3 × GII 2 eV/f.u.
Furthermore, we tentatively evaluated the instability for severely unstable structures. The instability of hypothetical cubic perovskite structure in Figure 1 was estimated from the differences in GII 2 between the minimum values in cubic structure (GII 2 cubic ) and values for  Table 1. Tolerance factor (t), ionic radius of A-and B-ions, GII for cubic (GII cubic ) and distorted perovskite structures (GII exp ), a 0 , and crystal structure at room temperature determined by experiment for selected ABo 3 -type oxides.  [39] determined by the structures at room temperature. On the other hand, total energies obtained calculated by DFT represent the internal energies at 0 K; the temperature dependence of internal energy and entropy term in the Helmholtz free energy are not taken into account. The changes of free energy introduced by lattice strain can be estimated quantitatively by performing phonon calculations under quasi-harmonic approximation if a target compound is dynamically stable [24,40]. Reflecting the thermal expansion of a crystal, the equilibrium volume at finite temperature slightly shifted from that at 0 K. The curvature of energy-volume curve (similar to k 2 in this manuscript) also changes as increasing the temperature. However, the changes of these quantities at room relationship of E total = 0.06(5) + 3.2(5) × GII 2 eV/f.u. The coefficient is about 2.5 times larger than those obtained by the analysis of stable cubic phases in Figure 3, but their scales are at the same level. Large-scale data analysis is needed to improve the total energy-GII relationship, but it is beyond the aim of the present study.
It should be noted that the GIIs were calculated using bond-valence parameters, which were empirically  temperature from those at 0 K are insignificant (order of 0.1%).
We verified adequacy of structure instability obtained from GII and E total calculations. It is known that stable structures with GII > 0.2 v.u. are rare. This GII value corresponds to energy instability by ~50 meV/f.u according to the relationship obtained from cubic perovskites (SrTiO 3 and BaZrO 3 ) in Figure 3, whereas ~200 meV/f.u. was obtained from the analysis of non-cubic perovskites (CaTiO 3 , BaTiO 3 , CaZrO 3 , SrZrO 3 , NaTaO 3 , and LaAlO 3 ) in Figure 4. To our knowledge, there are no established criteria to judge instability in total energy but energy differences among perovskite polymorphs are usually obtained by the order of several ten meV [41][42][43][44]. Therefore, GII calculations can reasonably exclude unstable phases at low cost.

Conclusions
We compared GII on bond-valence model and total energy by DFT calculations for cubic perovskite oxides. The obtained conversion coefficient is in roughly agreement with the energy scale in the previous DFT calculations, supporting the adequacy of GII as structure instability index. As well as ionic radii and Goldschmidt tolerance factor, GII can be adopted for structure-instability evaluation and prediction. For further study, integration of the above indices based on large-scale data analysis is expected to achieve fast and reliable materials design and exploration.

Disclosure statement
No potential conflict of interest was reported by the authors.