Plane-wave based electronic structure calculations for correlated materials using dynamical mean-field theory and projected local orbitals
B. Amadon, F. Lechermann, A. Georges, F. Jollet, T. O. Wehling, A. I. Lichtenstein
aa r X i v : . [ c ond - m a t . s t r- e l ] J a n Plane-wave based electronic structure calculations for correlated materialsusing dynamical mean-field theory and projected local orbitals
B. Amadon, F. Lechermann, A. Georges, F. Jollet, T. O. Wehling, and A. I. Lichtenstein CEA, D´epartement de Physique Th´eorique et Appliqu´ee,Bruy`eres-le-Chˆatel, 91297 Arpajon Cedex, France I. Institut f¨ur Theoretische Physik, Universit¨at Hamburg, D-20355 Hamburg, Germany Centre de Physique Th´eorique, ´Ecole Polytechnique, 91128 Palaiseau Cedex, France
The description of realistic strongly correlated systems has recently advanced through the combi-nation of density functional theory in the local density approximation (LDA) and dynamical meanfield theory (DMFT). This LDA+DMFT method is able to treat both strongly correlated insulatorsand metals. Several interfaces between LDA and DMFT have been used, such as (N-th order) LinearMuffin Tin Orbitals or Maximally localized Wannier Functions. Such schemes are however eithercomplex in use or additional simplifications are often performed (i.e., the atomic sphere approxi-mation). We present an alternative implementation of LDA+DMFT, which keeps the precision ofthe Wannier implementation, but which is lighter. It relies on the projection of localized orbitalsonto a restricted set of Kohn-Sham states to define the correlated subspace. The method is imple-mented within the Projector Augmented Wave (PAW) and within the Mixed Basis Pseudopotential(MBPP) frameworks. This opens the way to electronic structure calculations within LDA+DMFTfor more complex structures with the precision of an all-electron method. We present an applicationto two correlated systems, namely SrVO and β -NiS (a charge-transfer material), including ligandstates in the basis-set. The results are compared to calculations done with Maximally LocalizedWannier functions, and the physical features appearing in the orbitally resolved spectral functionsare discussed. PACS numbers: 71.10.Fd,71.30.+h,71.15.Ap,71.15.Mb
I. INTRODUCTION
The description of strong electronic correlations in arealistic framework has become an issue of major im-portance in current condensed matter research. Due tothe fast progress in the preparation of novel materials,especially in effectively reduced dimensions, and the ad-vances in experimental techniques in order to probe suchsystems, providing an adequate theoretical formalismthat can handle explicit many-body effects in a material-specific context is crucial. It is evident that standarddensity functional theory (DFT) in the local density ap-proximation (LDA) cannot meet those demands and, atleast, has to ally with manifest many-body techniques.In this respect, the combination of DFT-LDA with theDynamical Mean-Field Theory (DMFT), the so calledLDA+DMFT approach, has proven to be an invaluablemethod to face the challenge. Already numerous stud-ies within the LDA+DMFT framework have shown thatthis theory is capable of describing the effects of strongcorrelations in a realistic context, such as Mott transi-tions and volume-collapse transitions in d and f -electronsystems, effective mass enhancement, local moment for-mation and magnetism, spectral-weight transfers, orbitalphysics, etc...In view of these successes, it is crucial to push furtherthe range of applicability of LDA+DMFT methods. Tothis aim, implementing LDA+DMFT within those elec-tronic structure methods which are highly accurate andallow for the possible treatment of larger systems, is cer-tainly an important endeavor. One of the main aims of this paper is to report on such new implementations,within a projector augmented-wave (PAW) and a mixed-basis pseudopotential (MBPP) method. The formalismused in this paper is actually quite general and allows foran implementation of LDA+DMFT in a very large classof electronic structure methods. We present this formal-ism in an arbitrary basis set, which should make this taskeasier for other implementations. However, in our specificimplementation, we make an extensive use of the Blochbasis set of Kohn-Sham orbitals. In order to implementLDA+DMFT, we use local orbitals constructed by pro-jecting atomic-like orbitals onto a restricted set of Blochstates, a strategy similar to that introduced by Anisimov et al. .Another important route is to investigate how to op-timize the application of LDA+DMFT for a specific ma-terial, or class of materials of interest. Indeed, it is im-portant to realize that, while the results of course donot depend on the specific basis set used in the calcu-lation, they will actually depend on the specific set of local orbitals for which many-body effects will be takeninto account, or more precisely (following the terminol-ogy introduced in [2]) on the “correlated” subspace C ofthe full Hilbert space spanned by these orbitals. Indeed,the DMFT treatment applies local interaction terms tothose orbitals only, and furthermore neglects all non-localcomponents of the self-energy. This notion of locality isdefined with respect to the specific choice of the local or-bitals defining C , and the accuracy of the DMFT approx-imation cannot be expected to be identical for differentchoices.Early implementations of LDA+DMFT used thelinear muffin-tin orbital (LMTO) framework, and thecorrelated orbitals were frequently identified with a spe-cific subset of the LMTO basis functions, having d - or f -character. There is no specific reason to pick the localorbitals as members of the basis set however, and a setof atomic-like functions may prove to be better from thephysical point of view. Recently, several works have useddifferent kinds of Wannier functions (WFs) as correlatedorbitals, starting with the work of Pavarini et al. us-ing Nth-order muffin-tin orbitals (NMTO). Anisimov etal. used WFs constructed from a projection on a subsetof the Bloch functions , and Lechermann et al. usedmaximally-localized WFs and compared the results tothose using NMTOs.The energy window spanned by the basis functionswhich are retained in the implementation, and the spa-tial extension of the local orbitals defining C are impor-tant physical issues for the description of a given ma-terials. Those issues become particularly important forcharge-transfer materials (e.g. late transition-metal ox-ides, sulfides or selenides) in which the ligand states mustbe kept in the basis set in order to reach a satisfactoryphysical picture. In the present article, we address theseissues and provide explicit comparisons between calcula-tions performed with more spatially extended local or-bitals (hence spanning a smaller energy window whenprojected onto Bloch functions) and more localized localorbitals (hence spanning a wider energy window).This article is organized as follows. Sec. II is devotedto a presentation of the general theoretical framework.Secs. III and IV are devoted to applications to two com-pounds, which are considered as tests for the method: thetransition-metal oxide SrVO in Sec. III and the charge-transfer sulfide NiS in Sec. IV. II. THEORETICAL FRAMEWORKA. LDA+DMFT formalism in an arbitrary basisset
In order to implement DMFT within realistic elec-tronic structure calculations of correlated-electron ma-terials, one has to set up a formalism which keepstrack of the real-space (i.e., quantum chemical) and thereciprocal-space (i.e., solid-state) aspects on an equalfooting, while being computationally efficient. Follow-ing the work of Lechermann et al. , we therefore distin-guish between the complete basis set {| B k α i} in whichthe full electronic structure problem on a lattice is for-mulated (and accordingly the lattice Green’s function isrepresented), and local orbitals | χ R m i which span a “cor-related” subspace C of the total Hilbert space. Many-body corrections beyond LDA will be considered insidethis subspace. The index α labels the basis functions foreach wave vector k in the Brillouin zone (BZ). The in-dex R denotes the correlated atom within the primitive unit cell, around which the local orbital | χ R m i is centered,and m =1 , . . . , M is an orbital index within the correlatedsubset. Projection on that subset, for atom type R , isdone with the projection operator:ˆ P ( C ) R ≡ X m ∈C | χ R m ih χ R m | . (1)In Ref. 2, it was shown that the DMFT self-consistencycondition, which relates the impurity Green’s function G R imp to the Green’s function of the solid computed lo-cally on atom R , reads: G R , imp mm ′ ( iω n ) = X k X αα ′ h χ Rk m | B k α ih B k α ′ | χ Rk m ′ i ×× n [ iω n + µ − H KS ( k ) − ∆Σ ( k , iω n )] − o αα ′ , (2)In this expression: | χ Rk m i = X T e i k · ( T + R ) | χ RT m i (3)denotes the Bloch transform of the local orbitals whereby T denotes the Bravais lattice translation vectors. Notethat the object in the second line of Eq. (2) is, ofcourse, nothing else than the full lattice Green’s function G αα ′ ( k , iω n ) in the chosen {| B k α i} basis. The Kohn-Sham (KS) hamiltonian H KS ( k ) is obtained by solvingthe Kohn-Sham equations, which yields eigenvalues ε k ν and Bloch wave-functions | Ψ k ν i ( ν is the band index). Itcan be expressed in the {| B k α i} basis set as: H KS ,αα ′ ( k ) = X ν h B k α | Ψ k ν i ε k ν h Ψ k ν | B k α ′ i . (4)In order to obtain the self-energy for the full solid, one hasto promote (“upfold”) the DMFT impurity self-energyΣ imp mm ′ to the lattice via :∆Σ αα ′ ( k , iω n ) = X R X mm ′ h B k α | χ Rk m ih χ Rk m ′ | B k α ′ i ×× h Σ imp mm ′ ( iω n ) − Σ dc mm ′ i , (5)whereby a double-counting correction Σ dc mm ′ , taking careof correlation effects already accounted for in the LDAhamiltonian (see appendix B) has to be included in thegeneral case.The above equations form the general LDA+DMFTframework, in a general arbitrary basis set. Any specificimplementation must then make a choice for:i) The set of local orbitals {| χ R m i} spanning the cor-related subspaceii) The specific basis set {| B k α i} in which these equa-tions will be implementedIt is important to realize that the results will certainlydepend on the specific choice of {| χ R m i} : the quality ofthe DMFT approximation will indeed depend on howthe local orbitals are picked such as to minimize non-local contributions. In contrast, for a given choice of {| χ R m i} ’s, the results should be in principle independentof the basis set {| B k α i} which is chosen for the imple-mentation. However, in practice, considerations of nu-merical efficiency do limit the size of the basis whichcan be handled. Indeed, Eq. (2) involves inverting amatrix (at each k -point and for each frequency) of size N b × N b , in which N b is the number of basis functionswhich are eventually retained. Hence, in practice, onewill restrict the basis set to a certain set W of bands,as will be described in more details below. Regardingthe choice of local orbitals, we recall that, in Ref. [2],different kinds of Wannier functions were used and com-pared to one another, namely maximally-localized and N th-order linear-muffin-tin orbitals (NMTO) . The con-struction of such functions requires rather sophisticatedprocedures. In the present work we use a somewhatlighter implementation, with the same demands on accu-racy, by constructing the {| χ R m i} ’s out of entities that arealready existing in most of the common band-structurecodes, namely using the decomposition of local atomic-like orbitals onto the basis function retained in the set W . This is very similar in spirit to the construction pro-posed by Anisimov et al. in the LMTO framework. Thepresent construction, and the basic equations developedin this section, are fairly general however and can be im-plemented in an arbitrary electronic structure code. B. LDA+DMFT formalism in the Bloch basis set
In this section, we focus on a very simple choice for thebasis set, namely the Bloch basis itself {| B k α i} = {| Ψ k ν i} .This basis is most conveniently used, since it is a directoutput of any DFT-LDA calculation and furthermore di-agonalizes the KS hamiltonian H KS νν ′ ( k ) = δ νν ′ ε k ν .The basic LDA+DMFT equations in the Bloch basisset are easily written up, using the projection matrixelements of the local orbitals onto the Bloch functions,defined as: P R mν ( k ) ≡ h χ Rk m | Ψ k ν i , P R ∗ νm ( k ) ≡ h Ψ k ν | χ Rk m i . (6)Equations (2),(5) then read: G R , imp mm ′ ( iω n ) = X k ,νν ′ P R mν ( k ) G bl νν ′ ( k , iω n ) P R ∗ ν ′ m ′ ( k ) (7)∆Σ bl νν ′ ( k , iω n ) = X R X mm ′ P R ∗ νm ( k ) ∆Σ imp mm ′ ( iω n ) P R m ′ ν ′ ( k )(8)where G bl νν ′ ( iω n ) = n [( iω n + µ − ε k ν ) δ νν ′ − ∆Σ bl ( k , iω n )] − o νν ′ (9)∆Σ imp mm ′ ( iω n ) = Σ imp mm ′ ( iω n ) − Σ dc mm ′ . (10) C. Truncating the Bloch basis set, and choice oflocal orbitals
As pointed out above, it is computationally impossi-ble to implement these LDA+DMFT equations withoutrestricting oneself to a finite subspace of N b Bloch func-tions. Those Bloch functions span a certain energy win-dow, corresponding to a subspace W of the total Hilbertspace. Naturally, the local atomic-like orbitals | χ Rk m i will in general have a decomposition involving all Blochbands. The Bloch-transform of these local orbitals read: | χ Rk m i = X ν h Ψ k ν | χ Rk m i| Ψ k ν i . (11)Note that, starting from an orthonormalized set of localorbitals h χ R m k | χ R ′ m ′ k ′ i = δ mm ′ δ RR ′ δ kk ′ , it is easily checkedthat the matrix h Ψ k ν | χ Rk m i is unitary (from the complete-ness of the Bloch basis). Hence the χ m ’s can formally beviewed as Wannier functions associated with the com-plete basis set of all Bloch states. This property no longerholds, however, when the sum in (11) is restricted to thesubset W of Bloch band. Defining: | ˜ χ Rk m i ≡ X ν ∈W h Ψ k ν | χ Rk m i| Ψ k ν i (12)it is seen that the functions | ˜ χ Rk m i are not true Wannierfunctions associated with the subspace W since the trun-cated projection matrix is no longer unitary. However,these functions can be promoted to true Wannier func-tions | w Rk m i by orthonormalizing this set according to: | w Rk m i = X R ′ m ′ S RR ′ mm ′ ( k ) | ˜ χ R ′ k m ′ i , (13)where S RR ′ ( k ) is given by the inverse square root of theoverlap matrix between the Wannier-like orbitals, i.e., O RR ′ mm ′ ( k ) ≡ h ˜ χ Rk m | ˜ χ R ′ k m ′ i = X ν ∈W P R mν ( k ) P R ′ ∗ νm ′ ( k )(14) S RR ′ mm ′ ( k ) ≡ n ( O ( k )) − / o RR ′ mm ′ (15)Naturally, the functions w R m are more extended in spacethan the original atomic-like functions χ R m since they canbe decomposed on a smaller number of Bloch functions,spanning a restricted energy range.In the end, LDA+DMFT is implemented by taking for C the correlated subset generated by the set of functions | w Rk m i . Since those functions have a vanishing overlapwith all Bloch functions which do not belong to the set W , the LDA+DMFT equations can now be put in a com-putationally tractable form, involving only a N b × N b ma-trix inversion within the selected space W . Hence, theequations which are finally implemented read: G R , imp mm ′ ( iω n ) = X k , ( νν ′ ) ∈W ¯ P R mν ( k ) G bl νν ′ ( k , iω n ) ¯ P R ∗ ν ′ m ′ ( k ) , (16)∆Σ bl νν ′ ( k , iω n ) = X R X mm ′ ¯ P R ∗ νm ( k ) ∆Σ imp mm ′ ( iω n ) ¯ P R m ′ ν ′ ( k ) , (17)where ¯ P R mν ( k ) ≡ X R ′ m ′ S RR ′ mm ′ ( k ) P R ′ m ′ ν ( k ) , (18)¯ P R ∗ νm ( k ) ≡ X R ′ m ′ S RR ′ ∗ m ′ m ( k ) P R ′ ∗ m ′ ν ′ ( k ) . (19)It is important to realize that the truncation to a limitedset of Bloch functions was not reached by simply neglect-ing matrix elements between the local orbitals and Blochfunctions outside this set, but rather by constructing anew set of (more extended) local orbitals such that thedesired matrix elements automatically vanish, hence re-defining C accordingly. In this view, the choice of C and of W , although independent in principle, become actuallyinter-related.We also note that it is not compulsory to insist on form-ing true Wannier functions out of the (non-orthogonal)set | ˜ χ Rk m i . It is perfectly legitimate formally to choosethe correlated subspace C as generated by orbitals hav-ing a decomposition in W , but not necessarily unitarilyrelated to Bloch functions spanning W . Although or-thogonality of the χ m ’s is also not compulsory, several(but not all) impurity solvers used within DMFT do re-quire however that the χ ′ m s be orthogonal on a givenatomic site. One possibility, for example, is to orthonor-malize this set on identical unit-cells only, i.e requiringthat | w RT m i and | w R ′ T ′ m ′ i in real space are orthogonal for T = T ′ , but not in neighboring cells T = T ′ . Thisamounts to orthonormalize the | ˜ χ Rk m i set with respect tothe k -summed overlap matrix, instead of the one com-puted at each k -point.In our actual implementations, the wave-functionsspanning the correlated subspace C are obtained by fol-lowing the above orthonormalization procedure, startingfrom atomic-like orbitals χ R m centered on the atomic site R in the primitive unit-cell. These local orbitals are ei-ther all-electron atomic partial waves in the PAW frame-work, or pseudo-atomic wave functions when using theMBPP code. Since, in the present work, we are notdealing with full charge self-consistency including self-energy effects, the matrix elements P R mν ( k ) and the wave-functions | w Rk m i can be computed once and for all at thebeginning of the DMFT cycle . Details on the specificconstruction of the local orbitals used in this article, andthe corresponding calculation of (6) are summarized inappendix A. D. Physical considerations on the choice of thecorrelated subspace C and of the Wannier/Blochspace W Let us now discuss some physical considerations re-garding the choice of the truncated Bloch space W whenusing the LDA+DMFT framework to describe a givenmaterial. Operationally, this means that a certain num-ber of Bloch bands N b (spanning a certain energy win-dow) will be retained when solving Eqs. (16,17). As dis- cussed above, the choice of W also influences the actualdefinition of the correlated subspace C , since we requirethat the orbitals generating C can be expanded upon thebasis functions generating W .Let us consider, to be specific, the case of a transition-metal oxide, such as for example SrVO . This material,which is described in more details in the following section,has a set of three t g bands, well separated from boththe O p and e g bands, and containing nominally one d -electron. We can make two rather extreme choices whendescribing this material with LDA+DMFT:i) Focus only on a very limited set of low-energy Blochbands, such as the three t g bands, and generate W justfrom the three corresponding Bloch functions. In thiscase, we shall have also C = W , and the | w RT m i will beWannier functions unitarily related to these three Blochbands. Since these bands span a narrow energy win-dow, this also means that these Wannier functions willbe rather extended spatially: although centered on Vana-dium atoms, they will have a sizeable contribution onneighboring oxygen atoms as well. This kind of approachhas been emphasized and studied in details in Ref. [2].Of course, it is then out of the scope of such approachesto investigate the indirect effects of correlations on bandsother than the t g ones.ii) Alternatively, one may choose a large energy win-dow to define W , including in particular all Bloch bandscorresponding to O p , V- t g and V- e g . Then, the or-bitals | w RT m i defining the correlated subspace C may bechosen as having a component on Bloch states spanninga much larger energy range. As such, they will be morelocalized spatially, i.e closer to (vanadium) atomic-likeorbitals. When working with such an enlarged space W ,the physics of O p and e g states can also be addressed.One of the goals of the present paper is to present andcompare calculations done with such different choices of W and C . Of course, in a fully first-principle approach,the screened interaction matrix elements should also becalculated (e.g in a GW framework) in a manner whichis consistent with these choices of Hilbert spaces. This isleft for future investigations however, and in the presentwork, these matrix elements will be taken as parameters.Let us note that the local orbitals are constructed inthe present paper from an atomistic point of view. Hencethere is no “entangling problem” as the one encounteredwhen constructing maximally localized Wannier func-tions for strongly hybridized band complexes. Of course,in the present formalism, the Wannier functions obtainedby projection of atomic orbitals onto Bloch states be-longing to a narrow energy window are not maximally-localized in the sense of Ref. 11, but this feature does notbring essential differences in the results, as clear from theresults reported below.Finally, in order to relate the LDA+DMFT results toexperiments performed using e.g, photoemission spec-troscopy, the real-frequency spectral functions must beobtained. This can be done using e.g. a maximum en-tropy treatment of Monte Carlo data, but is importantto understand how to connect the calculated quantitiesto physically observable spectra. The most direct out-put of the LDA+DMFT calculation is the local impurityspectral function, obtained (for orbital m ) as: A imp m ( ω ) ≡ − π ImG impmm ( ω + i0 + ) (20)This also corresponds to matrix elements of the fullGreen’s function of the solid within atomic-like orbitals χ m ’s. This however, is not a quantity which can be eas-ily related to photoemission experiments, since the χ m ’s,when very localized spatially, extend over a large energyrange. When considering a certain energy, contributionsof other electrons with a different orbital character thanthe χ m ’s will contribute significantly to the photoemis-sion signal. Instead, if one is interested in the measuredspectral function in a given energy window, one mustconsider the matrix elements of the full Green’s functionwithin Bloch (or Wannier) functions spanning that en-ergy range, namely: A ν ( ω ) ≡ − π Im X k G bl νν ( k , ω + i0 + ) (21)This quantity can be obtained either by first reconstruct-ing the local self-energy on the real-frequency axis byanalytical continuation of the impurity Green’s func-tion, or by direct analytic continuation of the imaginary-frequency Bloch Green’s function G bl ( k , iω n ). III. APPLICATION: SRVO SrVO is a t g e g metal. It is a good test casefor LDA+DMFT calculations because it is cubic andnon magnetic and also the t g bands are isolated fromboth e g and oxygen p bands – in the LDA bandstruc-ture. Numerous calculations (including LDA+DMFT)and experiments have been done on this compound . From theses studies, itappears that there is a need to include correlations todescribe correctly this compound. This is thus an idealsystem to benchmark this new implementation. A. LDA
For the PAW calculations, semicore states of V andSr are treated in the valence. Valence states for Sr, Vand O thus include respectively 4 s p s , 3 s p s d and2 s p states. PAW matching radius are 1.52 a.u., 1.92a.u. and 2.35 a.u. respectively. The experimental cubiccrystal structure is used (space groupe Pm ¯3 m with lat-tice constant of 7.2605 a.u.). Atomic data are generatedusing ATOMPAW . Calculations are done with thePAW code ABINIT . The density of states and theLDA bandstructure are shown on Fig. 1 and 2. The pro-jection of the density of states on O- p , V- t g and V- e g R Γ X M Γ -8-6-4-202468 ε − ε f ( e V ) FIG. 1: LDA band structure for SrVO . -8 -6 -4 -2 0 2 4 6 8eV0246810 ( e l ec t r on s / e V ) totalO-pV-t V-e g FIG. 2: (Color online) LDA total and projected density ofstates for SrVO . and the character of the bands (see Fig. 3) show thatbands with Oxygen p and with Vanadium t g charactersare indeed isolated from the others. The hybridizationbetween Oxygen and Vanadium orbitals is neverthelessclearly seen. B. LDA+DMFT
In this section, we will present the results of our newscheme based on the projection of Bloch states upon Lo-cal Orbitals (Projection on Local Orbitals (PLO)). Aspreviously emphasized, the number of KS bands used forthe projection have to be chosen in a given range of en-ergy. The extension of the Wannier-like renormalized or-bitals (Eq. 13) will depend on this choice. We will use thefact that the band structure of SrVO is made of isolatedblocs. The choices of W are summarized in table III B.We will compare the results with LDA+DMFT calcula-tions done with Maximally Localized Wannier Functions(MLWF) (see appendix C). FIG. 3: (Color online) LDA band structure for SrVO com-puted in PAW, with ”Fatbands” to show the amplitude of theprojection of each band on a given atomic orbital (O- p , V-t g and V- e g ). N b N imp Energy range (eV)3 3 -1.5 → → → N b used within the PLO scheme of LDA+DMFT.The energy range spanned by these bands is also indicated. N b = 3 corresponds to only the three bands of mainly t g character. N b = 12 corresponds to the 9 Oxygen- p bandsand the 3 Vanadium- t g bands (see Fig. 3). N b = 21 corre-sponds to a large number of bands, including e g bands and 7bands above. N imp is the number of orbitals in the impuritymodel. When used, MLWF will be extracted with the sameparameters. Calculations have been carried out at T=0.1 eV. Adensity-density interaction vertex is used , with J=0.65eV, similarly to another study . The Hubbard parame-ters U=4 eV and U=6 eV have been used in the calcu-lations. The impurity problem of DMFT is solved withHirsch-Fye Quantum Monte Carlo with 128 time slices.The Around Mean Field formulation of the double count-ing is used in these calculations (see appendix B). t g bands in the basis for W ( N b = 3 ). Spectral functions obtained from the PLO and MLWFschemes are reproduced on figure 4.The impurity Green’s function obtained within the twoschemes are in very good agreement: they both show alower Hubbard band at -1.8 eV and a upper Hubbardband at 2.5 eV. This shows that the two schemes contain -4 -2 0 2 4 6 8 (eV) ( e V - ) MLWFWannier from PLO τ (eV -1 )-0.8-0.6-0.4-0.20 G ( τ ) FIG. 4: (Color online) Local impurity spectral function ofSrVO for U=4 eV within LDA+DMFT using the MLWFbasis and the PLO scheme ( N b = 3). Inset: Green’s functionin imaginary time. Note that when using such a small ( t g )energy window, the impurity spectral function and Bloch-resolved spectral function coincide. the same physical content and that PLOs constructedfrom a small energy window give very similar results to maximally localized Wannier functions .
2. Oxygen- p and Vanadium- t g bands in the basis for W ( N b = 12 ). Oxygen p states are now included in the calculation.MLWFs of t g symmetry will thus be more localized, andwill have less O- p character, because the MLWFs of Oxy-gen p symmetry will have mainly the weight on Oxygensatoms. In the PLO scheme also, the Wannier orbitals willinclude more Bloch states, and thus will be closer to tolocalized local orbitals χ m . In other words, because theWannier function of the local impurity problem is nowconstructed from a wider energy windows which includesoxygen p states, the corresponding spectral function willhave a non negligible weight in the energy area of theoxygen- p states. In the LDA case without Hubbard cor-rection, this spectral function would in fact be identicalto the projected t g density of states plotted on Fig 1.In order first to benchmark our calculation using PLOwith respect to calculations using MLWF, we plot on fig-ure 5 the comparison of the impurity Green’s functionin the two cases. The local impurity spectral functionshows a large hybridization band in the area of the oxy-gen p states. We emphasize that this band is mainly nota many-body feature: it is a manifestation of hybridiza-tion with oxygen states, and is readily visible at the pureLDA level. The lower Hubbard band is in fact partly con-tained in this spectral function as a shoulder in the hy-bridization band around -1.5 eV. However, disentanglingthe Hubbard band from the hybridization contributionis somewhat difficult to achieve with the Maximum en- -8 -6 -4 -2 0 2 4 6 8 (eV) ( e V - ) τ (eV -1 )-0.8-0.6-0.4-0.20 G ( τ ) FIG. 5: (Color online) Local impurity spectral function ofSrVO for U=6 eV within LDA+DMFT using the PLOscheme ( N b = 12). Inset: Green’s function in imaginary timefor both the PLO and the MLWF schemes. tropy method. The upper Hubbard band is visible in theimpurity spectral function, as a hump around 3.0 eV.These many-body features (lower and upper Hubbardband) are more clearly revealed, however, by looking atthe ( k -averaged) spectral function of t g Bloch states(21). It is also useful to look at the spectral functionsof the other Bloch states in the basis-set. The summa-tion of these spectral functions over bands belonging tothe same group (e.g Oxygen- p or Vanadium-t g states,ie bands with mostly Oxygen- p or Vanadium-t g char-acter ) enables us to have a clear view on the impactof correlation on LDA bands. The spectral function ofOxygen- p s tates, and Vanadium-t g s tates are plottedon figure 6. We emphasize that these spectral functionsdo not correspond to atomic-like orbitals with Oxygen- p character and Vanadium-t g character (the latter beingthe impurity spectral function plotted on Figure 5).The Hubbard band appears as a hump in the k -averaged Bloch spectral function corresponding to the t g states for U=4 eV. This hump is located between theO- p states and the t g -states quasiparticle peak. TheHubbard band is more clearly resolved for U=6 eV. Inthis case, however, it is hidden inside the Oxygen p band.We shall come back to this point at the end of this sec-tion.The fact that a higher value of U is necessary in thiscase (with respect to N b = 3.) to recover the lower Hub-bard band is consistent with the fact that Wannier func-tions are more localized for N b = 12.
3. Large number of bands in the basis for W ( N b = 21 ) In this case, the e g states are included in the calcu-lation. The impurity model is now solved with all five d − orbitals. The agreement between impurity Green’sfunctions computed in the MLWF and PLO schemes is (a) -10 -8 -6 -4 -2 0 2 4 (eV) ( e V - ) O pV t Total
U=4 eV, J=0.65 eV (b) -10 -8 -6 -4 -2 0 2 4 (eV) ( e V - ) O pV t Total
U=6 eV, J=0.65 eV
FIG. 6: (Color online) Total Bloch spectral function forSrVO , and spectral functions of Oxygen- p , and Vanadium-t g states for U=4 eV (a) and U=6 eV (b) withinLDA+DMFT using the PLO scheme ( N b = 12). These spec-tral functions are not the local orbitals projected spectralfunctions (see text): the t g local impurity spectral functionsare plotted on figure 5. shown on figure 7. Note again, that resolving the Hub-bard band from the impurity Green’s function is diffi-cult because this quantity is dominated by hybridizationeffects with oxygen states. Again, we have to turn toBloch-resolved spectral functions, plotted in Fig. 6 forOxygen- p , Vanadium-t g and Vanadium-e g Bloch states.The results are quite similar to the previous ones with N b = 12, with a lower Hubbard band clearly visible forthe t g states at about -2.0 eV in the spectral function.Note that here we use U = 6 eV. It shows that the basisof Kohn-Sham bands is adapted to the calculation: theconvergence of physical properties as a function of thenumber of bands is rather fast. -8 -6 -4 -2 0 2 4 6 8 (eV) ( e V - ) V-t V-e g τ (eV -1 )-0.8-0.6-0.4-0.20 G ( τ ) FIG. 7: (Color online) Local impurity spectral function ofSrVO for U=6 eV within LDA+DMFT using the PLOscheme ( N b = 21). Inset: Green’s function in imaginary timefor both the PLO and the MLWF schemes. -8 -6 -4 -2 0 2 4 6 8 (eV) ( e V - ) O pV t V e g Total
U=6 eV, J=0.65 eV
FIG. 8: (Color online) Total Bloch spectral function forSrVO , and spectral functions of Oxygen- p , Vanadium-t g and Vanadium-e g s tates for U=6 eV within LDA+DMFT us-ing the PLO scheme ( N b = 21). These spectral functions arenot the local orbitals projected spectral functions (see text):the t g and e g local impurity spectral functions are plottedon figure 7.
4. Discussion
The comparisons made above should not be simplythought of as a convergence study as a function of thesize of the Bloch basis W . Indeed, as N b is increased,we change the spatial extension of the local orbitals (Eq.13) spanning the correlated subspace, so that the DMFTtreatment does not apply to the same objects. Conver-gence studies with fixed χ m ’s could also be performed,but this is not the main scope of this article.Instead, we would like to emphasize some of the physi-cal issues when applying DMFT to different local orbitals -6 -4 -2 0 2 4 6 8(eV)00.511.52 ( e V - ) N b =3N b =12N b =21 FIG. 9: (Color online) Spectral function of t g Bloch states forSrVO within the PLO scheme of LDA+DMFT. For N b = 12and N b = 21: U=6eV. For N b = 3: U=4eV. with different degrees of spatial localization. As we haveseen, for more localized orbitals (larger N b ), we have toincrease the value of the on-site U on the vanadium siteto get consistent results. This is of course expected phys-ically. We see however that the main features of the spec-tral function (see Fig. 9) are almost identical for N b =12and N b =21 with similar values of U . The lower Hubbardband (between -1.7 eV and -2.0 eV), is not far from theposition found experimentaly ( between -1.5eV and -2.0 eV ). The Hubbard band (between 2.5 and3eV) is also in the range of experiments (2.5 eV).In the results that we have obtained, the lower Hub-bard band is rather systematically at an energy in whichoxygen states already give a sizeable contribution in thetotal density of states. Experiments, however, seem tosuggest a somewhat larger separation. We believe thatthis is due to the fact that the relative location of oxy-gen and V- t g state is not accurately obtained at theLDA level, and that a better starting point (such asGW) is required to handle this problem with better ac-curacy. Indeed, we have verified that the intensity of thelower Hubbard band, and especially, of the upper Hub-bard band (as revealed in the Bloch-resolved t g spec-tral function) is very sensitive to the precise value of thedouble-counting correction which is used. This is shownin Fig. 10, in which we have slightly shifted downwardsthe O − p states with respect to the V- t g ones by choos-ing purposedly a smaller value of the double countingcorrection (3.6 eV instead of 6.6 eV). The more isolatedthe t g states are, and the more prominent are correlationeffects within that band (for a given U ).In future work, calculations with extended basis setsshould naturally face the issue of calculating the on-siteCoulomb interaction from first-principles, but also of tak-ing into account Coulomb repulsion terms U pp on theoxygen sites, as well as inter-site repulsions U pd betweenvanadium and oxygen. We note that treating the lat-ter in the Hartree approximation precisely brings in acorrection to the relative position of oxygen and vana-dium states, of the same nature than the double-counting -10 -8 -6 -4 -2 0 2 4 6 8 10 (eV) ( e V - ) AMF double countingweaker DC (=shift of the oxygen bands)
FIG. 10: (Color online) Spectral functions of Oxygen- p , andVanadium-t g states for U=6 eV within LDA+DMFT usingthe PLO scheme ( N b = 12). Two values of the double count-ing shift are used: 6.6 eV (AMF) and 3.6 eV. terms.The general conclusion of this study of SrVO at thisstage is that our formalism is able to describe the mainfeatures of the experimental function within a generalformalism which can take into account all states in thebasis. IV. APPLICATION: β - NiS
The hexagonal form of nickel sulfide ( β -NiS) has at-tracted a lot of interest over the years since it exhibitsa first-order electronic phase transition at about 260K. Whereas the high-temperature phase may be classi-fied as a paramagnetic metal, below the transition β -NiSshows antiferromagnetic order and the resistivity behav-ior corresponds to characteristics of a semi metal ora degenerate semiconductor . Hence the term “anti-ferromagnetic nonmetal” is commonly used for the low-temperature phase. Note that hexagonal NiS is onlymetastable at room temperature, the true stable phaseis given by the millerite structure . The crystal struc-ture (Fig. 11) of β -NiS is of the NiAs-type (space group P /mmc ) with two unit cells in the primitive cell .In this rather simple structure the NiS octahedra shareedges within the ab plane and share faces along the c axis. There is a slight decrease of the cell parameters( δa/a ∼ δc/c ∼ β -NiS, since this wouldinvolve deeper considerations concerning the role of mag-netism above and below the transition. For an overviewof this, still controversially discussed, topic see e.g. thereview by Imada et al. .Here we mainly want to use the high-temperaturephase of β -NiS as an example for a correlated metal withstrongly hybridized Ni(3 d ) and S(3 p ) states, formally inthe charge-transfer regime of the Zaanen-Sawatzky-Allen FIG. 11: (Color online) Projected β -NiS structure with a localNi( e g ) orbital obtained from a MLWF construction using theblock of 16 Ni(3 d )/S(3 p ) bands. The c axis is identical to thevertical axis. classification scheme . In this regard it will become clearthat this compound may not easily be treated in the “tra-ditional” LDA+DMFT scheme via projecting solely ontolow-energy states close to the Fermi level. Thus the goalis to use the here outlined projection technique of in-terfacing LDA with DMFT in order to explore the im-portance of electronic correlations for the local spectralfunction. A. LDA investigation
Electronic structure calculations for β -NiS date backto the original work of Mattheis . Here the DFT-LDAcalculations are performed with a mixed-basis pseudopo-tential (MBPP) code described in appendix A 2. Forthe lattice parameters we use a =3.440 ˚A and c =5.351˚A from Ref. 36. Figure 12 shows the resulting bandstructure and density of states (DOS) from this compu-tation. It is seen that the common block of hybridizedNi(3 d ) and S(3 p ) bands is isolated from the remainingbands. The dominantly S(3 s )-like bands are ∼ , while the bands above starting withmainly Ni(4 s ) character are separated by ∼ ∼ d L peak at -1.5 eV (P ), a second at -4.4eV (P ) and the last at -6.5 eV (P ). Whereas P hasmixed Ni(3 d )/S(3 p ) character, P is a nearly pure S(3 p )peak and does not stem from bonding between Ni and S.In order to obtain more detailed information about theinvolvement of the different orbital sectors, the orbital-resolved DOS is additionally incorporated in Fig. 12b.Furthermore, the fatband representation of the decompo-sition of the Bloch bands is presented in Fig. 13. Becauseof the hexagonal symmetry the Ni(3 d ) multiplet splitsinto two degenerate e g levels, two degenerate E g levels,and one A g level per atom, respectively. As usual, the e g states hybridize more strongly with the S(3 p ) states0than the remaining ( E g , A g ) states. Due to the increas-ingly filled d shell of Ni, the ( E g , A g ) levels are more-over also nearly completely occupied in β -NiS. In num-bers, the LDA filling, including degeneracy, amounts to(2.9;3.9,1.9) for ( e g ; E g , A g ). The fatband representationclearly shows that the hybridization between the differ-ent orbitals are strong and its not quite obvious to singleout a low-energy regime for this compound within thecentral block of bands. Although the bands at the Fermilevel appear to be dominantly of e g character, the corre-spondig orbitals have strong weight at lower energy, too.Despite the strong filling of the ( E g , A g ) states, abandon-ing those orbitals in a minimal hamiltonian description of β -NiS is likely to fail. Fluctuations originating from hole-like states may be important in the end to reach a bet-ter understanding of the complex magnetic behavior ofthis material. Note also in this respect that the effectivebandwidth of the energy range where the e g character ismanifest is roughly twice as large as the correspondingenergy range for ( E g , A g ). Hence this reduced relativebandwidth for ( E g , A g ) within the Ni(3 d ) multiplet willalso have influence on the degree of the orbital-resolvedcorrelation effects. (a) -16-14-12-10-8-6-4-20246 ε - ε F ( e V ) Γ K M Γ A H L A (b) -16 -12 -8 -4 0 4
E-E F (eV) DO S ( / e V ) totale g E g A S(3p)
FIG. 12: (Color online) LDA data for β -NiS. (a) Band struc-ture. (b) DOS. For the local Ni(3 d )/S(3 p )-DOS the cutoffradius was half the nearest-neighbor distance, respectively. FIG. 13: (Color online) LDA fatband decomposition for β -NiS. (red/gray) Ni( E g ) and (blue/dark) Ni( A g ). The localNi(3 d )/S(3 p ) projection cutoff radius was half the nearest-neighbor distance, respectively. B. LDA+DMFT investigation
The influence of correlation effects in β -NiS has beeninvestigated, experimentally and theoretically, by severalauthors . It is generally believed that thiscompound is moderately correlated in view of the stronghybridization of the ∼ / e g states with sulfur. Avalue for the Hubbard U of the order of 4 − . Theseexperiments in the metallic phase are in rough agreementwith the shown LDA DOS below the Fermi level in so faras they also reveal three peaks with comparable relativeintensity as the theoretical set (P , P , P ). However, aneffective single-particle picture appears to be insufficientto understand those peaks, especially when going to thenonmetallic phase . It is theoretically expected that alower Hubbard band, i.e., satellite, originating from theNi(3 d ) states should be located within the Bloch statesdeep in energy (starting around 6 eV) and dominantlycharacterized by S(3 p ). This idea relies on the fact that β -NiS may be viewed as belonging to the charge-transfercategory of transition-metal chalcogenides, yet not be-ing that strongly correlated for the lower Hubbard bandappearing completely below the S(3 p ) states.The strong hybridization between Ni(3 d ) and S(3 p )renders the PLO version of the LDA+DMFT methodmost suitable for this compound. Concerning the valueof U we chose a pragmatic approach and performed thecalculations for two possibly reasonable choices, i.e., U =4eV and U =5 eV. The value of J is certainly less materialsdependent and was fixed to J =0.7 eV. In order to takecare of the double counting, we used the formalism offixing the local total charge (see appendix B). Note thatfor the present crystal structure there are two symmetry-equivalent Ni atoms in the primitive unit cell, i.e., twocorrelated sites R . The local Green’s function and self-energy were thus computed by symmetrizing the site-resolved quantities. For the projection onto local orbitals1we limited the number of bands to N b =16, i.e., all bandsof the central block around the Fermi level are used.This renders the corresponding correlated subspace al-ready rather localized, e.g., the effective Ni(3 d )-like WFsfrom this set of bands are not expected to leak much tothe sulfur sites (see also Fig. 11).By explicitly including the correlation effects, theorbital-resolved fillings in the 3 d shell of Ni do not changerelevantly, thus effects due to changes in orbital popula-tions induced by correlations are not expected to play animportant role for this compound. Figure 14 exhibits theresulting local spectral function for the different values of U at inverse temperature β =10 eV − . It is seen that theinfluence of correlation effects on the metallic spectralfunction are indeed rather subtle for this compound. Alower Hubbard band appears to show up for U > p ) energy regime. It is however -12 -8 -4 0012345 ρ ( ω ) ( e V ) U = 0 -12 -8 -4 0 ω (eV) U = 4 eV -12 -8 -4 0 4 total10 upper bands6 lower bands
U = 5 eV
FIG. 14: (Color online) The LDA+DMFT local spectral func-tion for β -NiS derived from the local Green’s function in Blochbasis at β =10 eV − . The contribution from the upper 10bands of the Ni(3 d )/S(3 p ) block was encoded red, the onefrom the lower 6 bands was encoded blue. This guidance tothe eye should roughly seperate dominant Ni(3 d ) from domi-nant S(3 p ) character. -12 -8 -4 00.00.51.01.5 ρ ( ω ) ( / e V ) e g A E g U = 0 -12 -8 -4 0 ω (eV) U = 4 eV
FIG. 15: (Color online) The LDA+DMFT local impurityspectral function for β -NiS at β =15 eV − . hard to extract this atomic-like excitation from the to-tal spectral function (as in experiment ). The plotof the impurity spectral function in Fig. 15 reveals thatthe ( E g , A g ) orbitals yield indeed effective bands withsmaller bandwidth than the e g orbitals. It seems thatthe contribution to the lower Hubbard band then alsostems more significantly from the former set of orbitals.In the end, the total spectral function computed for U =5 eV shows close resemblance to recent experimen-tal curves obtained from photoemission . Note howeverthat it appears tricky to disentangle the lower Hubbardband from the P peak of in principle pure S(3 p ) content.Our calculation suggests that the correlation effects havealso impact on the lowest S(3 p )-like states. Further stud-ies are necessary to clarify this issue. V. CONCLUSIONS
We have implemented an effective and flexibleLDA+DMFT scheme into two electronic structure for-malisms based on a plane-wave description, i.e., pro-jector augmented-wave and mixed-basis pseudopotential.The orbitals defining the correlated subspace C , in whichmany-body effects are included, are constructed by pro-jecting local atomic-like orbitals onto a restricted set W of Kohn-Sham Bloch states, similar to the procedureadopted in Ref. 1. Orthonormalization of the projectedorbitals yields effective Wannier functions spanning C .These WFs are not unique: they depend on the energywindow covered by the Bloch functions in W (the largerthe window, the more localized the WFs are). Althoughmore sophisticated constructions of explicit WFs willsurely remain an important tool in the LDA+DMFT con-text, we feel that the more straightforward projectiontechnique will render future developments in this areaeasier and more efficient.Using this method, we have investigated SrVO andcompared different implementations involving differentchoices for the set of Bloch states W and for the cor-responding local orbitals spanning C . We have shownthat the basic physical findings for this compound areconsistent with previous LDA+DMFT treatments, inde-pendently of the chosen implementation. However, thepresent study allows for an explicit treatment of ligandstates, which raises several issues which should be thesubject of further studies. One of these issues is a first-principle determination of the on-site Coulomb matrix el-ements in the different choices of local orbitals and basis-set. Our study supports the expected fact that a largervalue of U dd has to be taken when a larger energy range(and more localized orbitals) is considered. Furthermore,a proper treatment of the on-site repulsion on oxygensites U pp , as well as of oxygen-transition metal inter-siterepulsion U pd should certainly be considered. This is-sue is tightly connected to the choice of double-countingcorrection, which, as we have shown, plays a significantrole.2We have also presented a first LDA+DMFT study ofmetallic β − NiS, a charge-transfer compound for whichthe inclusion of ligand states is crucial. More detailedstudies, extending also into the nonmetallic regime, areneeded in order to gain more insight into this “traditionalcompanion” of the more famous NiO compound wheresuch investigations recently took place .In general, the explicit inclusion of the ligand states inthe LDA+DMFT description, based on generic highly-accurate electronic structure codes, will open the doorto new possibilities for the investigation of strong-correlation effects in real materials.
APPENDIX A: CALCULATION OF THEPROJECTION P R mν In this appendix, we describe the details of the calcula-tion of the projection P R mν from KS Bloch states onto lo-cal orbitals in both the PAW and the MBPP method.
1. Projected Augmented-Wave (PAW) method
In the projected augmented-wave formalism, two kindsof atomic functions are used: the pseudowavefunction e ϕ and the true wavefunction ϕ . They are used to recoverthe correct nodal structure of the wavefunction near thenucleus (see Eq. (9) of Bl¨ochl et al ). A direct cal-culation of the projection P R mν ( k ) would thus require tocompute the three terms resulting from this equation,namely ( e p i is the projector n i for angular momenta l i ,and its projection m i ): P R mν ( k ) = h χ R m | e Ψ k ν i + X i h e p i | e Ψ k ν i ( h χ m | ϕ i i − h χ R m | e ϕ i i )(A1)However, atomic d or f wavefunctions are mainly local-ized inside spheres. As a result, we can compute the pro-jection only inside sphere, in the spirit of current LDA+Uimplementation in PAW. In this particular case, thefirst and third terms of Equation A1, cancel with eachothers. This cancellation is only exact for a complete setof projectors in the energy range of the calculation. Thiscompleteness can be easily tested during the construc-tion of projectors and partial waves. The projection thuswrites as a sum over projectors e p n i : P R mν ( k ) = X n i h e p n i | e Ψ k ν ih χ R m | ϕ n i i (A2)We note that integrals in this equation are computedonly inside sphere. It implies that the projection is doneon an unnormalized wavefunction χ . Nevertheless, thelater normalization of the projection P R mν ( k ) will make anormalization of χ redundant. This implementation hasbeen made with the code ABINIT .Additionally, we implemented the proposedLDA+DMFT scheme in a somehow simplified interfac-ing with the PAW-based Vienna Ab Initio Simulation Package (VASP) by using only the PAW projectorsfrom the standard output for P R mν ( k ) . Althoughapproximate, this latter identification yielded (after, ofcourse, proper normalization) very similar results incomparison with the more well-defined explicit codingdescribed above.
2. Mixed-Basis pseudopotential (MBPP) method
The combination of normconserving pseudopotentialswith a mixed basis of plane waves and localized orbitalsin order to represent the pseudo crystal wave function isthe main ingredient of our employed mixed-basis pseu-dopotential code . More concretely, the MBPP bandstructure code uses the following representation for theKS pseudo wave functions: | Ψ k ν i = X G ψ Gk ν | k + G i + X γlm β γl k νm | φ γl k m i , (A3)with h r | k + G i = 1 √ Ω c e i ( k + G ) · r (A4) h r | φ γl k m i = X T e i k · ( T + R γ ) φ γlm ( r − T − R γ ) . (A5)In these eqns., γ denotes the atom in the unit cell, Ω c thevolume of the unit cell, ν the band and l , m are angularmomentum and azimutal quantum number. The planewaves in this basis extend up to a chosen energy cutoff E pw . The analytical form of the local orbitals φ γlm reads φ γlm ( r ′ ) = i l f γl K lm (ˆ r ′ ) , r ′ = r − T − R γ , (A6)whereby K lm descibes a cubic harmonic and the radialfunction f γl is an atomic pseudo wave function modifiedwith a proper decay and cutoff function in order for f γl tovanish beyond some chose radial cutoff r ( γl ) c . Note thatthe local orbitals are orthonormal and by definition donot overlap between neighboring sites, i.e., Z r ( γl ) c dr | f γl ( r ) | r = 1 (A7) h φ γ ′ l ′ m ′ | φ γlm i = δ γ ′ γ δ l ′ l δ m ′ m . (A8)By identifying the local orbitals φ γlm as the objects toproject onto (i.e. the χ m in the outlined LDA+DMFTformalism), the projection P R mν ( k ) may be written withthe help of (A3) and (A8) as P R γ mν ( k ) = h φ γl k m | Ψ k ν i = X G ψ Gk ν h φ γl k m | k + G i + β γl k νm . (A9)3 APPENDIX B: DOUBLE COUNTING
The double-counting shift (used in eq. 5) is neces-sary in the formalism in order to correct the fact thatelectronic correlations are already treated in an aver-age way within LDA. Since there is no unique way toextract this double counting from LDA, two alternativeformulations were used in this work. First the so-calledAround Mean Field (AMF) method, which is usu-ally adequate for metals. Second, a different formalismin which the double counting is such that the electroniccharge computed from the local noninteracting Green’sfunction and the one computed from the interacting im-purity Green’s function is constraint to be identical:Tr G imp ( iω n ) = Tr G (0)imp ( iω n ) . (B1)Here G (0)imp is generally given by G R , imp , (0) mm ′ ( iω n ) = X k X αα ′ h χ Rk m | B k α ih B k α ′ | χ Rk m ′ i ×× n [ iω n + µ − H KS ( k )] − o αα ′ . (B2)Alternatively one may want to use the Weiss field G ( iω n )instead of the local noninteracting Green’s function ineq. (B1). Note that in any case this formalism asks foran additional convergence parameter δµ dc to fulfill (B1).Hence this parameter has to be found iteratively togetherwith the total chemical potential µ in the DMFT cycle.For SrVO3 (with N b = 12), the respective double-counting shifts are 6.6 eV (AMF) and 5.5 eV (fixing localcharge). However these two calculations give very similarresults. APPENDIX C: MAXIMALLY LOCALIZEDWANNIER FUNCTIONS IN PAW
To define maximally-localized WFs, we need the fol-lowing quantity: M ( k , b ) m,n = h Ψ n, k | e − i b . r | Ψ m, k + b i (C1)In the PAW framework, this quantity can be expressedas a function of the pseudo wave function e Ψ, the projec-tors p , the atomic wave function φ and the pseudo atomic wave function e φ . We use the expression for an operator A within the PAW formalism (Eq. (11) of Bl¨ochl et al ).We thus obtain: M ( k , b ) m,n = h e Ψ n, k | e − i b . r | e Ψ m, k + b i + X i,j h e Ψ k ,m | e p i ih e p j | e Ψ k + b ,m i (C2) × (cid:16) h φ i | e − i b . r | φ j i − h e φ i | e − i b . r | e φ j i (cid:17) A similar expression has been used for Ultra-Soft Pseu-dopotentials by Ferretti et al . The first term is com-puted in the Fourier basis, whereas the second is com-puted in the radial grid. We use the expansion of e − i b . r in spherical harmonics and bessel functions in order tocompute its expectation value over atomic wavefunctions.The overlap matrix M ( k , b ) m,n is used to minimize thespread and localize the WFs. This is done with the pub-licly available Wannier90 code .The code is tested on t g MLWF for SrVO3 using onlythe three KS t g bands. The spread of one of the WFsfor a 8 × × . In excellent agreementwith the FLAPW result (6.96 a.u in Ref. 2). ACKNOWLEDGMENTS
We are grateful to O.K. Andersen, F. Aryasetiawan,S. Biermann, T. Miyake, A. Poteryaev, L. Pourovskii,J. Tomczak and M. Torrent for collaborations and dis-cussions related to the subject of this article. We ac-knowledge support from CNRS, Ecole Polytechnique,the Agence Nationale de la Recherche (under con-tract
ETSF ) and the Deutsche Forschungsgemeinschaft(DFG) within the scope of the SFB 668.
Note added:
As the writing of this work was beingcompleted, we became aware of the related work of Ko-rotin et al. (arXiv:0801.3500) also reporting on an imple-mentation of LDA+DMFT in a pseudopotential plane-wave framework. V. I. Anisimov, D. E. Kondakov, A. V. Kozhevnikov, I. A.Nekrasov, Z. V. Pchelkina, J. W. Allen, S.-K. Mo, H.-D.Kim, P. Metcalf, S. Suga, et al., Phys. Rev. B , 125119(2005). F. Lechermann, A. Georges, A. Poteryaev, S. Biermann,M. Posternak, A. Yamasaki, and O. K. Andersen, Phys.Rev. B , 125120 (2006). V. I. Anisimov, A. I. Poteryaev, M. A. Korotin, A. O. Anokhin, and G. Kotliar, J. Phys. Cond. Matter , 7359(1997). A. I. Lichtenstein and M. I. Katsnelson, Phys. Rev. B ,6884 (1998). S. Y. Savrasov, G. Kotliar, and E. Abrahams, Nature ,793 (2001). O. K. Andersen, Phys. Rev. B , 3060 (1975). E. Pavarini, S. Biermann, A. Poteryaev, A. I. Lichtenstein, A. Georges, and O. K. Andersen, Phys. Rev. Lett. ,176403 (2004). O. K. Andersen and T. Saha-Dasgupta, Phys. Rev. B ,16219 (2000). J. des Cloizeaux, Phys. Rev. , 554 (1963). W. Ku, H. Rosner, W. E. Pickett, and R. T. Scalettar,Phys. Rev. Lett. , 167204 (2002). N. Marzari and D. Vanderbilt, Phys. Rev. B , 12847(1997). M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys. , 1039 (1998). A. Fujimori, I. Hase, H. Namatame, Y. Fujishima,Y. Tokura, H. Eisaki, S. Uchida, K. Takegahara, andF. M. F. de Groot, Phys. Rev. Lett. , 1796 (1992). K. Maiti, D. D. Sarma, M. Rozenberg, I. Inoue, H. Makino,O. Goto, M. Pedio, and R. Cimino, Europhys. Lett. , 246(2001). K. Maiti, Ph.D. thesis, IISC, Bangalore (1997). I. H. Inoue, I. Hase, Y. Aiura, A. Fujimori, Y. Haruyama,T. Maruyama, and Y. Nishihara, Phys. Rev. Lett. , 2539(1995). A. Sekiyama, H. Fujiwara, S. Imada, S. Suga, H. Eisaki,S. I. Uchida, K. Takegahara, H. Harima, Y. Saitoh, I. A.Nekrasov, et al., Phys. Rev. Lett. , 156402 (2004). A. Liebsch, Phys. Rev. Lett. , 096401 (2003). I. A. Nekrasov, G. Keller, D. E. Kondakov, , A. V.Kozhevnikov, T. Pruschke, K. Held, D. Vollhardt, and V. I.Anisimov, Phys. Rev. B , 155106 (2005). E. Pavarini, S. Biermann, A. Poteryaev, A. I. Lichtenstein,A. Georges, and O. K. Andersen, Phys. Rev. Lett. ,176403 (2004). T. Yoshida, K. Tanaka, H. Yagi, A. Ino, H. Eisaki, A. Fuji-mori, and Z.-X. Shen, Phys. Rev. Lett. , 146404 (2005). H. Wadati, T. Yoshida, A. Chikamatsu, H. Kumigashira,M. Oshima, H. Eisaki, Z. X. Shen, T. Mizokawa, and A. Fu-jimori, cond-mat/0603642 (2006). I. V. Solovyev, Phys. Rev. B , 155117 (2006). I. A. Nekrasov, K. Held, G. Keller, D. E. Kondakov, T. Pr-uschke, M. Kollar, O. K. Andersen, V. I. Anisimov, andD. Vollhardt, Phys. Rev. B , 155112 (2006). R. Eguchi, T. Kiss, S. Tsuda, T. Shimojima, T. Mizokami,T. Yokoya, A. Chainani, S. Shin, I. H. Inoue, T. Togashi,et al., Phys. Rev. Lett. , 076402 (2006). A. R. Tackett, N. A. W. Holzwarth, and G. E. Matthews,Comp. Phys. Com. (2001). N. A. W. Holzwarth, M. Torrent, and F. Jollet, ATOM-PAW (http://pwpaw.wfu.edu/) (2007). X. Gonze, J.-M. Beuken, R. Caracas, F. Detraux,M. Fuchs, G.-M. Rignanese, L. Sindic, M. Verstraete,G. Zerah, F. Jollet, et al., Comp. Mat. Science , 478(2002). M. Torrent, F. jollet, F. Bottin, G. Zerah, and X. Gonze, Comput. Mater. Sci. (2007 (in press)), doi:10.1016/j.commatsci.2007.07.020. R. Fr´esard and G. Kotliar, Phys. Rev. B , 12909 (1997). K. Maiti, U. Manju, S. Ray, P. Mahadevan, and D. D. S.I. H. Inoue, C. Carbone, cond-mat/0509643 (2005). J. T. Sparks and T. Komoto, Phys. Lett. , 398 (1967). R. M. White and N. F. Mott, Philos. Mag. , 845 (1971). T. Ohtani, K. Kosuge, and S. Kachi, J. Phys. Soc. Jpn. , 1588 (1970). R. Benoit, J. Chem. Phys. , 119 (1955). J. Trahan, R. G. Goodrich, and S. F. Watkins, Phys. Rev.B , 2859 (1970). J. Zaanen, G. A. Sawatzky, and J. W. Allen, Phys. Rev.Lett. , 418 (1985). L. Mattheiss, Phys. Rev. B , 995 (1974). B. Meyer, C. Els¨asser, F. Lechermann, and M. F¨ahnle,
FORTRAN 90 Program for Mixed-Basis-PseudopotentialCalculations for Crystals , Max-Planck-Institut f¨ur Metall-forschung, Stuttgart (unpublished). S. H¨ufner, T. Riesterer, and F. Hulliger, Solid State Comm. , 689 (1985). A. Fujimori, K. Terakura, M. Taniguchi, S. Ogawa,S. Suga, M. Matoba, and S. Anzai, Phys. Rev. B , 3109(1988). V. I. Anisimov, J. Zaanen, and O. K. Andersen, Phys. Rev.B , 943 (1991). M. Nakamura, A. Sekiyama, H. Namatame, H. Kino, andA. Fujimori, Phys. Rev. Lett. , 2891 (1994). M. Usuda and N. Hamada, J. Phys. Chem. , 744 (2000). S. R. Krishnakumar, N. Shanthi, P. Mahadevan, and D. D.Sarma, Phys. Rev. B , 16370 (2000). J. Kuneˇs, V. I. Anisimov, A. V. Lukoyanov, and D. Voll-hardt, Phys. Rev. B , 165115 (2007). J. Kuneˇs, V. I. Anisimov, S. L. Skornyakov, A. V. Lukoy-anov, and D. Vollhardt, Phys. Rev. Lett. , 156404(2007). P. E. Bl¨ochl, Phys. Rev. B , 17953 (1994). S. G. Louie, K. M. Ho, and M. L. Cohen, Phys. Rev. B ,1774 (1979). O. Bengone, M. Alouani, P. Bl¨ochl, and J. Hugel, Phys.Rev. B , 16392 (2000). B. Amadon, F. Jollet, and M. Torrent (2007). G. Kresse and J. Hafner, J. Phys.: Condens. Matter ,8245 (1994). M. T. Czy˙zyk and G. A. Sawatzky, Phys. Rev. B ,14211 (1994). B. B. A Ferretti, A Calzolari and R. D. Felice, Journal ofPhysics: Condensed Matter , 036215 (16pp) (2007).55