Complete phase diagram of rare-earth nickelates from first-principles

The structural, electronic and magnetic properties of AMO3 perovskite oxides, where M is a 3d transition metal, are highly sensitive to the geometry of the bonds between the metal-d and oxygen-p ions (through octahedra rotations and distortions) and to their level of covalence. This is particularly true in rare-earth nickelates RNiO3 that display a metal–insulator transition with complex spin orders tunable by the rare-earth size, and are on the border line between dominantly ionic (lighter elements) and covalent characters (heavier elements). Accordingly, computing their ground state is challenging and a complete theoretical description of their rich phase diagram is still missing. Here, using first-principles simulations, we successfully describe the electronic and magnetic experimental ground state of nickelates. We show that the insulating phase is characterized by a split of the electronic states of the two Ni sites (i.e., resembling low-spin 4+ and high-spin 2+) with a concomitant shift of the oxygen-2p orbitals toward the depleted Ni cations. Therefore, from the point of view of the charge, the two Ni sites appear nearly identical whereas they are in fact distinct. Performing such calculations for several nickelates, we built a theoretical phase diagram that reproduces all their key features, namely a systematic dependence of the metal–insulator transition with the rare-earth size and the crossover between a second to first order transition for R = Pr and Nd. Finally, our results hint at strategies to control the electronic and magnetic phases of perovskite oxides by fine tuning of the level of covalence. A new theoretical approach provides a complete phase diagram of rare-earth nickelates, reproducing the key features seen in experiments. Transition metal oxides with a pervoskite crystal structure exhibit a broad range of behaviours due to a complex the interplay between lattice, electronic and magnetic degrees of freedom. Rare-earth nickelates are a particularly interesting class of perovskite oxide that undergo a highly tunable (and potentially exploitable) metal-insulator transition. Theoretically describing the insulating phase of nickelates, however, is far from trivial and two seemingly distinct descriptions have emerged. A team of researchers led by Manuel Bibes at Unité Mixte de Physique CNRS/Thales use first-principles simulations to somewhat reconcile these conflicting visions, fully describing the electronic and magnetic ground state of nickelates, as well as hinting at strategies for tuning these fascinating materials.


INTRODUCTION
Transition metal oxides with an AMO 3 perovskite structure have attracted widespread interest over the last decades, both from academic and application points of view. This can be ascribed to their wide range of functionalities that originates from the interplay between lattice, electronic, and magnetic degrees of freedom. 1 Among all perovskites, rare-earth nickelates R 3+ Ni 3+ O 3 (R = Lu-La, Y) might be considered as a prototypical case because they posses almost all possible degrees of freedom present in these materials. Nickelates were intensively studied during the nineties 2, 3 and have regained interest in the last few years due to their great potential for engineering novel electronic and magnetic states. [4][5][6][7][8][9][10][11] Except for R = La, all rare-earth nickelates undergo a metal-insulator phase transition (MIT) at a temperature T MI , accompanied by a symmetry lowering from Pbnm to P2 1 /n yielding two different Ni sites. 2,3 In the insulating phase, a subsequent breathing distortion of the NiO 6 cage octahedra is observed, producing a rock-salt-like pattern of large and small NiO 6 groups. In the following, we will use the notation Ni S and Ni L to refer to the Ni cations belonging to the small and large NiO 6 groups, respectively. At T N ≤ T MI , they undergo an antiferromagnetic (AFM) phase transition leading to a quadrupling of the magnetic unit cell ( k ! ¼ ð 1 2 ; 0; 1 2 Þ with respect to the Pbnm primitive 20-atom unit cell) and possible collinear or non collinear spin orderings. [12][13][14][15][16] The electronic structure is also characterized by strong overlaps between O-2p and Ni-3d states inducing large covalent effects. 2 As a consequence, external stimuli, such as temperature, or chemical or hydrostatic pressure, can modify the electronic bandwidth and influence the MIT. [17][18][19][20] Nickelates have thus been used as a platform to search for novel electronic phases through strain engineering or confinement. 10,11,[21][22][23] Nevertheless, the description of the insulating phase of nickelates-and its origin-is still under debate. Originally, the MIT has been ascribed to ionic charge disproportionation effects of Ni cations from 2Ni 3+ to Ni ð3þδÞþ S + Ni ð3ÀδÞþ L , [24][25][26] that result in the rock-salt-like pattern of large and small NiO 6 groups in the low symmetry phase. 25 This presence of the breathing distortion has been argued to be an alternative to Jahn-Teller distortions-Ni 3+ cations are expected to be nominally Jahn-Teller active-, by removing electronic degeneracies between Ni sites and rendering a localized charge-ordered state opening a band gap. 27 In other words, through such a distortion, the e g levels of large NiO 6 groups are stabilized with respect to those of small NiO 6 groups. This scenario was first suggested by several neutron diffraction studies, 16,24,25,28,29 for which a charge ordering of Ni sites was derived from Ni-O bond lengths analysis based on the empirical valence bond model of Brown. 30 X-ray absorption spectroscopy (XAS) measurements 26 on several elements present also characteristics of a Ni valence splitting at the MIT, with a magnitude modulated by the rare-earth. Finally, the strongest evidence is provided by Mössbauer spectroscopy showing a clear-cut split of the Ni charges when entering the insulating phase. 29 Most recently, an alternative covalent vision has emerged and appears to be somehow contradictory with this original scenario. 31,32 The starting point of this theory is the ability to transfer charges from the ligands to the Ni cations due to the strong covalent character of nickelates. It follows that Ni cations adopt a Ni 2+ 3d 8 L electronic configuration in the metallic phasethe notation L defines a ligand hole that is created on the surrounding oxygens-rather than the expected Ni 3+ 3d 7 configuration. In the insulating phase, this configuration transforms into a Ni L 3d 8 and a Ni S 3d 8 L 2 electronic configurations through a transfer of charge from the oxygens from the small NiO 6 groupsthe notation L 2 labels two holes created on the oxygens-, which is compatible with the strong breathing mode observed in this phase. 31 Consequently, this mechanism is rather labeled as a bond disproportionation effect and implies similar valence states for both Ni sites and de facto an absence of charge ordering in the insulating phase. However, these theoretical observations are at odds with the clear experimental evidences of charge disproportionation extracted through several techniques. 26,29 Here we performed density functional theory (DFT) calculations and a Wannier analysis on various representative nickelates in order to reconcile these two apparently conflicting visions of the insulating phase. Our simulations reveal the co-existence of strong ionic effects, with a clear-cut split of the electronic structure of the inequivalent Ni sites resembling a high-spin Ni 2þ L and a low-spin Ni 4þ S , and of covalent features, with O-2p electrons shared with the charge-depleted Ni S cations leaving virtually no charge ordering. Our results therefore propose a unified picture for the description of the insulating phase. By performing calculations on several species, we build a theoretical phase diagram in which all the key features of nickelates can be recognized, and most notably the correct AFM states and the crossover from a second to a first-order metal-insulator phase transition for R = Nd and Pr. Finally, we unveil a new pathway to control electronic and magnetic phases in perovskites by tuning the level of covalence.

Structural properties
First, we performed full geometry relaxations considering 80-atom supercells of both possible Pbnm and P2 1 /n structures with different magnetic orderings: ferromagnetic (FM) as well as complex E-type, S-type, and T-type AFM orderings 5 based on ↑↑↓↓ spin chains in the (ab)-plane with different stackings along the c ! axis (see Fig. 1a-c). The results are presented for spins treated at the collinear level, although we have carefully checked that allowing the spin-orbit interaction and non collinear spin arrangements only yields metastable magnetic phases and does not change the main conclusions of the study (see Supplementary Material). We employed the PBEsol functional 33 in combination with a U correction of 2 eV on Ni-3d states as implemented by Dudarev et al. 34 in order to account for electronic correlations. Several nickelates (R = Y, Dy, Tb, Gd, Eu, Sm, Nd, and Pr) were considered, covering the phase diagram as a function of rare-earth radius. All nickelates relax to a P2 1 /n insulating ground state with complex AFM structures (S-type or T-type depending on the rareearth) and band gaps compatible with experiments 35, 36 (see Table 1). All our Pbnm phases favor a metallic FM solution, while classical AFM orders-A, C and G type-yield solutions that are higher in energy by at least 180 meV per 80 atoms unit cell (see supplementary Material). We checked the reliability of our DFT + U calculations by changing the U correction to either 0 eV or 5 eV in SmNiO 3 . While the ground state is unchanged when no Ucorrection is applied-the gap is eventually decreased to 91 meV-, imposing U = 5 eV yields a P2 1 /n FM and insulating solution that is much more stable than the considered complex AFM orderings (ΔE≃160 meV per 80-atom unit cell). This further supports our choice of a relatively small Hubbard correction for the Ni-3d electrons.
Our optimized ground state structures are characterized by three main lattice distortions (the optimized lattice parameters are provided in the Supplementary Material). First, they feature two antiferrodistortive (AFD) modes that can be described, respectively, as a − a − c 0 and a 0 a 0 c + patterns using Glazer's notation. 37 These AFD modes are the main features of the phase with Pbnm symmetry. Second, we have a breathing of the O 6 octahedra, B oc (see Fig. 2a-c). The breathing mode only appears in the P2 1 /n symmetry and produces the rock-salt-like pattern of small and large NiO 6 groups, automatically resulting in two different Ni sites (see Fig. 2c).
As usual in perovskites, the magnitude of the metal-oxygenmetal bond angles associated with the O 6 rotations is governed by steric effects (see Fig. 2d, e), and nickelates with low tolerance factors (i.e., smaller R cations) 38 are more distorted. The alternating expansion/contraction pattern of the oxygen cage associated with the B oc breathing mode also appears to be modulated by the rare earth (see Fig. 2f), as smaller R cations yield larger distortions. Finally, we observe a Jahn-Teller distortion in the ground state that is one to two orders of magnitude smaller than the breathing mode or the two AFD motions, and it is also much smaller than the Jahn-Teller motion appearing in prototype compounds like LaMnO 3 39 (amplitudes of distortions of our nickelates ground state structures are reported in the Supplementary Material). Hence, the relaxed structures indicate that there is no significant orbital order in these systems, although the 3d 7 -t 6 2g e 1 g electronic configuration of Ni 3+ in the high temperature Pbnm phase is nominally Jahn-Teller active. 27 Disproportionation signatures The electronic structure of the optimized ground states is characterized by strong hybridizations between O-2p and Ni-3d Fig. 1 Sketch of the three complex colinear AFM orderings used in the calculations. a E-type AFM ordering. b S-type AFM ordering. c T-type AFM ordering. The A cations are not displayed for clarity. The dashed lines represent the size of the crystallographic unit cell levels, as inferred from the projected density of states (pDOS, see Fig. 3a for the representative case of SmNiO 3 ). Comparing the pDOS corresponding to the 3d levels of the two different Ni sublattices reveals some small differences, likely reflecting weak disproportionation effects and a small charge ordering (see Fig. 3b). Although atomic charges are not uniquely defined in DFT calculations and there are no unambiguous ways to extract them, 40,41 sphere integrations around the Ni cations can provide some insight into the possible charge ordering, keeping in mind that they might not reflect the real ionization state 40 (We emphasize that the number of electrons extracted from the sphere integrations is not reminiscent of the real site occupancy. We have performed calculations on BaNiO 3 -a non magnetic insulator adopting a P6 3 cm hexagonal structure and in which the Ni cations are formally in a 4+ oxidation state-and La 2 TiNiO 6 -an insulator adopting a double perovskite P2 1 /n phase and in which Ni cations are formally in a 2+ oxidation state. Those calculations reveal that the total Ni site occupancy increases with increasing the oxidation state of Ni cations-Q Ni = 9.631 e in BaNiO 3 , Q Ni = 9.376 e in SmNiO 3 (Pbnm phase) and Q Ni = 9.176 e in La 2 TiNiO 6 . At the same time, the extracted d occupancy on Ni sites remains contants for the three valence states-Q Ni,d = 8.257 e for BaNiO 3 , Q Ni,d = 8.276 e for SmNiO 3 (Pbnm phase) and Q Ni,d = 8.272 e for La 2 TiNiO 6 ). Figure 3c reports the occupancy of both Ni sites as a function of the rare earth. A weak and rather constant level of charge ordering is observed between the Ni sites, going from δ = 0.13 (YNiO 3 ) to δ = 0.11 (PrNiO 3 ), but the sign of δ is opposite to what is expected: the Ni S cations, sitting at the center of the smallest O 6 octahedra, appear to hold more electrons than the Ni L cations, located in the largest oxygen cages. Since the breathing mode B oc enhances the crystal field splitting at the small NiO 6 groups, the e g levels of Ni S lie higher in energy than those of Ni L 27, 42 and therefore Ni S should have fewer electrons than Ni L associated to it.
Let us now consider better defined-and experimentally measurable-quantities, such as Born effective charges (BECs) that measure the amount of charge displaced upon the movement of individual atoms. Figure 3d reports the average of the diagonal components of the tensor for the different nickelates (see the Supplementary Material for the full tensors). In the representative case of SmNiO 3 , we obtain Z NiL % þ2:5, which is not far from the nominal oxidation state of 2+ that this Ni site is associated with in the complete-charge-disproportionation picture. However, we find a similar Z NiS % þ2:1, which sharply deviates from the expectation value (4+) in the chargedisproportionation picture. As shown in Fig. 3d, we observe the same behavior across the various studied compounds; an approximately constant difference of Born charges, of about 0.4 electrons, is found across the whole series. Hence, from the point of view of the effective charges, the two Ni sites behave in a rather similar fashion, the disproportionation effects being weak and not complying with the usual picture.
However, our computed magnetic moments on both Ni sites appear to be in contradiction with the conclusion of the charge analysis. Indeed, as shown in Fig. 3e, we observe a large difference between Ni L -with a moment estimated between 1.25 μ B and 1.15 μ B when going from R = Y to R = Pr-and Ni S -for which the  magnetic moment is null -, which suggests two very different electronic states. Finally, Fig. 3f reports the charge disproportionation extracted from the empirical Brown bond valence model, 30 which is just an alternative way to monitor the breathing mode magnitude. This model suggests a differentiation of the Ni valences in our optimized insulating phases, with a disproportionation δ modulated by the rare-earth, but the result of this phenomenological model is at odds with the charges-statics or dynamicsextracted directly from our simulations. However, this simple quantification yields values that are very similar to experiments. 16,24,25,28,29 It is worth emphasizing that this latter model predicts the expected sign for the charge ordering.

Wannier analysis
The conclusion of the previous discussion is that, in the RNiO 3 compounds, all Ni atoms seem to display a similar oxidation state. Yet, the presence of a significant breathing distortion and of a strong difference in the local magnetic moments clearly suggest two markedly different electronic states. Note that similar results have been reported in previous theoretical works using a variety of methods. 23,31,43,44 Here, to gain insights into the real site occupancies, we ran a maximally-localized Wannier function (MLWF) analysis of our first-principles results, which allowed us to resolve this pending issue.
We used the Wannier90 package [45][46][47] to determine the MLWFs that reproduce the occupied electronic manifold. More precisely, our purpose was to count how many occupied MLWFs are centered at the different Ni cations, and how many at the surrounding oxygen anions, and to characterize them. Further, we wanted to run our analysis without having to make any assumption on the precise character of the occupied Ni and O orbitals, which complicated the choice of the initial atomic orbitals onto which the Kohn-Sham states are projected in order to extract the initial gauge matrix for the localization procedure. Nevertheless, we found the following robust strategy to proceed. We considered the whole occupied manifold and sought to extract from it (i.e., to disentangle) a set of 2 × (144 + 10) MLWFs, where 2 × 144 = 2 × 3 × 48 is the total number of O-2p orbitals available in our 80-atoms supercell and 2 × 10 = 2 × 5 × 2 is the number of Ni-3d orbitals corresponding to two specific Ni atoms. (Note that we ran separate MLWFs optimizations for the spin-up and spin-down channels.) Hence, we used the following initial atomic orbitals: three generic p orbitals centered at each O anion and 5 generic d orbitals centered at two neighboring Ni cations; this couple of Ni cations were chosen to be first-nearest neighbors, so that we considered one Ni L and one Ni S . The basic qualitative results of this optimization were the same for all the nickelates considered, and thus the following discussion is not compound specific.
Our optimization renders 2 × 3 MLWFs centered at each oxygen anion (i.e., three spin-up MLWFs and three-very similar-spindown MLWFs), suggesting that all oxygens in our nickelates are in a 2-oxidation state. The oxygen-centered MLWFs have a clear p character, as can be appreciated in Fig. 4c, d. We also obtained 2 × 3 t 2g -like MLWFs centered at each of the two considered Ni atoms, indicating that the t 2g states are fully occupied and there is no magnetic moment associated to them. Further, we obtained 2 e glike spin-up MLWFs centered at the Ni L site (see Fig. 4a, b), indicating that this cation is in a 2+ oxidation state and has a significant magnetic moment associated to it. Finally, as regards the other initial atomic orbitals centered at the chosen Ni L (two spin-down d orbitals) and Ni S (2 × 2 d orbitals) atoms, they did not lead to any MLWF centered at those sites. Instead, the maximallocalization procedure resulted in MLWFs centered at Ni and R cations in the vicinity of the considered Ni L -Ni S pair. Thus, in particular, it was impossible to localize any e g -like MLWFs at a Ni S site, which strongly indicates that these Ni cations are close to a 4 + oxidation state. We note that a Ni close to a 4+ oxidation state was experimentally observed for some Ni compounds through 61 Ni Mössbauer spectroscopy. 48 These conclusions were ratified by considering larger clusters of Ni sites for the MLWFs optimization, as well as individual Ni's and/or optimizations in which the oxygen bands were not included. Hence, the Wannier analysis yields a picture of strong charge disproportionation between the Ni sites, which is clearly at odds with the quantitatively similar behavior discussed in the section above. To resolve this apparent contradiction, we now inspect in more detail the oxygen-centered MLWFs.  average. Thus, while the Ni S cations appear to be in a 4+ state when we count how many MLWFs are centered at them, they also receive a significant fraction of electrons coming from the surrounding oxygens, with which they are strongly hybridized. Therefore, this is the explanation why quantitative measures of the charge around the Ni S ions render results that are similar to those of the Ni L and suggest valence state much more reduced than the expected 4+. Across the series, we observe that O-2p centers get closer to the Ni S cations for smaller rare earths, increasing the localization of the charge on Ni S . Additionally, the O-2p MLWFs behave similarly irrespective of their spin polarization, so that their shifting does not result in any magnetic moment at the Ni S sites.
We have run additional MLWF optimizations by initially projecting the Kohn-Sham states on all the "missing" e g states of the small Ni S sites and by restricting the localization procedure to the 16 first unoccupied states (2 e g -like MLWFs × 8 Ni S = 16 MLWFs per spin channel). We are able to obtain the unnocupied e g -like MLWFs on Ni S (see Supplementary Material); these MLWFs exhibit a certain O-2p character, which gives further evidence for the hybridization between oxygen and Ni S states.
The Wannier analysis therefore leads to the conclusion that a strong disproportionation occurs in the system, going toward clearly distinct Ni 4þ S (low-spin, non-magnetic) and Ni 2þ L (high-spin, magnetic) sites. Simultaneously, the O-2p MLWFs approach the Ni 4þ S sites, ultimately yielding Ni S and Ni L that are nearly equivalent from the point of view of the charge (static and dynamic) in their vicinity.

DISCUSSION
While the charge-static or dynamic-analysis yields very similar Ni sites occupancies, the Wannier analysis leads to a picture with a strong disproportionation of Ni oxidation states going toward Ni 4þ S and Ni 2þ L sites in the insulating phase. Our results thus appear to be compatible with the ionic disproportionation effects originally proposed to occur at the MIT. They are indeed reminiscent of this picture, with the observation of a subsequent breathing distortion modulated by the rare-earth (see Fig. 2f)-which results in the varying ionization state according to the phenomenological Brown valence bond model (see Fig. 3f)-and a charge disproportionation between Ni sites complying with numerous experimental works, and most notably with XAS and Mössbauer measurements.
They are also compatible with the model of Mizokawa et al. 49 as well as with recent dynamical mean field theory studies and exact Hartree-Fock calculations, 31,32,[50][51][52][53] proposing a ligand-hole structure in rare-earth nickelates. Indeed, our Wannier analysis indicates that we have a 3d 8 electronic configuration for the Ni L cations. More importantly, it is also compatible with the 3d 8 L 2 configuration proposed for the Ni S sites. 31 According to our MLWF analysis, such a situation would correspond to having all six oxygens around Ni S sharing the 2p electrons that occupy orbitals along the Ni-O-Ni bonds. Further, given that our integrated and dynamical Born charges suggest that the Ni S and Ni L sites host a similar number of electrons, our results do in fact point to a situation in which each O 6 cage shares approximately two electrons with the Ni S at its center, exactly as in the ligand-hole picture. More generally, our results are consistent with the selfregulation response upon a change of the oxidation state of transition metals 40, 54 -we note here that this mechanism explains that similar sphere integrated charges around transition metals might be evaluated irrespective of a change of the oxidation state of the elements.
As regards the magnetic moments at the different Ni sites, our results are also clear and compatible with both proposed pictures: Ni L bears a magnetic moment greater than 1 μ B , as it corresponds to having Ni 2+ in a high-spin configuration-we note here that while the magnetic moment associated with Ni L should be around 2 μ B in a purely ionic vision, screening effects from oxygens reduce the magnetic moment beared by Ni L . Then, Ni S has no magnetic moment associated with it, as it would correspond to a nominal Ni 4+ low-spin configuration. The latter result is partly a consequence of the fact that the electrons shared by Ni S and its surrounding oxygens are spin paired, which can be interpreted as a ligand-hole screening.
Our first-principles calculations, based on the DFT framework, shed light on the nature of the ground state of rare-earth nickelates and agree with the different visions proposed to occur at the MIT. The electronic structure is summarized in Fig. 4e and is based on having (i) a strong charge disproportionation between Ni L and Ni S sites accompanied by a breathing mode and (ii) nearly two electrons from surrounding oxygens "shared" with the depleted Ni S site, leaving the impression of two similar Ni 2+ sites in the system. This electronic structure therefore yields a Ni L site with a dominant ionic character and a Ni S site with a dominant covalent character, as it was already suggested by a structural analysis. 55 Although the electronic structure is similar between all nickelates, the rare-earth cations do have an impact on the level of covalence of the system through their induced lattice distortions. An example is their influence on the magnetic moment of the Ni L cations. Note that our computed Ni L magnetic moments are always far from the nominal value of 2 μ B , as a consequence of the hybridization between Ni-3d and O-2p states. This can be seen in Fig. 4c, d where slightly larger overlaps between O-2p orbitals over the Ni L site arise in the spin-up channel (this particular Ni L site has a net spin up). From Fig. 3e, the level of covalence clearly increases with the tolerance factor (a calculation assuming FM order yields similar conclusions, see Supplementary Material). This is further confirmed by the increase of the average BECs on both Ni S and Ni L , following the interpretation of such effects in other perovskite oxides. 56 We also observe increasing O-2p overlaps with Ni L sites going from R = Y to Pr. Therefore, the level of covalence increases with the tolerance factor of the nickelate, i.e. increasing with the rare-earth ionic radius.
This level of covalence seems to correlate with the stability of the magnetic ordering and the insulating phase (see Table 1 and Fig. 5c right scale). On one hand, with increasing covalence, the energy difference between the AFM (S or T) and FM solutions increases strongly. On the other hand, the insulating gap decreases. In order to further corroborate these trends, we performed additional calculations on SmNiO 3 by applying a hydrostatic pressure of ±8% on the ground state volume and relaxed the atomic positions for each magnetic ordering (see methods). Under compression, we find a sizable enhancement of the covalent character as μ NiL decreases to 1.103 μ B -this is accompanied by a strong reduction of the breathing mode with a variation of the NiO 6 octahedra volume of roughly 3.25% and by a reduction of the two a − a − c 0 and a 0 a 0 c + oxygen cage rotations. The stability of the S-type antiferromagnetic ordering with respect to the FM solution is therefore doubled (ΔE = −293 meV per 80-atom cell). Eventually, the band gap is decreased by 0.11 eV with respect to the ground state value. Under expansion, we observe a weakening of the covalent character of the system with μ NiL increasing to 1.308 μ B -this is accompanied by a strengthening of the breathing mode reaching 17% of NiO 6 volume variation and by an increase of the a − a − c 0 and a 0 a 0 c + oxygen cage rotations; simultaneously, the stability of the complex AFM ordering is roughly reduced by a factor of 2 (ΔE = −67 meV per 80-atom cell).
Finally, the rare-earth atom is known to play a key role in the nature of the MIT. Experimentally it is observed that T MI is different from the magnetic-ordering transition temperature T N for all nickelates except for those in which the rare earth is Pr or Nd. Interestingly, our calculations reflect this differentiated behavior.
Complete phase diagram of rare-earth nickelates from firstprinciples J Varignon et al.
As already mentioned, we obtain an insulating solution for the AFM monoclinic ground state of all considered nickelates. Then, when we consider the P2 1 /n structure with a FM spin arrangement, we also obtain an insulating phase for all R cations ranging from Y to Sm; the corresponding band gaps range from 77 to 41 meV and we observe a relatively large energy gain with respect to the orthorhombic phase (see Fig. 5a, b). Thus, the breathing distortion and disproportionation effects seem sufficient to open the band gap in these compounds, irrespective of the spin arrangement. As a consequence, our results indicate that these nickelates can potentially present an insulating, spin-disordered phase, as they indeed do experimentally for R = Lu to Sm. In contrast, for the monoclinic phase of NdNiO 3 and PrNiO 3 , the FM spin configuration is found to be metallic; further, the stability of the low-symmetry structure with respect to the orthorhombic one drastically decreases. The complex AFM ordering therefore appears to be a necessary condition for the MIT to occur in these two compounds.
Considering the case of SmNiO 3 under 8% of compression, we unveil a similar behavior to bulk PrNiO 3 and NdNiO 3 , since the level of covalence increases. These results suggest that it is possible to control the electronic and magnetic structures of these compounds by tuning the level of covalence in the system.
Our first-principles results can be summarized by the computed phase diagram pictured in Fig. 5c. We have extrapolated the metal-insulator transition temperature (left scale of Fig. 5c) from our simulations by following the formula established in reference: 57 where Δ is the charge-transfer energy (assumed to be constant along the series), ϕ is the average tilting angle of all NiO 6 octahedra, d Ni−O is the average Ni-O bond lengths on the two Ni sites and A is a proportionality constant related to the bare bandwidth. By imposing our computed T MI to match those of the end members YNiO 3 and PrNiO 3 , we fit a charge-transfer energy Δ of 0.893 eV, in very good agreement with experiments and Hartree-Fock methods, 23,58 and A = 8.741 eV Å 7/2 , again very close to recent reports. 23 Remarkably, as is clear from Fig. 5c, when we evaluate Eq. 1 for the other nickelates considered here, we obtain a behavior for T MI that resembles very closely the experimental one. 2, 3 On the right scale of Fig. 5c, we report the energy difference between the FM and the complex AFM orderings for R = Y to Sm. This quantity can be expected to correlate with the Néel temperature for these compound and, indeed we do observe that it exhibits a trend very close to experiments. 2,3 Our computed phase diagram further supports our choice of a small Hubbard U correction on Ni-3d sites and demonstrate that our DFT + U calculations capture all the key physical properties of rareearth nickelates. More generally, this constitutes a good example of the efficiency and reliability of DFT to study strongly correlated systems.
In summary, we have used first-principles methods to investigate the ground state electronic structure and phase diagram of rare-earth nickelates. Our DFT simulations using a moderate Hubbard correction on Ni-3d states reproduce all features reported from experiments (insulating character, structural trends, complex AFM structures, disproportionation effects, covalence, apparent absence of charge ordering). In particular, we show that the insulating phase is characterized by a clear-cut split of the electronic states of the two Ni sites, which can be strictly described as resembling low-spin 4+ and high-spin 2+. At the same time, our simulations reveal a shift of the oxygen-2p orbitals toward the depleted Ni S cations, so that, ultimately, from the point of view of the static and dynamical charge, the two Ni sites appear to be nearly identical. These findings are clearly related to the various pictures proposed in the literature to explain the ground state of these compounds, and they further highlight the efficiency of DFT to study the electronic properties of complex oxides. Performing a systematic study of various materials, we have been able to build a theoretical phase diagram of the rareearth nickelates, underlining the subtle interplay between the level of covalence and the nature of the metal-insulator phase transition as well as the stability of the complex AFM ordering. These latter results hint at strategies to control electronic and magnetic phases of late transition-metal oxide perovskites.

METHODS
First-principles calculations were performed with the Vienna ab-initio simulation package package. 59,60 The results are presented for spin arrangements treated only at the collinear level, although the role of spin-orbit interaction and non collinear spin arrangements have been carefully checked (see Supplementary Material). We used a 3 × 6 × 2 Γcentered K-point mesh. The cut-off was set to 500 eV. We used projector augmented wave pseudopotentials 61 with the following electron configurations: 4 s 2 3d 8 (Ni), 2 s 2 2p 4 (O), 4 s 2 4p 6 5 s 2 4d 1 (Y), 4 s 2 4p 6 5 s 2 (Pr, Nd, Sm, the extra 4f electrons are frozen in the pseudopotentials), 4p 6 5 s 2 (Gd, Eu, Tb, Dy, the extra 4f electrons are frozen in the pseudopotentials). We did not treat explicitly the f electrons as they order at very low temperature 62 and they were included in the pseudopotential. Full geometry relaxations were treated until forces were lower than 0.001 eV/Å and energy was converged to 1 × 10 −7 eV. BECs were computed using density functional perturbation theory. 63 Symmetry adapted modes allowing the extraction of lattice distortion amplitudes were performed using the Bilbao crystallographic server. 64,65 The results presented in Fig. 2d-f have been obtained by freezing single lattice distortion in a cubic reference that has the same lattice parameters for all reported compounds. This allows to extract quantities that are strain independent for all nickelates. Wannier functions reported in Fig. 4a, b are plotted for isosurfaces equal to 2. Wannier functions reported in Fig. 4c, d are plotted for isosurfaces equal to 9. Parameters for the Brown valence bond model analysis can be found in reference. 16 The calculations of SmNiO 3 under hydrostatic pressure have been carried by considering the ground state Bluefilled circles and red-filled squares correspond to second order (T MI ≠ T N ) and first order (T MI = T N ) metal-insulator phase transitions, respectively. Right scale: energy difference in meV per formula unit between the monoclinic FM and AFM solutions for R = Y to Sm unit cell lattice parameters with a homogeneous contraction or expansion of them (±2%). We then relaxed atomic positions by keeping lattice parameters fixed while considering the different magnetic orderings.