{'fraction_non_alphanumeric': 0.042885102463734746, 'fraction_numerical': 0.024752475247524754, 'mean_word_length': 4.509990485252141, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Theory of the magnetic and metal-insulator transitions in RNiO3 bulk and layered structures

A slave rotor-Hartree Fock formalism is presented for studying the properties of the p-d model describing perovskite transition metal oxides, and a flexible and efficient numerical formalism is developed for its solution. The methodology is shown to yield, within an unified formulation, the significant aspects of the rare earth nickelate phase diagram, including the paramagnetic metal state observed for the LaNiO3 and the correct ground-state magnetic order of insulating compounds. It is then used to elucidate ground state changes occurring as morphology is varied from bulk to strained and un-strained thin-film form. For ultrathin films, epitaxial strain and charge-transfer to the apical out-of-plane oxygen sites are shown to have significant impact on the phase diagram.

Understanding the unusual electronic behavior of transition-metal oxides has been a long-standing question in condensed matter physics, 1 and interest has intensified following the demonstration 2 that the materials could be used as components of atomically precise oxide heterostructures. [3][4][5] The theoretical challenge posed by the materials is to treat simultaneously the strong local correlations in the transition metal d-orbitals and their substantial hybridization with oxygen p orbitals. In some systems, the p orbitals can be integrated out and the physics represented in terms of a (possibly multi band) Hubbard model representing the d orbitals only, for which many theoretical methods are available 1 . However, in many cases the charge transfer between p and d orbitals is large enough that the p states cannot be neglected. This "negative charge transfer energy" regime 6 has been less extensively studied. While much useful information has been provided by density functional theory (DFT) 7 and extensions such as DFT+U [8][9][10][11][12][13][14] , hybrid functionals 15,16 , and DFT+DMFT 17,18 , these methods are computationally intensive, so that the large supercells required for long-period ordered phases are difficult to study. Furthermore, the variety of experimental bulk and superlattice configurations and of many-body phenomena emphasizes the need for a model-system treatment that encapsulates the essential physics so the importance of different contributions can be disentangled.
The rare earth nickelates 19 , ReN iO 3 , are an important case in point. According to standard valence counting arguments, Ni is in the 3+ oxidation state with a d 7 valence configuration. However, photoemission experiments 21 reported that the dominant GS configuration is d 8 L, in agreement with unrestricted Hartree-Fock (UHF) 20,21 and DFT+U [12][13][14] , and DFT+DMFT 17 approaches. Additionally, O 1s x-ray absorption experiments have found significant hole concentration on the oxygen p orbitals [21][22][23] . This large degree of charge transfer places the material near the negative charge-transfer regime 6,20 . As Re is varied across the lanathanum row of the periodic table, the ground state (GS) of bulk materials changes from paramagnetic metal (PM-M) to cor- Nd .. Ho related insulator 1,19,24 . The correlated insulator phases exhibit a rock-salt-pattern lattice distortion 25-29 and a nontrivial long-period magnetic ordering [28][29][30][31][32][33] . In ultrathin films a metal-insulator transition which is apparently unaccompanied by rock-salt-pattern lattice distortion occurs as film thickness and strain are varied [34][35][36][37][38][39][40] .
In this paper, we develop a mean-field (MF) approach based on a combination of slave-rotor (SR) [44][45][46][47][48] and Hartree-Fock (HF) methods, along with an efficient and flexible numerical formalism for solving the MF equations. The method is powerful enough to permit the examination of the large supercells needed to investigate long-period ordering patterns in the context of realistic crystal structures. The results reconcile the bulk and film phase diagrams and allows us to identify the key role played by the oxygen degrees of freedom and the lattice distortions. Figure 1 summarizes our key new results. It presents phase diagrams in the space of Ni charging energy U and amplitude δd pd of rocksalt-type lattice distortion (defined more precisely below). The upper panel shows that, in the absence of oxygen breathing distortion, bulk materials are PM-M at any value of U , consistent with experiment 19 . As the lattice is distorted, transitions occur, first to a metallic magnetic phase and then to a magnetic insulator. The nontrivial ordering wavevector found in experiment [28][29][30][31][32][33] is correctly obtained as a 3D 16-formula-unit magnetic pattern of the type ⇑0⇓0, denoting the z-projection of magnetic moment detailed below. The 4 ⇑ and 4 ⇓ Ni sites are surrounded by 6 spin-0 Ni sites, which in our interpretation form singlets with its 6 oxygen neighbors. Each of the 8 spin-0 Ni sites are connected to 3 ⇑ and 3 ⇓ Ni sites via oxygen bonds. States with the ↑↑↓↓ pattern are unstable against PM or ⇑ 0 ⇓ 0 solutions, and indeed are not observed in experiments 30 . In agreement with symmetry-based arguments 43 , slight charge ordering always accompanies the spin ordering. For reasonable parameter values the change in electronic structure across the Re series is properly accounted for, with one exception: the magnetic metallic phase found in theory is not observed in experiment. This is discussed in more detail below. The lower panels show the evolution of the phase diagram with film thickness M and (bottom panel) with applied strain, revealing that in ultra-thin films an insulating phase can occur even in the absence of rocksalt-type lattice distortions.
Our theoretical approach begins from a p-d lattice model of the form H = k H pd k + l (U l + J l ), with k and l sum over crystal momenta and Ni sites, respectively. The parameters of this model may be obtained from, e.g. maximially localized Wannier fits to DFT results 49 , but the precise form is not important for this paper (see also Ref. 50 for another example of the insensitivity of results to the precise p-d model band parameters). We adopt a simplified scenario in which we retain the Ni e g and O 2pσ orbitals, with nearest-neighbor p-d and p-p hopping. The difference between the bare Ni-d and O-2p energies ε d -ε p is important, as discussed below.
The structure of the ReNiO 3 materials is derived from the ideal ABO 3 cubic perovskite structure which is a lattice of corner-sharing oxygen octahedra, each containing a Ni site at the center. The crystal structure of the actual materials is distorted from this structure by rotations of the octahedra which are not important for our purposes and, in the insulating cases, by a twosublattice distortion in which adjacent Ni's have signif-icant different Ni-O bond lengths [26][27][28][29] . To incorporate the bond disproportionation in the Hamiltonian under the cubic approximation, we scale the hopping according to the Harrison rule 51 , t pd = t 0 pd . An additional effect may occur in layered structures. Liu et al. 34 showed that the presence of Al at the interface would deplete holes on the out-of-layer oxygen sites linking Al and Ni, raising the charge-transfer energy from those apical sites by ∼1eV. We model the Mlayer 2D structures using supercells with M NiO 3 units in the z-direction, terminated on both ends with apical oxygen sites whose e p are shifted by -1eV. For reference we also performed some calculations without such shift.
The interaction terms U l and J l contain the on-site repulsion and the Hund's interactions in the rotationallyinvariant Slater-Kanamori form 1 .
In this expression, J l differentiates only states with the same occupancy but inequivalent configurations.
To treat the on-site Coulomb interactions, we adopt the slave rotor (SR) approach [44][45][46][47][48] . The method was originally applied to the Hubbard model and then to other d-only models with Hund's-like interactions 47,48 under the j << U approximation. We apply it here to the p-d model, noting that for RNiO 3 j∼1eV is much smaller than either the d-d repulsion and the electron bandwidth. For each e g site, the approach introduces an auxiliary SR field, θ l ∈ [0, 2π), and decomposes electron operators as The consistency of SR state and doccupancy is enforced by the constrain The spectrum of U (θ) l ∝ L 2 l . We treat the large U l with θ and the small j-scale J l with f † using the weak coupling Hartree-Fock approximation. That is, we solve the p-d model with the single-Slater-determinant ansatz |M F = |p, f |θ . We follow previous SR applications using Lagrange multipliers, h l . Along with constraint Eq. 3 and up to a constant, the system of equations reads where V 0 ǫα is given by the product of t 0 pd and the Slater-Koster orbital symmetry factors in Table I of Ref. 52.
We solve the MF equations without any restriction other than the single-Slater-determinant assumption. A 3D Bravais lattice of Ni 16 O 48 supercells is required for unrestricted modeling of the (1/2, 0, 1/2) pattern with respect to the orthorhombic unit cell. To systematically capture the effects of dimensionality, we model the M -layer structures using a 2D Bravais lattice of Ni 4M O 12M+4 supercells, which is connected to the 3D Ni 16 O 48 supercell for M → ∞. The system is solved by the standard T=0 iterative procedure for up to 65536 supercells. Without much optimization, the worst case initial condition with bi-partite charge-, orbital-, and magnetic-order would converge into a PM-M or (1/2, 0, 1/2)-ordered insulator within 30 cpu hours on an Opteron-2350 cluster.
We perform bulk calculations for (u,e d ) such that n eg ∼ 2 as found in UHF with parameters fitted to photoemission spectrum 20 , DFT 12-14 , and DMFT+DFT 17 . For the reference set of U = 5.3, j = 1, t pσ,d x 2 −y 2 = 1.5, t pp = 0.5, ∆ pd = e d − e p + U = 0.3, the undistorted 3D lattice has n eg = 1.95, slightly less than those of the above calculations. To isolate dimensionality effects, we use the bulk's (u,e d ) pairs, in addition to the aforementioned apical e p shift, for layered calculations.
The left panel of Fig. 2 shows representative results for the energies of different locally stable phases computed at interaction (U,j)=(5.3,1) as a function of (δd pd ) 2 . We see that in agreement with experiment, the undistorted lattice is a PM-M. As oxygen breathing disproportionation δd pd is increased a transition to a period-4 metallic magnetic state occurs. At yet larger distortion a metalinsulator transition occurs. The ferromagnetic states found in UHF 20,21 , DFT+U 12-14 and DFT+DMFT 17 become stable only at larger δd pd and gain less energy than does the experimentally observed ordering pattern. Also, E∼−δd 2 pd implies an equilibrium distortion ultimately determined by the aharmonic lattice restoration force beyond the scope of this study.
The dashed and dotted lines in the right panel of on the two sublattices of the distorted structure. We see that the Coulomb renormalizations are not large, and are only weakly dependent on sublattice. This renormalization has a moderate effect on physical properties, e.g. Fermi velocity renormalization at (U,j)=(5.3,1) is v 5.3,1 /v 0,0 ∼0.65. The relatively modest SR renormalizations confirms that near the "negative charge transfer" region of the phase diagram 6 the on-site Coulomb correlations do not drive the phase transition in an orthodox way, in particular the physics is far from the e iθ → 0 Brinkman-Rice insulator. Even then, the renormalized kinetic terms, combined with lattice distortion and the magnetic order opens a gap in the density of states. We stress that SR renormalization e iθ due to U is required for insulation within physical parameter range. Parameters scaled for different Re bulk materials 26-29 are marked on the phase diagram in top panel of Fig. 1. It captures all experimentally observed phases over the Re series, with the exception of Re=Pr, which is predicted to be a metal instead of an insulator, albeit with the correct charge and magnetic order. This can be understood by noting that our approach underestimates the insulating gap, which is 160meV at δd pd d 0 pd =2.5% compared to ∼200meV in a 2-site DFT+DMFT study 17 , and that Re=Pr is also extrapolated to be a metal from those Re=Lu results.
In agreement with experiments, the ⇑0⇓0 insulator (top panel of Fig. 1) has a (1/2,0,1/2) magnetic structure with respect to the orthorhombic unit cell, a rock-salt charge ordering pattern, but without orbital ordering. Let us identify the Ni sites with longer/shorter Ni-O bond length as Ni long and Ni short and consider the from the corresponding δt pd /t 0 pd ∼10% hopping modulation. This shows that Ni short is strongly hybridized with the oxygen bands while Ni long retains its d 8 characteristic. The magnetic moments are m long =1.1, which agrees with experiment 28,29 , and m short =0. For charge-transfer insulators, the effective physics of strong p-d hybridization is a strong antiferromagnetic S d ·S p spin correlation. In this context, having m short = 0 in the single-Slaterdeterminant approximation is evident of a d 8 -high-spin state forming a singlet with its oxygen neighbors. The formation of p-d singlets has also been reported in a twosite DFT+DMFT study which predicted FM GS 17 .
We now discuss the effects of dimensional confinement. The bulk and M ≥ 3 layers structures have similar phase diagrams, with increasing ordered and insulating region for decreasing M and increasing U . In particular, the lower panel of Fig. 3 shows that PM-M exists for M ≥ 3. The phase diagram of 2-, and 1-layer ultra-thin films are shown in the lower panels of Fig 1. We see that, for M ≤ 2 layers, the ground state is always magnetically ordered, even in the absence of the bond disproportionation. Further, the distortion amplitude required for the transition to the insulating phase decreases with increasing u and decreasing M . At (u,j)=(5.3,1), the metalinsulator boundary ( δd pd d 0 pd ) MIT is reduced from 1.9% to 0.96% and 0.44% for 2-and 1-layer, respectively. A comparison between the solid and dotted lines in Fig. 1 shows that charge transfer to the out-of-plane orbitals 34 greatly enlarge the insulating regions.
To study the effects of epitaxial strain, we focus on 1-layer structure for concreteness. We include the strain effects by scaling hopping integrals according to the Harrison rule 51 for the geometry specified in Ref. 39, which introduced compressive (tensile) strain with LaSrAlO 4 (SrTiO 3 ) such that the lattice parameters are a=b=3.769 (3.853)Å and c = 3.853 (3.790)Å . Breathing-mode distortion is then introduced as discussed above. In agreement with experiments on ultra-thin layers 35,36,39 , Fig. 1 shows that compressive (tensile) epitaxial strain enlarges (shrinks) the insulating region. Under realistic values of U and compressive stain, ordered insulator is possible for δd pd =0, agreeing with experiments which found breathing-mode distortion only under tensile strain 40 . Note that further shift in apical oxygen e p also decreases δd pd required for insulation. Figure 3 demonstrates an orbital polarization of ∼5%. We also found that the polarization decreases (increases) with compressive (tensile) strain, in agreement with a DFT+U study 11,12 . For example, LaSrAlO 4 (SrTiO 3 ) substrate changes the polarization to ∼4 (6)%. While DFT+U predicted FM order, our M =1,2 GS have an AFM in-plane order of ⇑↑⇓↓, which differs slightly from the bulk's ⇑0⇓0 pattern by having small m short =0 which decreases with increasing distortion. The 2-layer structure has an additional magnetic transition between FM and AFM ordering in the z-direction.
Lastly, we compare our results with other theoretical works. In our formulation, the bulk PM-M develops into ⇑0⇓0-I for δd pd d 0 pd >1.9% compared to UHF's δd pd d 0 pd >7.5% transition from ferromagnetic to MCO state with similar pattern 20 . The pattern is also similar to the j/U > 1 S-SDW state in the d-only model 42 , but we have j/U < 1/5 and a bond and charge disproportionation with d 8+δ Ld 8−δ L involving strong pd hybridization. Unlike the d-only model, our results are sensitive to dimensional confinement. Even though different DFT implementations have been employed to capture the bulk LaNiO 3 PM-M 15 and LaNiO 3 /LaAlO 3 ordered insulating layer 16 , the demonstration of all bulk and layered phases by a single formulation had been elusive to the best of our knowledge.
Our formulation bridges the gap between HF-like approaches and the expensive cluster DMFT. The results connect all ReNiO 3 delicate phases in bulk and layer form and also provide detailed insights about them. This demonstrates a viable pathway of treating systems with different superstructures as well as compounds such as Fe and Co oxides with strong p-d hybridization and large number of partially filled strongly correlated orbitals. For example, Sr 2 FeO 4 has also been shown to exhibit strong hybridization and non-trivial magnetic order 53 .
We thank G. A. Sawatzky, H. Chen, H. T. Dang, R. Fernandez, S. Park, and D. Zgid for helpful discussions. This effort is supported by US National Science Foundation under grant NSF-DMR-1006282.