DFT+DMFT study of spin-charge-lattice coupling in covalent LaCoO$_3$

TL;DR

This study uses DFT+DMFT to explore how structural changes influence the spin, charge, and lattice interactions in LaCoO$_3$, revealing the importance of covalency and charge ordering in its spin-state transition.

Contribution

It demonstrates the effectiveness of DFT+DMFT in capturing the multiconfigurational and fluctuating nature of spin and charge states during the transition, surpassing DFT+U results.

Findings

Homogeneous spin state shows multiconfigurational nature with fluctuating HS and LS states.

Charge ordering arises from spin-state ordering in the bond-disproportionated structure.

Structural changes significantly influence the energetics and electronic structure of LaCoO$_3$.

Abstract

We study energetics and the nature of both homogeneous and mixed spin (MS) states in LaCoO$_3$ incorporating structural changes of the crystal volume expansion and the Co-O bond disproportionation (BD) during the spin-state transition using the density functional theory plus dynamical mean field theory (DFT+DMFT) method. DFT+DMFT predicts that energetics of both excited spin states are almost the same while DFT+U calculations of the same structures energetically favor the MS states and produce various metastable solutions whose energetics depend sensitively on final spin states. Within DFT+DMFT, the homogeneous spin state in the expanded crystal volume shows the multiconfigurational nature with non-negligible occupancy probabilities of both high spin (HS) and low spin (LS) states along with $d^6$ and $d^7$ charge configurations indicating the dynamically fluctuating nature of spin and charge states due to the Co-O covalency. The nature of the MS state under the BD structure reveals that Co sites with the long Co-O bonds develop a Mott insulating state and favor HS with a $d^6$ configuration, while more covalent Co sites with the short Co-O bonds occupy more LS states with a $d^7$ configuration and behave as a band insulator, as a result, charge ordering is induced in the BD structure from the spin-state ordering. We also find that both energetics and electronic structure sensitively depend on the Co-O covalency effect, which can be tuned by changing the double counting potential and the resulting $d-$occupancy ($N_d$), and $N_d$ close to 6.7 is consistent with the nature of the spin-state transition. Our results show that structural changes during the spin-state transition can play an important role in understanding energetics and electronic structure of LaCoO$_3$.

Tables

Table 1. Table 1: The structural information of LaCoO 3 relaxed structures obtained using DFT+U with U=6eV and J=0.9eV. Three different structures are obtained by relaxing with different magnetic moments and they are denoted as S1 (relaxed with LS), S2 (relaxed with IS), and S3 (relaxed with MS). The space group (S.G.), the crystal volume per formula unit, the Co-O bond length, magnetic moments, and the d − limit-from 𝑑 d- occupancy ( N d subscript 𝑁 𝑑 N_{d} ) are given in this table.

Structures	S.G.	Vol.[Å3]	Co-O[Å]	Mom.[${\mu }_B$]	$N_d$
S1(LS)	$R\overline{3}c$	56.40	1.95	0	7.1
S2(IS)	$R\overline{3}c$	57.98	1.97	2.2	7.1
S3(MS) Co1	$R\overline{3}$	57.98	1.99	3.0	6.9
S3(MS) Co2	$R\overline{3}$	57.98	1.94	0.3	7.2

Equations

$$E$$

$$E^{DC}=\frac{U}{2}\cdot \overline{N_d}\cdot (\overline{N_d}-1)-\frac{J}{4}\cdot \overline{N_d}\cdot (\overline{N_d}-2)$$

$$V^{DC}=U\cdot (\overline{N_d}-\frac{1}{2})-\frac{J}{2}\cdot (\overline{N_d}-1)$$

$${\hat{H}}_{S1}^{dd}=\left(\begin{array}{ccccc}6.912& 0.006& -0.05& -0.00& -0.01\\ 0.006& 6.236& 0.019& -0.01& -0.02\\ -0.05& 0.019& 6.241& -0.02& -0.02\\ -0.00& -0.01& -0.02& 6.917& -0.01\\ -0.01& -0.02& -0.02& -0.01& 6.237\end{array}\right)$$

$${\hat{H}}_{S2}^{dd}=\left(\begin{array}{ccccc}6.575& -0.01& -0.050& -0.00& 0.009\\ -0.01& 5.948& -0.02& 0.013& -0.02\\ -0.050& -0.02& 5.954& -0.02& 0.020\\ -0.00& 0.013& -0.02& 6.579& 0.010\\ 0.009& -0.02& 0.020& 0.010& 5.949\end{array}\right)$$

$${\hat{H}}_{S3,LB}^{dd}=\left(\begin{array}{ccccc}6.484& 0.003& 0.054& -0.00& 0.005\\ 0.003& 5.913& -0.02& -0.00& 0.018\\ 0.054& -0.02& 5.921& 0.027& -0.02\\ -0.00& -0.00& 0.027& 6.489& 0.002\\ 0.005& 0.018& -0.02& 0.002& 5.915\end{array}\right)$$

$${\hat{H}}_{S3,SB}^{dd}=\left(\begin{array}{ccccc}6.583& 0.006& 0.049& -0.00& 0.014\\ 0.006& 5.909& -0.02& -0.02& 0.022\\ 0.049& -0.02& 5.914& 0.013& -0.02\\ -0.00& -0.02& 0.013& 6.587& 0.012\\ 0.014& 0.022& -0.02& 0.012& 5.910\end{array}\right)$$

$${\hat{H}}_{S1}^{dd}=\left(\begin{array}{ccccc}4.332& -0.01& -0.01& -0.00& 0.014\\ -0.01& 3.737& -0.01& 0.022& 0.006\\ -0.01& -0.01& 3.737& -0.02& 0.004\\ -0.00& 0.022& -0.02& 4.333& 0.022\\ 0.014& 0.006& 0.004& 0.022& 3.735\end{array}\right)$$

$${\hat{H}}_{S2}^{dd}=\left(\begin{array}{ccccc}3.767& 0.013& -0.01& -0.00& -0.02\\ 0.013& 3.158& 0.006& -0.02& 0.007\\ -0.01& 0.006& 3.157& -0.02& -0.01\\ -0.00& -0.02& -0.02& 3.767& -0.02\\ -0.02& 0.007& -0.01& -0.02& 3.156\end{array}\right)$$

$${\hat{H}}_{S3,LB}^{dd}=\left(\begin{array}{ccccc}3.150& -0.02& 0.016& -0.00& -0.02\\ -0.02& 2.490& 0.009& 0.037& -0.01\\ 0.02& 0.009& 2.488& 0.033& 0.006\\ -0.00& 0.037& 0.033& 3.148& -0.04\\ -0.02& -0.01& 0.006& -0.04& 2.487\end{array}\right)$$

$${\hat{H}}_{S3,SB}^{dd}=\left(\begin{array}{ccccc}4.519& -0.01& 0.006& -0.00& -0.01\\ -0.01& 3.984& 0.006& 0.017& -0.01\\ 0.006& 0.006& 3.984& 0.014& 0.005\\ -0.00& 0.017& 0.014& 4.520& -0.02\\ -0.01& -0.01& 0.005& -0.02& 3.983\end{array}\right)$$

$$E^{DC}=\frac{U}{2}\cdot N_d\cdot (N_d-1)-\frac{J}{4}\cdot N_d\cdot (N_d-2)$$

$$V^{DC}=\frac{\partial E^{DC}}{\partial N_d}=U\cdot (N_d-\frac{1}{2})-\frac{J}{2}\cdot (N_d-1)$$

$$V^{DC}=U\cdot (N_d^0-\frac{1}{2})-\frac{J}{2}\cdot (N_d^0-1)$$

$$E^{DC}=\frac{U}{2}\cdot \overline{N_d}\cdot (\overline{N_d}-1)-\frac{J}{4}\cdot \overline{N_d}\cdot (\overline{N_d}-2)$$

$$V^{DC}=U\cdot (\overline{N_d}-\frac{1}{2})-\frac{J}{2}\cdot (\overline{N_d}-1)$$

$$E^{DC}=\frac{U^{\prime }}{2}\cdot N_d\cdot (N_d-1)-\frac{J}{4}\cdot N_d\cdot (N_d-2)$$

$$V^{DC}=U^{\prime }\cdot (N_d-\frac{1}{2})-\frac{J}{2}\cdot (N_d-1)$$

Full text

I Introduction

Transition metal oxides exhibit complex and rich phase diagrams arising from the strongly correlated nature of spin, charge, orbital, and lattice degrees of freedom Imada et al. (1998). LaCoO3 has been known for the spin-state transition of partially filled $d$ orbitals in a Co ion. At very low temperatures, LaCoO3 is a non-magnetic insulator with the low spin (LS) state. As the temperature is elevated above 90K, the magnetic susceptibility changes to a Curie-Weiss form indicating that paramagnetism dominates with higher spin states while retaining an insulating behavior Ivanova et al. (2009). The spin-state transition can be explained based on the atomic multiplet structure of the Co $d$ orbital, namely from the $|S_{z}|$=0 LS ground state to $|S_{z}|$=1 intermediate spin (IS) or $|S_{z}|$=2 high spin (HS) state. Various experimental results have been used to interpret the higher spin state as either IS Saitoh et al. (1997); Abbate et al. (1993); Klie et al. (2007); Zobel et al. (2002) or HS Haverkort et al. (2006); Podlesnyak et al. (2006). The MS of LS and HS has been also suggested to explain other experimental measurements Raccah and Goodenough (1967); Bari and Sivardière (1972); Señarı́s-Rodrı́guez and Goodenough (1995); Kyômen et al. (2005); Yaroslavtsev et al. (2016). Despite extensive experimental works, the nature of excited spin-states has not been clarified yet.

The spin-state transition in LaCoO3 occurs since the Hund’s coupling tends to maximize the spin and excite electrons from $t_{2g}$ to $e_{g}$ orbitals by overcoming the crystal field splitting between them. Here, the interplay between electron and lattice degrees of freedom plays an important role as the $e_{g}$ orbital occupation increases the Co-O bond-length to reduce the repulsive Coulomb interaction of electrons between Co and O ions. This structural change also gives the positive feedback since the reduction of the crystal-field splitting can favor the spin-state transition. This strong electron-lattice coupling has been measured experimentally by the anomalous lattice expansion arising due to the Co-O bond-length elongation at the spin-state transition Asai et al. (1994); Berggold et al. (2008); Knízek, K. et al. (2005); Radaelli and Cheong (2002). Several scenarios of local structural distortions due to the spin excitation have been proposed although clear experimental evidences have not been given yet. The Co-O bond-disproportionation (BD) with alternating the long bond (LB) site and the short bond (SB) site Raccah and Goodenough (1967); Yan et al. (2004) was suggested to accommodate MS with HS and LS. The Jahn-Teller distortion Zobel et al. (2002); Gnezdilov et al. (2006) was also discussed possibly due to the IS state. The strong electron-lattice coupling has been also shown in the tensile-strained LaCoO3 film promoting various competing orders including spin Fuchs et al. (2007); Freeland et al. (2008); Choi et al. (2012), charge Sterbinsky et al. (2018), and orbital orderings Fujioka et al. (2013).

Alongside experimental measurements, various theoretical scenarios based on first-principle calculations have been proposed to address this long-standing problem of the spin-state transition in LaCoO3. Density functional theory (DFT)+U calculations have been predicting that excited spin-states including IS Korotin et al. (1996); Knížek et al. (2005) and MS Knížek et al. (2006, 2009) states can be energetically stable. Dynamical mean field theory (DMFT) can capture the multi-configurational nature of a paramagnetic state fluctuating dynamically beyond DFT+U Georges et al. (1996); Kotliar et al. (2006). Early DFT+DMFT studies computed the single-particle spectra with spin-state crossover Craco and Müller-Hartmann (2008) comparable to the experimental X-ray absorption spectra Abbate et al. (1993) and also studied effects of the pressure Zhang et al. (2012) and the Co-O covalency Křápek et al. (2012) on the spin-state transition. Both the homogeneous spin excitation including the electronic entropy Chakrabarti et al. (2017) and the mixed LS and HS solution without any structural distortions Kuneš and Křápek (2011); Karolak et al. (2015) have been discussed as the possible origin of the spin-state transition within DFT+DMFT. However, no DFT+DMFT studies have addressed yet both the nature and energetics of paramagnetic states with all possible structural changes in this material.

In this paper, we adopt the DFT+DMFT method to study both homogeneous and MS states by incorporating structural changes during the spin-state transition. We show that the DFT+DMFT energetics treating the multi-configurational nature of a paramagnetic state can be noticeably different from the static DFT+U solutions in which various meta-stable spin-states are possible. We also find that both structural changes and the Co-O covalency effect tuned by the double counting potential can strongly affect the electronic structure and energetics of spin-states in LaCoO3. In the expanded crystal volume, the occupation probabilities of higher spin-states in the Co ion increase as Co $d$ orbitals become more correlated. As temperature increases, DFT+DMFT can reproduce the insulator-to-metal transition, consistent with experiment, as the correlation in the Co ion becomes weaker. In the BD structure, the Co ion in the Co-O LB behaves as a Mott insulator with HS while the Co ion in the SB becomes a band insulator favoring LS. The charge ordering is also induced in the BD structure as the LB Co ion favors a $d^{6}$ configuration while the SB Co ion occupies more $d^{7}$ states. The charge-self-consistency in DFT+DMFT plays a role to reduce the charge ordering in the BD structure.

This paper is organized as follows. In Sec. II, the method used in this paper is explained in details. Sec. III.1 discusses the energetics of possible spin-state transition in LaCoO3 comparing DFT+DMFT and DFT+U. The nature of the paramagnetic state treated in DFT+DMFT is shown in Sec. III.2 by computing the occupation probabilities of different spin-states. The origin of charge ordering driven by spin-state ordering in the BD structure is revealed in Sec. III.3. Results of the density of states computed using DFT and DFT+DMFT for different structures and temperatures are shown in Sec. III.4. The self-energy data are displayed in Sec. III.5 and the effect of the charge-self-consistency in DFT+DMFT is also discussed in Sec. III.6. We summarize our paper with conclusions in Sec. IV.

II Method

II.1 Structural relaxation

To study possible structural distortions during the spin-state transition, we first perform the structural relaxation by adopting the DFT+U method as implemented in the Vienna ab-initio simulation package (VASP) code Kresse and Hafner (1993, 1994); Kresse and Furthmuller (1996a, b) using different spin-states as initial guess. The Perdew-Burke-Ernzerhof (PBE) functional Perdew et al. (1996) is used for the exchange-correlation functional within DFT. The DFT+U convergence is achieved using the plane-wave energy cut-off of 600eV and the $k-$point mesh of 8$\times$8$\times$8. The convergence of the structural relaxation is achieved if the atomic forces of all ions are smaller than 0.01eV/Å. Within DFT+U, we use the on-site Hubbard interaction U=6eV, which is obtained using the constraint DFT method Křápek et al. (2012), and the Hund’s coupling J=0.9eV. Although the ground-state of LaCoO3 is paramagnetic, DFT+U relaxation calculations are performed with the ferromagnetic configuration since correlations are included by imposing a long-range magnetic order in DFT+U.

LaCoO3 is a rhombohedral structure with the $R\bar{3}c$ symmetry containing two Co ions per unit cell. The experimental crystal volume $V$ is 56.0Å3 per formula unit at low temperatures Radaelli and Cheong (2002). We find that the DFT+U relaxation calculation with the PBE functional converging to LS produces the volume $V$ of 56.40Å3 per formula unit with the Co-O bond-length $a\sim 1.95$Å  (S1 structure in Table. 1), while the IS state (magnetic moment=2.2$\mu B$) results in the $2.7$% volume expansion resulting $V\sim 57.98$Å3 per formula unit and $a\sim 1.97$Å(S2 structure in Table. 1). The HS structure converged to a more expanded volume but the total energy is much higher than either LS or IS one, therefore we do not consider the HS structure in this paper. Interestingly, the MS imposing HS to one Co ion and LS to the other Co ion within the unit-cell produces the BD structure by lowering the crystal symmetry from $R\bar{3}c$ to $R\bar{3}$ (S3 structure in Table. 1). The crystal volume $V$ is similar to the IS volume ($V\sim 57.98$Å3) and the HS site (magnetic moment=3$\mu B$) becomes a Co-O LB site with $a\sim 1.99$Å  and the LS site (magnetic moment=0.3$\mu B$) is the SB with $a\sim 1.94$Å, resulting in the bond-length difference $\delta a\sim 0.05$Å. The $d-$occupancy ($N_{d}$) for relaxation results in all structures are close to $d^{7}$ (leaving one hole in surrounding oxygen atoms) although the nominal $d-$occupancy is $d^{6}$ for LaCoO3. This is due to the strongly covalent nature of the Co-O bonding in LaCoO3. The summary of the relaxed structure information is given in Table 1.

II.2 DFT+DMFT

Using LS, IS, and MS structures (denoted as S1, S2, and S3 in Table  1) obtained from DFT+U relaxations, we employ a charge-self-consistent DFT+DMFT method Park et al. (2014a); Singh et al. (2020) to study the nature and energetics of the spin-state transition. The DFT+DMFT method is implemented using the maximally localized Wannier functions (MLWFs) Mostofi et al. (2014); Marzari et al. (2012) as localized orbitals. First, we solve the non-spin-polarized Kohn-Sham (KS) equation within DFT using the VASP code. Then, we construct Co 3$d$ and O 2$p$ orbitals using MLWFs to represent the hybridization subspace for solving DMFT equations. The Co-O covalency effect can be treated within DMFT by including both $d$ and $p$ orbitals in the hybridization subspace. The $p-d$ Hamiltonian of the MLWF basis is constructed from the DFT bands in the hybridization energy window of 11eV. Then, the correlated subspace of Co 3$d$ orbitals is treated using the continuous time quantum Monte Carlo (CTQMC) Haule (2007); Gull et al. (2011) impurity solver by solving the DMFT self-consistent equations. The Hubbard U is 6eV and the Hund’s coupling J is 0.9eV within the Co 3$d$ shell while we compute the local self-energy for Co 3$d$ orbitals within DMFT calculations for the study of spin-state transitions. Within DMFT, we use temperatures from 100K to 1000K to study the temperature effect on the spectral function while 300K is used for most calculations unless specified otherwise. Here, both U and J are parameterized by Slater integrals ($F^{0}$, $F^{2}$, and $F^{4}$), namely U=$F^{0}$ and J=$(F^{2}+F^{4})/14$. The Coulomb interaction matrix elements with only density-density types are considered in CTQMC while the spin-flip and pair-hopping terms are neglected. Also, we chose the cartesian axes of Wannier orbitals to be aligned along the Co-O bonds so that the off-diagonal terms in the $d-$Hamiltonian can be much smaller (even zeros) than the on-site terms (see Appendix A). Therefore, the off-diagonal terms in the DMFT hybridization function can be also ignored within CTQMC.

The charge-self-consistency in DFT+DMFT can be achieved by updating the charge density using the DMFT local Green’s function while the DMFT loop is converged. The KS equation is solved again using the updated charge density and the new $p-d$ Hamiltonian is constructed using the updated MLWFs obtained from the KS solutions. The full charge-self-consistent DFT+DMFT loop is continued until both the charge density, the DMFT Green’s function $G^{loc}$, and the DMFT self-energy $\Sigma^{loc}$ are converged. More details of the DFT+DMFT implementation combining the projected augmented wave method in the DFT part and the formula for the charge update can be found in Ref. Park et al., 2014a. For the precise convergence, the energy cut-off of 600eV and the $k-$point mesh of 8$\times$8$\times$8 are used for DFT loops and a more dense $k-$point mesh of 30$\times$30$\times$30 is used within the hybridization window of Wannier orbitals for DMFT loops.

Once the DFT+DMFT self-consistency loop is converged, the total energy $E$ is computed using the following formula:

[TABLE]

where $E^{DFT}$ is the DFT energy computed using the charge density $\rho$ obtained within DFT+DMFT, $\epsilon_{i\mathbf{k}}$ is the DFT eigenvalues, $n_{i\mathbf{k}}$ is the diagonal part of the DMFT occupancy matrix element with the KS band index $i$ and the momentum $\mathbf{k}$, $W$ is the energy window for the hybridization subspace, and $n^{0}_{i\mathbf{k}}$ is the DFT occupancy matrix element with the KS band $i$ and the momentum $\mathbf{k}$. The potential energy $E^{POT}$ within DMFT is given by the Migdal-Galiski formula: $E^{pot}=\frac{1}{2}$Tr$[\Sigma^{loc}(i\omega)\cdot G^{loc}(i\omega)]$.

The double counting (DC) energy, $E^{DC}$ needs to be defined for beyond-DFT methods such as DFT+DMFT and DFT+U since the potential energy treated in the correlated subspace is already accounted as the part of the DFT energy and it needs to be subtracted from the total energy formula. Various DFT+DMFT calculations suggest that the DC potential, $V^{DC}(=\partial E^{DC}/\partial N_{d})$, smaller than the conventionally used fully-localized-limit (FLL) form Anisimov et al. (1991) can produce better agreements of energetics Park et al. (2014b, a), the metal-insulator transition Haule et al. (2014); Wang et al. (2012); Dang et al. (2014), and the $p-d$ orbital splitting Haule et al. (2014); Dang et al. (2014); Haule (2015) of oxides compared to experiments. Recently, it has been also shown that the more exact form of $V^{DC}$ within DFT+DMFT can be derived Haule (2015) and the exact $V^{DC}$ value can be close to the nominal DC form, in which the $d-$occupancy, $N_{d}$, in the FLL formulae is replaced to the nominal $d-$occupancy, $N_{d}^{0}$ ($d^{6}$ in the LaCoO3 case). More detailed discussions about different formula of double counting corrections are given in Appendix B.

In this paper, we propose the following form of $V^{DC}$ to allow the change between the FLL form and the nominal DC form by replacing $N_{d}$ to $\overline{N_{d}}$:

[TABLE]

where U and J are the same paramters which are defined above for the Slater-type interaction, $\overline{N_{d}}=N_{d}-\alpha$ where $N_{d}$ is the $d-$occupancy obtained self-consistently at each correlated site, and $\alpha$ is a parameter which can be tuned for obtaining different $E^{DC}$ and $V^{DC}$ values from the conventional FLL form. Our $V^{DC}$ formula can be derived from $E^{DC}$ ($V^{DC}=\partial E^{DC}/\partial N_{d}$) and allow site-dependent DC potentials. Our modified $V^{DC}$ form can recover the conventional FLL DC form by setting $\alpha$=0. By increasing $\alpha$, $V^{DC}$ approaches to the nominal $V^{DC}$ value as $\overline{N_{d}}$ becomes $N_{d}^{0}$ ($\alpha=N_{d}-N_{d}^{0}$). One should note that changing $V^{DC}$ with different $\alpha$ values can also tune the $p-d$ covalency effect by effectively shifting the $d$ orbital level. For example, a smaller $V^{DC}$ potential than the FLL DC potential will make the $d$ orbital level higher and the covalency effect weaker, resulting in a reduced $N_{d}$ value. The $p-d$ energy separation predicted by DFT or DFT+U with the FLL DC can be usually overbound and the physical role of $\alpha$ is to avoid this overbinding effectively by increasing the $p-d$ energy separation. In this paper, We studied the effect of the Co-O covalency on energetics and the nature of spin states in LaCoO3 by changing $V^{DC}$ potentials using different $\alpha$ values.

III Results

III.1 Energetics of spin-state transition: DFT+DMFT vs DFT+U

Here, we begin by showing the effect of different DC potential ($V^{DC}$) values on the energetics of spin-states in LaCoO3 computed using DFT+DMFT (Fig. 1 top panel) and DFT+U (Fig. 1 bottom panel). The $x-$axis shows $N_{d}$ values obtained by changing $V^{dc}$ in Eq. 3 using different $\alpha$ values, namely $\alpha$=0, 0.2, 0.4, 0.6, and 0.7. Two vertical dashed lines indicate the $N_{d}$ values obtained using the FLL DC formulae (the right line) and the nominal DC formulae (the left line). The calculations with the FLL DC ($\alpha$=0) converge to $N_{d}\sim 7.2$ for DFT+DMFT and 7.1 for DFT+U while those with the nominal DC result in $N_{d}\sim 6.68$ for DFT+DMFT and 6.65 for DFT+U. The $y-$axis indicates the total energy difference between the expanded volume structures (S2 and S3; $V\sim 57.98$Å3) accompanying excited spin-states and the S1 structure ($V\sim 56.4$Å3). Here, DFT+U energies are computed by adopting the MLWFs as correlated orbitals consistently with DFT+DMFT calculations, therefore the difference between DFT+DMFT and DFT+U results is attributed purely to the dynamical correlation effect beyond DFT+U. Also, the previous study of DFT+U calculations using different choices of orbitals (MLWFs vs projectors) has shown that results of the structural phase diagram of nickelates are almost the same as long as the same U and J values are used Park et al. (2015). Therefore we expect that the energetics results reported here will not depend much on the choice of correlated orbitals.

The energetics obtained using two methods show noticeable differences depending on $N_{d}$ (the Co-O covalency effect). The DFT+DMFT energy difference between two expanded volume structures (S2($\delta$a=0Å)and S3($\delta$a=0.05Å)) is much smaller (maximally 20meV) than the DFT+U energy difference while DFT+U solutions converge to various meta-stable states for a given structure (LS, IS, and HS for S2 and LS-HS and LS-IS for S3) and the DFT+U energy depends sensitively on $N_{d}$ and spin-states. This difference in energetics arises since the spin-state within DFT+DMFT is described as a multi-configurational state with a mixture of various spin states showing smooth crossover (see Fig. 2) while the DFT+U solution is based on a single-determinant form and produces meta-stable states depending on given structure and $N_{d}$.

Within DFT+DMFT, the S2 structure (triangular dot) is more stable when $N_{d}>7.0$ and the LS state is dominant while the S3 BD structure (circular dot) becomes stable when the spin-state transition to LS-HS occurs ($N_{d}<7.0$). Two structures become energetically almost the same when $N_{d}$ is further reduced ($N_{d}\sim 6.7$) since the energy of the S2 structure gets lowered as HS is more excited at smaller $N_{d}$. The dependence of DFT+U energetics on $N_{d}$ favors the MS (LS-HS) state although it behaves qualitatively similar to DFT+DMFT as LS with the S2 structure is stable when $N_{d}>6.8$ and LS-HS with the S3 BD structure becomes rapidly stable when $N_{d}$ is further reduced. The energy of HS with the S2 structure is also rapidly decreasing at smaller $N_{d}$ as the Hund’s coupling lowers the energy in the HS state. Within DFT+U, energetics of other meta-stable states including LS and LS-IS do not depend much on $N_{d}$ while IS is not favored as $N_{d}$ is reduced. Our calculations show that the spin-state transition from the S1 (LS) structure to excited spin-states with the expanded volume $(E-E(S1)<0)$ occurs when $N_{d}$ becomes smaller ($<6.9$ for DFT+DMFT and $<6.8$ for DFT+U) than the FLL DC result. This is consistent with the nominal $V^{DC}$ result (the left vertical dashed lines; $N_{d}\sim 6.68$ for DFT+DMFT) and due to the fact that higher spin states are more excited at smaller $N_{d}$ as will be shown in Fig. 2.

III.2 Nature of paramagnetic states

Now we turn to the nature of paramagnetic states in LaCoO3 obtained within DFT+DMFT. Fig. 2 displays the spin $|S_{z}|$ probabilities sampled using CTQMC as a function of $N_{d}$ for different structures. The multi-configurational nature of the paramagnetic state treated in DFT+DMFT means that various $|S_{z}|$ states ($S_{z}=0.0\sim 2.0$) contribute to the solution. In all structures, the charge $d^{7}$ state with $|S_{z}|=$ 0.5 and 1.5 is not negligible in addition to the nominal $d^{6}$ state with $|S_{z}|=$ 0, 1, and 2, therefore LaCoO3 is strongly covalent with the dynamically fluctuating nature of spin and charge states. When $N_{d}>7.0$, LS with $|S_{z}|=$ 0 and 0.5 has the highest probability for all structures, consistently with the fact that the S1 structure ($V=56.40$Å3) is energetically stable. As $N_{d}$ is reduced, LS is still dominant for the S1 structure while HS with $|S_{z}|=$ 1.5 and 2 becomes more excited for the S2 structure ($V=57.98$Å3). The nature of paramagnetic state in the S2 structure, which is energetically stable when $N_{d}<6.7$, is characterized by a mixture of both LS and HS, fluctuating dynamically with non-negligible probabilities. The IS state ($|S_{z}|=$1.0) is strongly suppressed for all structures. In the BD structure, HS with $|S_{z}|=$ 1.5 and 2 becomes rapidly favored in the LB site as $N_{d}$ is reduced while LS with $|S_{z}|=$ 0 and 0.5 is always dominant in the SB site. These DFT+DMFT results of paramagnetic states exhibiting the mixture of different spin states are distinct from the DFT+U results where meta-stable solutions of spin states are found.

III.3 Charge ordering induced by spin-state ordering

The strongly coupled spin and charge degrees of freedom also produce an intriguing rock-salt type charge ordering state induced from the MS state in the S3 structure. Fig. 3 shows that the DFT+DMFT charge in the LB site (triangular dots) gets smaller more rapidly than the SB charge (circular dots) as the average $\langle N_{d}\rangle$ becomes reduced (the overall Co-O covalency is reduced). This is because HS in the LB site favors the $|S_{z}|=2$ with $d^{6}$ state while the $|S_{z}|=0.5$ with $d^{7}$ state is dominant in the SB site as the Co-O covalency remains strong. This enhancement of charge ordering when $N_{d}<$6.9 is also consistent with the rapid increase of the HS probability in the LB site as shown in Fig. 2. The experimental evidence of charge ordering with the concomitant spin-state ordering in LaCoO3 has been reported in the tensile-strained film Sterbinsky et al. (2018). The ground-state of this charge ordered LaCoO3 is insulating and the nature of the LB state is a Mott insulator while the SB state is a band insulator due to the strong Co-O hybridization (see Fig. 8). This MS state in the S3 (BD) structure is somewhat reminiscent of the site-selective Mott physics occurring in nickelates in which the LB site is Mott insulating with the $d^{8}$ state while the SB site is a covalent insulator hybridized with O hole states as $d^{8}L^{2}$. Park et al. (2012) Although charge ordering between two Ni sites in nickelates is not important to induce the insulating state in nickelates, charge ordering ($\sim 0.2$) in LaCoO3 is naturally induced from the spin-state ordering. Also, nickelates are negative charge-transfer insulators meaning that almost one electron is donated to each Ni ion from surrounding O ions while the Co ion in LaCoO3 favors a mixed-valence state with $N_{d}=6.6\sim 6.7$ as the charge transfer from the O hole is smaller than nickelates.

III.4 Density of states

In this subsection, we study the correlation, structure, $p-d$ covalency, and temperature effects on the spectral functions of LaCoO3. First, we show the density of states (DOS) computed using DFT paramagnetic calculations (no spin polarization) in Fig. 4. The DOS computed for the different structures shows very similar features. All structures are metallic without correlations and the S1 structure has slightly larger band-width than other structures due to the smaller crystal volume. Co $t_{2g}$ states are all occupied below the Fermi energy Co $e_{g}$ states are widely spread due to the strong mixing with the O $p$ states. The LB and SB sites in the S3 structures exhibit the similar DOS without charge ordering ($N_{d}\sim 7.2$).

Now, we include strong correlation effects in three structures (S1, S2, and S3) within DFT+DMFT and compute the DOS using the paramagnetic spin configuration in Fig. 5. We see dramatic changes of electronic structures due to correlation effects as well as the structural changes from the volume expansion and the Co-O bond disproportionation. The DOS in the S1 structure (Fig. 5 a) shows that the band gap is almost 0.6eV, which is consistent with the optical gap measurement Chainani et al. (1992). The $t_{2g}$ state is almost occupied while the $e_{g}$ state is mostly unoccupied, as expected for the LS state. As the volume is expanded, the spin-state transition to higher spins occurs continuously and the $t_{2g}$ state begins to be unoccupied while more $e_{g}$ orbitals are occupied (Fig. 5 b). In the S1 structure, the major optical transition occurs from the Co $t_{2g}$ valence band peak to the Co $e_{g}$ conduction band peak while the transition between the same $t_{2g}$ bands dominates in the S2 structure. Therefore, the position of the first major peak in the optical conductivity will be reduced as the crystal volume expands at higher temperature, which is consistent with the optical conductivity measurement in experiment. Moreover, our calculation shows that the S3 structure with the MS state produces the smaller gap than the S2 structure, therefore the optical gap can be further reduced as temperature is raised presumably with more populated MS states. The nature of the paramagnetic insulating state in the S2 structure is a strongly correlated band insulator driven by both the Co-O hybridization due to the covalency and the electron localization of the Co d orbitals due to increased HS states. This physics is captured in the imaginary part of the self-energy (see Fig. 8 b) showing the emergence of the sharp pole outside the hybridization gap near the Fermi energy. This paramagnetic insulating state is represented as a mixture of fluctuating HS and LS states (see Fig. 2) and it is distinct from the normal Mott insulator driven by the pure electron localization which is expected from the HS state in a d6 configuration. In the S3 structure, the LB site becomes higher spin state with more unoccupied $t_{2g}$ states while the SB site remains as LS with a similar gap size as the homogeneous LS gap. Here, the $\alpha$ value in $V^{DC}$ (Eq. 3) is set to 0.7 resulting $N_{d}\sim 6.68$ for all structures. Different $V^{DC}$ values will change the relative position of the O $p$ peak from the Fermi energy as $N_{d}$ is also varied. Our O $p$ top peak position is located at -2eV when $N_{d}\sim 6.68$ and this peak position is consistent with the measured X-ray photo-emission spectra Abbate et al. (1993) validating the $\alpha$ value we used.

To clarify the role of the Co-O covalency (parametrized by $N_{d}$) on electronic structure, we show the DOS obtained using different $\alpha$ values for the S2 structure in Fig. 6. As the Co-O covalency is enhanced ($N_{d}$ is increased), the O $p$ peaks move closer to Co $d$ states near the Fermi energy. As a result, the spectral gap becomes smaller due to the larger Co-O hybridization and the ground-state is eventually metallic when $N_{d}\sim 7.1$. Therefore, the insulator-to-metal transition occurs as the correlation in the Co ion is reduced due to the increased Co-O covalency and, at the same time, the spin states are less excited as Co $t_{2g}$ states are more occupied and $e_{g}$ states are more unoccupied. This DOS result is also consistent with the occupation probability data of the S2 structure (Fig. 2 upper right panel) showing that the HS states ($|S_{z}|=1.5$ and 2.0) are suppressed and the LS and IS states ($|S_{z}|=0.0$, 0.5, and 1.0) are gradually increasing as $N_{d}$ increases. Our results suggest that the Co-O hybridization due to the covalency can play an important role in explaining the metal-insulator transition and the spin-state transition in LaCoO3.

Until now, our DFT+DMFT calculations have been performed at the fixed temperature $(\sim 300K)$. Experimentally, LaCoO3 also exhibits the insulator-to-metal transition as temperature is raised above near $T=400K$ but the nature of this transition has not been clarified. To better understand the role of temperature on the metal-insulator transtion in LaCoO3, we plot the DOS of LaCoO3 at the fixed S2 structure and different temperatures (Fig. 7). The DOS data computed at both $100K$ and $300K$ show similar features although the spectral gap at $300K$ is slightly larger than one at $100K$. The similarity of the DOS at low temperatures indicates that our spin-state calculations at $300K$ can be similarly reproduced at $100K$ (the experimental $T_{c}$ of the spin-state transition). More importantly, the strong variation of the DOS and spin-states depending on structures in Fig. 5 means that the structural changes should be incorporated for the better description of electronic structure while temperatures are varied. Nevertheless, the spectral gap becomes smaller as temperature is raised above $500K$ and the Co $t_{2g}$ state becomes the incoherent metallic state at around $1000K$.

III.5 Self-energy data

To study the nature of insulating states occuring in different stuctures, we show both the real part ($Re\Sigma$; Fig. 8 top panel) and the imaginary part ($Im\Sigma$; Fig. 8 bottom panel) of self-energies on the real axis, which are used to compute the DOS in Fig. 5. The S1 structure (Fig. 8a) exhibits rather small $Im\Sigma$ values in both $e_{g}$ and $t_{2g}$ orbitals when $\omega\sim 0$ while a sharp pole develops at $\omega\sim 300meV$ in the $t_{2g}$ orbital for both real and imaginary parts. This diverging nature of the self-energy indicates Co $d$ orbitals are still correlated even at the S1 structure with LS. This self-energy nature of LS is different from the typical band insulating nature of LS which is expected in a $d^{6}$ configuration. The nature of our insulating state is the correlated band insulator driven by both the Co-O hybridization and the electron correlation. The electron correlation is encoded in the sharp and narrow pole structure of the imaginary part of the self-energy developed outside the hybridization gap while the imaginary part of the self-energy is still zero at the Fermi energy. As the structure changes from S1 to S2 along with the volume expansion, the strength of poles becomes enhanced and the positions of the poles get close to the Fermi energy in both $e_{g}$ and $t_{2g}$ orbitals. This clearly shows that correlations are enhanced due to the increased higher spin probabilities (see Fig. 2) as the volume is expanded in LaCoO3. The insulating nature of the S3 structure shows the “site-selective” Mott physics, as the Co-O LB site undergoes a Mott transition with the diverging self-energies at the Fermi energy while the SB site behaves as a band insulator with the small imaginary part and the flat real part of self-energies originated from the strong Co-O covalency.

To understand the nature of the insulator-to-metal transition in LaCoO3 as temperature is raised, we also plot the $Im\Sigma$ computed for the S2 structure as a function of temperature in Fig. 9. Results of $Im\Sigma$ at $T$=100K show that poles are developed due to correlations at $\omega=$ -450meV for the $e_{g}$ orbital and at $\omega=$ 150meV for the $t_{2g}$ orbital. The strengths of these poles are enhanced at $T$=300K as the S2 structure can occupy higher spin states due to the spin-state transition and electron localization can be enhanced at this temperature. This physics is different from the typical Mott insulator without any spin-state transition, in which the electron localization is usually stronger at lower temperature. As the temperature is raised even above 300K, the pole strengths are reduced and correlations become weaker for both $e_{g}$ and $t_{2g}$ orbitals resulting the metallic phase obtained in Fig. 7.

III.6 Charge-self-consistency effect

Finally, we show the charge-self-consistency effect in DFT+DMFT on the energetics and electronic structure in LaCoO3. Fig. 10 shows the energetics of spin-state transition in LaCoO3 as a function of $N_{d}$ comparing charge-self-consistent (top panel) and non-charge-self-consistent (bottom panel) DFT+DMFT calculations. Here, the non-charge-self-consistent calculation means that the charge density ($\rho$ in Eq. 1) is fixed to the DFT one while the DMFT local Green’s function ($G^{loc}$ in Eq. 1) is obtained by converging DMFT self-consistent equations. As a result, the $V^{DC}$ potential is fixed during the DFT+DMFT loop since the charge-density is not updated and $V^{DC}$ is a function of the charge-density ($N_{d}$). Therefore, different $N_{d}$ results in non-charge-self-consistent DFT+DMFT (bottom panel) are obtained by shifting the $V^{DC}$ potential as data points are changed. However, in charge-self-consistent DFT+DMFT (top panel), $V^{DC}$ is computed using the Eq. 3 with the self-consistently determined charge density ($N_{d}$) and the corresponding $N_{d}$ is obtained from the DMFT Green’s function.

The energetics of the S2 structure (homogeneous spin states) show the very similar behavior as a function of $N_{d}$ between charge-self-consistent and non-charge-self-consistent calculations. In the case of the MS state for the S3 structure, non-charge-self-consistent DFT+DMFT energetics overestimate the tendency toward the spin-state ordering as the energy difference between the excited spin-state and the low-spin (S1) state becomes much lower ($\sim$-150meV) than the charge-self-consistent energetics ($\sim$-30meV) when $N_{d}\sim$6.7. The energetics of the MS state also depend on $N_{d}$ much sensitively in the non-charge-self-consistent case. Therefore, the main effect of charge-self-consistency in DFT+DMFT is to reduce the spin-state ordering effect in the S3 structure and the energetics between homogeneous and MS spin-states in the S2 and S3 structures become very close.

To further investigate the effect of charge-self-consistency on electronic structure of LaCoO3, we display the $d-$occupancy, $N_{d}$ in Fig. 11 obtained for both LB and SB sites in the S3 structure computed using charge-self-consistent DFT+DMFT (solid lines) and non-charge-self-consistent DFT+DMFT (dashed lines) as a function of the average $N_{d}$ between two sites. Without the charge update (dashed lines), the $N_{d}$ difference between two sites is more enhanced since the SB site occupies more $d-$orbitals while the LB site takes even less $N_{d}$ compared to the charge-self-consistent results (solid lines) across different $\langle N_{d}\rangle$ values. Therefore, charge-self-consistency within DFT+DMFT reduces the tendency toward both spin-state and charge orderings between correlated Co sites, as a result, the energy difference between different excited spin-states is also much decreased within the charge-self-consistent calculation.

Here, we fix the U and J values during the DFT+DMFT calculations, therefore the main difference between charge-self-consistent and non-charge-self-consistent results is originated from the change of the Wannier Hamiltonian due to the charge update (see Appendix A). The main role of the charge update is to decrease the crystal field splitting in the S1 and S2 structures and to promote the homogeneous spin-state transition. At the same time, the spin-state and charge ordering effects in the S3 structure has been reduced as the crystal field splitting in the LB site has been increased while it has been decreased in the SB site acting counterintuitively on the structural effect. Moreover, the overall $d-$orbital level in the LB site is substantially lowered than the level in the SB site due to the charge update, again decreasing the charge ordering effect. Therefore, the charge-self-consistent effect compensates for the spin-state and charge ordering effects in the S3 structure while it favors the homogeneous spin-state transition, and as a result, the energy difference between the spin-state ordering and the homogeneous spin state has been much reduced within the charge-self-consistent DFT+DMFT calculations.

IV Conclusion

In conclusion, we adopt the charge-self-consistent DFT+DMFT method to study the nature and energetics of both homogeneous and MS states in covalent LaCoO3. We find that structural changes during the spin-state transition are important to understand both energetics and electronic structure in LaCoO3. As the crystal volume is expanded, the occupation probability of higher spin state in a Co ion increases and a sharp and narrow peak in the imaginary part of the self-energy develops outside the Co-O hybridization gap due to the enhanced electron correlation. This paramagnetic insulating state also exhibit a multi-configurational mixture of both HS and LS states with strong spin and charge fluctuations. This DMFT result is different from the static DFT+U result showing various meta-stable solutions. The MS state in the BD structure accompanies the HS state in a LB Co site and the LS state in a SB site. Charge ordering is also induced from this spin-state ordering since HS favors the $d^{6}$ charge state while LS occupies more $d^{7}$ states as covalently bonded with O ions. The LB site with HS becomes a Mott insulator while the SB site with LS behaves as a band insulator. Our DFT+DMFT calculation reveals that energetics between homogeneous and MS states are very close while DFT+U energetics depend sensitively on spin-states and structures favoring the MS state.

We also find that the Co-O covalency plays a crucial role in electronic structure of LaCoO3. Changing the double-counting potential in DFT+DMFT can lead to a different $N_{d}$ value tuning the covalency effectively and both energetics and one-particle spectra are more consistent with experiments when $N_{d}\sim 6.7$. Increasing the Co-O covalency produces the insulator-to-metal transition faving the LS state. Increasing temperature beyond 500K also reduces the correlation effect and drives the insulator-to-metal transition, similarly as the experiment. The charge-self-consistency effect within DFT+DMFT reduces the tendency toward spin-state and charge orderings in LaCoO3 producing close energetics between different structures compared to the non-charge-self-consistent calculation.

Acknowledgments

H. Park acknowledges helpful discussions with Andrew Millis. H. Park and R. Nanguneri are supported by the U.S. Deptartment of Energy (DOE), Office of Science, Basic Energy Sciences (BES), Materials Sciences and Engineering Division. A. Ngo is supported by the Vehicle Technologies Office (VTO), Department of Energy (DOE), USA. We also gratefully acknowledge the computing resources provided on Bebop, a high-performance computing cluster operated by the Laboratory Computing Resource Center at Argonne National Laboratory.

Appendix A: The $d-$orbital Hamiltonian

In this Appendix, we show the matrix elements of the multi-orbital Hamiltonian, $\hat{H}^{dd}$, in the Co 3$d$ shell at each correlated site for three structures (S1, S2, S3 LB, and S3 SB) used in the DFT+DMFT calculations of this paper. First, we provide $\hat{H}^{dd}$ (in units of eV) represented using MLWF $d-$orbitals obtained from the DFT solution, which is relevant to the non-charge-self-consistent DFT+DMFT calculation. The orbital order representing $\hat{H}^{dd}$ is $d^{3z^{2}-r^{2}}$, $d^{xz}$, $d^{yz}$, $d^{x^{2}-y^{2}}$, and $d^{xy}$.

[TABLE]

[TABLE]

[TABLE]

[TABLE]

First of all, the local axis for initial projections of MLWFs are chosen to be aligned to the Co-O octahedron axis, therefore the off-diagonal elements are close to zeros. The crystal-field splittings between $e_{g}$ and $t_{2g}$ orbitals are $\sim 0.67$eV for S1, $\sim 0.63$eV for S2, $\sim 0.57$eV for S3 LB, and $\sim 0.68$eV for S3 SB. As expected, the Co-O bond length can tune the crystal-field splitting (the energy difference between $e_{g}$ and $t_{2g}$ orbitals). Namely, the long-bond Co ion has the smaller splitting and favors higher spin-states.

The charge-self-consistency effect in DFT+DMFT produces a new charge density $\rho$ which is different from the original DFT $\rho$. As a result, $\hat{H}^{dd}$ is computed from the updated $\rho$ and is also changed. Here, we show below $\hat{H}^{dd}$ obtained from charge-self-consistent DFT+DMFT calculations using U=6eV, J=0.9eV, and the DC parameter $\alpha$=0.6. In this case, the crystal-field splittings are $\sim 0.6$eV for S1 and S2, $\sim 0.66$eV for S3 LB, and $\sim 0.54$eV for S3 SB. The main effect of the charge-self-consistency on $\hat{H}^{dd}$ is to reduce the crystal-field splittings for S1 and S2 structures compared to the non-charge-self-consistency case, therefore it promotes the homogeneous spin-state transition. However, the charge-self-consistency effect also compensates for the spin-state and charge orderings as the crystal-field splitting has been increased for the LB site while it is decreased for the SB site in the S3 structure. The average $d-$orbital level for the LB site becomes also much lower than the SB site, again compensating for charge ordering.

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Appendix B: Double counting correction

A frequently used expression of the double-counting energy $E^{DC}$ is the fully localized limit (FLL) form Anisimov et al. (1991) which has been adopted frequently in DFT+U.

[TABLE]

where $N_{d}$ is the occupancy of the correlated site and can be obtained as the result of self-consistent DFT+DMFT or DFT+U calculations. Therefore, the $V^{DC}$ potential depends on the correlated site since $N_{d}$ is site-dependent.

Recently, it has been shown that the exact form of $V^{DC}$ within DFT+DMFT Haule (2015) can be computed and the formulae should be close to the nominal DC form, where $N_{d}$ in Eq. 13 is replaced to the nominal $d-$occupancy in the atomic limit, $N_{d}^{0}$, which is site-independent.

[TABLE]

The hybridization of $d$ and $p$ orbitals in transition metal oxides means that the resulting $d-$occupancy $N_{d}$ will be larger than the nominal value $N_{d}^{0}$ ($N_{d}>N_{d}^{0}$), therefore the nominal $V^{DC}$ will be always smaller than the FLL $V^{DC}$. This smaller $V^{DC}$ potential reduces the covalency effect between $d$ and $p$ orbitals.

In this paper, we use the following forms of $E^{DC}$ and $V^{DC}$ by modifying the FLL forms to allow the tuning of $V^{DC}$ for changing the Co-O covalency effect.

[TABLE]

where $\overline{N_{d}}=N_{d}-\alpha$ with a parameter $\alpha$. Our $V^{DC}$ formula can be derived from $E^{DC}$ ($V^{DC}=\partial E^{DC}/\partial N_{d}$) and allow site-dependent potentials similarly as the FLL form. The conventional FLL DC form is recovered by setting $\alpha$=0. By increasing $\alpha$, $V^{DC}$ can be close to the nominal $V^{DC}$ value as $\overline{N_{d}}$ approaches to $N_{d}^{0}$ ($\alpha=N_{d}-N_{d}^{0}$).

Another modified form of $V^{DC}$ for DFT+DMFT was also suggested as below since the $U$ value used in the FLL form can be smaller than the Hubbard $U$ ($U^{\prime}<U$) to allow the smaller $V^{DC}$ potential than the FLL one:

[TABLE]

where $U^{\prime}=U-\alpha$ with a parameter $\alpha$. Here, the role of $\alpha$ is the same as the one in Eq. 16. Namely, the covalency effect can be reduced by increasing the $\alpha$ value. It has been shown that using $U^{\prime}=U-0.2$eV can successfully reproduce the structural and electronic phase diagram of rare-earth nickelates Park et al. (2014b, a).