Correlated Band Structure of a Transition Metal Oxide ZnO Obtained from a Many-Body Wave Function Theory

0.44 MB·

by  Masayuki Ochi

Correlated Band Structure of a Transition Metal Oxide ZnO Obtained from a
Many-Body Wave Function Theory
Masayuki Ochi1,2, Ryotaro Arita2, and Shinji Tsuneyuki3,4
1Department of Physics, Osaka University, Machikaneyama-cho, Toyonaka, Osaka 560-0043, Japan
2RIKEN Center for Emergent Matter Science (CEMS), Wako, Saitama 351-0198, Japan
3Department of Physics, The University of Tokyo, Hongo, Bunkyo-ku, Tokyo 113-0033, Japan and
4Institute for Solid State Physics, The University of Tokyo, Kashiwa, Chiba 277-8581, Japan
(Dated: January 13, 2017)
Obtaining accurate band structures of correlated solids has been one of the most important and
7
challenging problems in first-principles electronic structure calculation. There have been promis-
1
0 ing recent active developments of wave function theory for condensed matter, but its application
2 to band-structure calculation remains computationally expensive. In this Letter, we report the
first application of the bi-orthogonal transcorrelated (BiTC) method: self-consistent, free from ad-
n
justableparameters, andsystematically improvablemany-bodywavefunction theory,tosolid-state
a
calculations with d electrons: wurtzite ZnO. We find that the BiTC band structure better repro-
J
ducestheexperimentalvalues of thegaps between thebandswith different characters than several
2
otherconventionalmethods. Thisstudypavesthewayforreliablefirst-principlescalculationsofthe
1
properties of strongly correlated materials.

]
l PACSnumbers: 71.10.-w,71.15.-m,71.20.-b
e
-
r
t To reveal fertile and nontrivial physics in condensed called post-Hartree-Fock (post-HF) methods, have also
s
. matter, first-principles electronic-structure calculation beenappliedtocondensedmatterinrecentyears[15–19].
t
a has established itself as an indispensable tool in recent However,their targetsarein mostcaseslimited to solids
m studies. For this purpose, density functional theory with small unit cells, owing to their expensive computa-
- (DFT) [1, 2] has played a leading role and has been ap- tional cost. In addition, the correlated band structure,
d
pliedtovariousmaterials;however,thelimitationsofthis which is quite useful in various kinds of theoreticalanal-
n
theoreticalframeworkhavecometolight. Oneofthema- yses, is not easily obtained in many WFTs. For exam-
o
c jor problems is an inaccurate description of strong elec- ple,calculationofthebandstructureintheframeworkof
[ tron correlations, e.g., in transition metal oxides. The VMCorDMCrequiresalargenumberofsingle-pointcal-

GW method [3–5] is a promising way to ameliorate the culationsofthe excitedstates,whichis a cleardifference
2
v inaccuracy of the band structures and has been applied from a mean-field-like approach such as DFT, whereby
2 to several solids, including d-electron systems. However, the whole band structure is obtained at once.
1 becausetheGW methodisoftenappliedwithoutsatisfy- From this viewpoint, the transcorrelated (TC)
4
ingself-consistency,anontrivialdependenceontheinitial method [20–23] is a fascinating WFT that can be ap-
2
0 DFTcalculationsisintroduced. Ithasalsobeenreported plied to solids with reliable accuracyandmoderate com-
. thattheGW methodsometimesexhibitsseveredifficulty putational cost [24–29]. The TC method adopts the so-
1
in obtaining convergedresults [6], owing to its perturba- called Jastrow ansatz, which is based on a promising
0
7 tive nature. Another possible choice is to construct ef- strategy often adopted in several WFTs such as QMC
1 fective models from DFT that include correlation terms methods to describe strong electron correlations; i.e.,
v: such as Hubbard U and to solve the models using elab- theelectron-electrondistanceisincludedintomany-body
i orated methodologies [7] such as dynamical mean-field wavefunctions. Explicitlycorrelatedelectronic-structure
X
theory [8–10]. However, the exact correspondence be- theory [30] in quantum chemistry also adopts this strat-
r tweentheeffectivemodelsandthefirst-principlesHamil- egy. Infact,the Gutzwiller- andJastrow-correlationfac-
a
tonian is a nontrivial problem. torshaveoftenbeenusedtodescribestrongelectroncor-
Recently, wave function theory (WFT), which had relation, including the Mott physics in systems such as
been mainly applied to molecular systems and es- the Hubbard model [31–34]. It is also important to note
tablished itself as the gold standard in theoretical that, unlike several WFTs, the whole band structure is
chemistry [11], has become a promising alternative to obtained at once by solving a one-body self-consistent-
DFT for accurate descriptions of electron correlation in field (SCF) equation in the TC method, as described
solids [12]. Among the most powerful frameworks in later in this Letter. Moreover, the TC method is de-
WFT are first-principles quantum Monte Carlo (QMC) terministic, i.e., free fromthe statisticalerror,unlike the
methods [13, 14], such as the variational Monte Carlo QMC methods. Accurate calculations for the Hubbard
(VMC), diffusion Monte Carlo (DMC), auxiliary-field model[35,36]andmolecularsystems[37–39]werealsore-
quantumMonteCarlo(AFQMC),andfull-configuration- portedusingtheTCmethodorothertheoriesthathavea
interaction (FCI) QMC methods. Other kinds of WFT, closerelationshipwiththeTCmethod. However,insofar
2
assolid-statecalculationsareconcerned,the TCmethod respectively. As is evident, the HF method can be re-
hassofarbeenappliedonlytoweaklycorrelatedsystems. garded as the TC method with u=0. Owing to this ef-
In this Letter, we present the first application of the fectiveone-bodypicture,itispossibletotreatthemany-
TC method to the band-structure calculation of a d- body correlation with moderate computational cost. In
electron system: wurtzite ZnO. 3d transition metal ox- addition,onecanobtainthe bandstructureofthequasi-
ides have posed theoretical challenges for first-principles particles by using the real part of the eigenvalues of the
band-structurecalculations,asitiswellknownthatpop- ǫ matrix on the right-hand side of Eq. (2). One of the
ular approximations such as the local-density approxi- authorsprovedinpriorwork[23]thatsuchauseoftheǫ
mation (LDA) fail to provide their band structures ac- matrix as quasiparticleenergies is consistentwith Koop-
curately, as we shall see later. We find that the TC mans’theorem. Wenotethatonecansystematicallyim-
method with the bi-orthogonal formulation (the BiTC prove the accuracy of the TC method by utilizing quan-
method) [29, 40] successfully reproduces the experimen- tum chemical methodologies such as the coupled-cluster
tal band structure of ZnO. We also clarify how the Jas- and configuration interaction methods to go beyond a
trowfactorimprovesthefirst-principlesdescriptionofthe single Slater determinant (e.g. Refs. [22, 29, 40]).
correlated electronic states through comparison of the Here, the Jastrow function, u(x,x′), is set to the fol-
band structures and electron densities among the BiTC lowing simple form without adjustable parameters: [13,
and other methods. 24, 42, 43]
ThecentralconceptoftheTCmethodistomakeuseof
tnhiaensiwmiitlhartihtyetJraasntsrfoowrmfaacttioonr oFf=theexmpa(ny-bodyu(Hxa,mxil)t)o,- u(x,x′)= r Ar′ (1−exp(−|r−r′|/Cσ,σ′)), (5)
− i,j i j | − |
P
Ψ=EΨ Φ=EΦ ( =F−1 F), (1)
TC TC whereA= V/(4πN)(N is the numberofvalence elec-
H ⇔H H H
where the correlated wave function is represented as trons in thepsimulation cell, V is the volume of the sim-
Ψ=FΦand TC iscalledtheTCHamiltonian. Here, x ulation cell) and Cσ,σ′ = √2A (spin parallel: σ = σ′),
denotes a pairHof space and spin coordinates: x=(r,σ). √A (spin antiparallel: σ = σ′). The long-range asymp-
6
By adopting the so-called Slater-Jastrow ansatz, Φ be- totic form of this function describes the screening effect
comesaSlaterdeterminantconsistingofone-electronor- of the electron-electron Coulomb interaction [44]. The
bitals, φ(x): Φ =det[φ (x )], and Eq. (1) yields an SCF short-rangebehavioroftheexactJastrowfunctionshould
i j
equation for one-electronorbitals that experience the ef- obey the cusp condition [45–47]. The Jastrow ansatz
fective interaction described with the TC Hamiltonian: adopted here works well for state-of-the-artQMC meth-
ods[13,14]. AlthoughourchoiceoftheJastrowfunction
1
2+v (x ) φ (x ) is rather simple, we shall see that, nevertheless, it works
(cid:18)−2∇1 ext 1 (cid:19) i 1
well not only for weakly correlated systems [24, 25] but
N
+ dx φ∗(x )v (x ,x )det φi(x1)φi(x2) also for the 3d-electron system. Of course, it is possible
Xj=1Z 2 j 2 2body 1 2 (cid:20)φj(x1)φj(x2)(cid:21) tuosiinmgparocvoemaplqicuaatleitdyJoafstthroewmfaancyto-br,odbyutwwaeveafduonpctteidonthbiys
1 N N simple trial wave function to realize moderate computa-
∗ ∗
dx dx φ (x )φ (x )v (x ,x ,x )
− 2 Z 2 3 j 2 k 3 3body 1 2 3 tional cost.
Xj=1Xk=1 In this study [49], we adopted the BiTC method [29,
φi(x1)φi(x2)φi(x3) N 40], in which the left one-electron orbitals, χ(x) in the
×detφj(x1)φj(x2)φj(x3)= ǫijφj(x1), (2) left Slater determinant X = det[χi(xj)], replace the
φk(x1)φk(x2)φk(x3) Xj=1 bra orbitals in the SCF equation (2), while the ket or-
 
bitals remain φ(x). Becausethe TC Hamiltonian is non-
where v (x ), v (x ,x ), and v (x ,x ,x ) are
ext 1 2body 1 2 3body 1 2 3
Hermitian, φ(x) and χ(x) become different. We do not
the external potential including the nucleus-electron in-
showthebandstructurecalculatedusingtheTCmethod
teraction [41] and the two- and three-body effective in-
without the bi-orthogonalextension here, because of the
teractions in the TC Hamiltonian, defined as
large imaginary part of the eigenvalues [57].
2
1 1 Figure 1 presents the band structures of ZnO calcu-
v (x ,x ) + [ 2u(x ,x )
2body 1 2 ≡ r r 2 ∇i 1 2 lated with the LDA, all-electron G0W0 starting from
| 1− 2| Xi=1
LDAusingtheLAPW[58]method,BiTC,andHFmeth-
( u(x ,x ))2+2 u(x ,x ) ], (3)
i 1 2 i 1 2 i ods. The characteristic energy values in these band
− ∇ ∇ ·∇
and structures, as evaluated at the Γ point, are listed in
Table I. Table I also lists calculated values with vari-
v (x ,x ,x ) u(x ,x ) u(x ,x )
3body 1 2 3 1 1 2 1 1 3 ous other methods. By focusing on the gaps between
≡∇ ·∇
+ u(x ,x ) u(x ,x )+ u(x ,x ) u(x ,x ), the bands with different characters as listed in Table I,
2 2 1 2 2 3 3 3 1 3 3 2
∇ ·∇ ∇ ·∇
(4) the BiTC band structure exhibits better accuracy than
3
Bandgap[Error] O-2pbottom[Error] Zn-3daveraged Zn-3dbottom
DFT LDA 0.7 [−2.7] − − − −5.8
HSE03a 2.1 [−1.3] −4.9 [+0.4] − −6.5
GW G0W0(LDA)b 2.4 [−1.0] −5.2 [+0.1] − −6.5
G0W0(HSE03) 3.2a,3.46c [−0.2,+0.06] − − −6.21c −7.2a
G0W0 (GGA+U)c 2.94 [−0.46] −5.6 [−0.3] −6.33 −7.1
(U−J =6eV)
G0W0+Vd (GGA+U)c 3.30 [−0.1] −5.5 [−0.2] −7.45 −8.0
(U−J =6eV,Vd=1.5eV)
QSGWd 3.87 [+0.47] −5.3 [±0] − −7.2
scGW (RPA)e 3.8 [+0.4] − − −6.4 −
scGW (e-h)e 3.2 [−0.2] − − −6.7 −
WFT AFQMCf 3.26(16) [−0.14] − − − −
DMCg 3.8(2) [+0.4] − − − −
BiTC 3.1 [−0.3] −5.1 [+0.2] −9.1 −9.7
HF 11.4 [+8.0] −5.7 [−0.4] −9.3 −9.9
Expt. 3.4h −5.3h,−5.2(3)i −7.5c,−7.5(2)j,−8.5(4)k, −
−8.6(2)k,−8.81(15)i
TABLE I: Some characteristic values in the band structures of ZnO, as calculated by several methods [74]. The bottoms of
theZn-3d and O-2p bandsare evaluated at theΓ point where thevalence-band top is set to zero for each method. The errors
betweenthecalculatedandexperimentalvaluesshowninRef.[62]arepresentedinparentheses. ForLDA,theO-2pbottomand
Zn-3d-averagedlevelsarenotpresentedherebecauseoftheoverlapamongtheO-2pandZn-3dbands. ForG0W0 (HSE03),the
Zn-3d bottom position was read from thedensity of states presented in Ref. [64]. All values are in eV. a Ref. [64]. b Ref. [60].
c Ref. [63]. d Ref. [65]. e Ref. [66]. f Ref. [67]. g Ref. [68]. h Ref. [62]. i Ref. [69]. j Ref. [70]. k Ref. [71].
consistentGW methodssuchastheQSGW method[59].
It is noteworthy that both the GW and BiTC methods
yield such successful results despite being based on con-
ceptually differentformulations. Itis importantthatthe
BiTC method is basedonthe self-consistentformulation
andthusisindependentofDFTcalculations,whereasthe
G0W0 method strongly depends upon the unperturbed
DFTcalculations,asseeninTableI.Weshouldalsonote
that, while some methods employ parameters U and V ,
d
which are difficult to determine in the ab initio way, the
BiTC method does not use such parameters.
Consistencyofthe calculatedbandgapsamongWFTs
that use similar trial wave functions (AFQMC, DMC,
and BiTC) is also remarkable. We again stress that the
wholebandstructureisnoteasilyobtainedinQMCsim-
ulationsandrequiresmanysingle-pointexcited-statecal-
culations. Wecanseethatthebandstructurescalculated
with the BiTC and HF methods are very similar, except
FIG. 1: Calculated band structures using the LDA (solid
lines), all-electron G0W0 [60] (green circles), BiTC (solid forthebandgap. Therefore,themainroleoftheJastrow
factor used in this study on the band structure seems
line), DMC [68] (orange circles), and HF (solid line) meth-
ods. Experimental data taken from Ref. [61] are shown with to be improvement of the size of the band gap through
black dots in the middle figure. Blue broken lines show the the screening effect of the electron–electron interaction,
other experimental data for the positions of the conduction- which is described with the long-ranged asymptotic be-
band minimum, the O-2p bottom [62], and the Zn-3d peak havior of the Jastrow factor. To obtain a more accurate
position [69], which might correspond to the averaged posi-
bandstructure,e.g.,withrespecttothe depthofthe Zn-
tion of theZn-3d bands. The valence-bandmaximum energy
3d bands, more elaborated Jastrow factors, as used in
is set to zero.
QMCstudies [13,14], willbe necessary. Becausethe Zn-
3d bands are almost flat, a key point might be accurate
description of the atomic states, which is an important
many other conventional methods, including the G0W0 issue for future investigation.
method starting fromLDA, and accuracycomparableto For the GW method, it was pointed out that the
that for the most accurate varieties of the GW scheme, shallow Zn-d bands can strengthen the p-d hybridiza-
i.e.,theG0W0 methodstartingfromHSE03andtheself- tion, and thus can result in underestimation of the band
4
solelywithone-bodydegreesoffreedom. Thisisacharac-
teristic advantage of the Slater-Jastrow-type wave func-
tionforthedescriptionoflocalizedelectronicstates,such
as in strongly correlated systems. It is noteworthy that
the atomic calculations of the TC method also exhibit a
similar tendency for localization [23].
Finally, we mention the computational effort required
for the BiTC calculation. Computation takes place on
time scales given by (N2N2N logN ), where N ,
N , and N are thOe nukmbberspwof k-ppowints, occupiekd
b pw
bands, and plane waves, respectively [76]. This is the
same order as that for the HF or hybrid DFT calcula-
tions with a prefactor about 20 to 40 [77]. The BiTC
calculation involves neither the frequency index nor the
FIG.2: Calculatedelectrondensitiesonthebluelineshownin
convergence with respect to the number of conduction
the crystal structure are presented for the LDA, BiTC, and
bands, unlike some perturbative methodologies such as
HF methods with solid black, broken blue, and broken red
the GW method. As can be seen from the fact that
lines, respectively. The crystal structure was depicted using
theVESTA software [75]. the hybrid DFT calculations have now been applied to
various periodic systems, the computational cost of the
BiTC method is reasonable for solid-state calculations.
OneremainingobstacleforwideapplicationoftheBiTC
gap [63, 66]. A similar situation might also be realized methodisthatweusethenorm-conservingpseudopoten-
in the BiTC band structure,whereas the Zn-dbands are tialwithaveryhighcutoffenergytohandlethesemicore
rather deep [72] and so the band gap can be overesti- states at present. However, this problem is not inherent
mated. However, the BiTC band gap is also affected by to the BiTC method and can be overcome in principle
the A parameter in the Jastrow factor, as mentioned in by the development of a pseudopotential formalism such
the above comparison between the BiTC and HF band as the PAW method [78, 79] adaptedto the TC method,
structures. Because the A parameter used in this study whichisanimportantfutureissue[80]. Wealsonotethat
was determined by RPA analysis of the uniform elec- the deformationof the electrondensity nearatomic sites
tron gas,it can cause overscreeningin the insulator [26], shown in Fig. 2 also implies the importance of careful
thereby decreasing the band gap. One possibility is that treatment for the core states. This can also be a com-
thesetwofactorsarecanceledhere,butmoredetailedin- mon problem for QMC calculations using the Jastrow
vestigationonothermaterialsisalsoanimportantfuture correlationfactor.
issue [73]. Toconclude,weapplythe bi-orthogonalversionofthe
Figure2presentstheelectrondensitiescalculatedwith TC method to wurtzite ZnO and find that it well re-
theLDA,HF,andBiTCmethodsonthelineshowninthe produces the experimental band structure. Our study
crystal structure. The electron density obtained by the encourages further investigation of other strongly corre-
BiTC method is defined as n(r) = Re[ N χ∗(r)φ (r)], lated materials using the BiTC method.
i=1 i i
where the condition dr n(r) = NPis satisfied due A part of the calculations was performed using su-
to the bi-orthonormalRization condition χi φj = δij. percomputers at the Supercomputer Center, Institute of
h | i
We can see that the electron densities of these meth- SolidStatePhysics,TheUniversityofTokyo. Thisstudy
ods are almost the same. However, a slight increase of was supported by Grant-in-Aid for young scientists (B)
the electron density at the atomic sites is observed for (No. 15K17724) from the Japan Society for the Promo-
the BiTC method compared with the others. Such a tion of Science, MEXT Element Strategy Initiative to
tendency is consistent with the fact that the strong di- Form Core Research Center, CREST, JST, and Compu-
vergence of the electron-electron Coulomb repulsion is tational Materials Science Initiative, Japan.
alleviated for the effective two-body interaction in the
similarity-transformedTCHamiltonian,becausetheJas-
trowfactorsatisfiesthe cuspcondition. Moreconcretely,
the 2u(x ,x ) and 2u(x ,x ) terms in Eq. (3) yield
1/r∇1 r 1 di2vergenc∇e2with1 a2different sign than the [1] P. Hohenberg and W. Kohn, Phys. Rev. 136, B864
1 2
| − | (1964).
electron-electron Coulomb repulsion. As can be seen in
[2] W.KohnandL.J.Sham,Phys.Rev.140,A1133(1965).
the proof of the cusp condition [45, 46], the true many-
[3] L. Hedin, Phys.Rev.139, A796 (1965).
bodywavefunctionshouldexhibitdeformationdescribed
[4] M. S. Hybertsen and S. G. Louie, Phys. Rev. Lett. 55,
with the two-body degrees of freedom near the electron- 1418 (1985).
electron coalescence point, which cannot be represented [5] M.S.HybertsenandS.G.Louie,Phys.Rev.B34,5390
5
(1986). 425 (2010).
[6] B.-C. Shih, Y. Xue, P. Zhang, M. L. Cohen, and S. G. [38] H. Luo, J. Chem. Phys.133, 154109 (2010).
Louie, Phys. Rev.Lett. 105, 146401 (2010). [39] T. Yanai and T. Shiozaki, J. Chem. Phys. 136, 084107
[7] M. Imadaand T. Miyake, J. Phys.Soc. Jpn. 79, 112001 (2012).
(2010). [40] O. Hino, Y. Tanimura, and S. Ten-no, J. Chem. Phys.
[8] A.Georges,G.Kotliar,W.Krauth,andM.J.Rozenberg, 115, 7865 (2001).
Rev.Mod. Phys. 68, 13 (1996). [41] Note that the external potential vext does not include
[9] G. Kotliar, S. Y. Savrasov, K. Haule, V. S. Oudovenko, the DFT potential unlike some other post-Kohn-Sham
O.Parcollet, andC.A.Marianetti,Rev.Mod.Phys.78, treatments.
865 (2006). [42] D. M. Ceperley, Phys.Rev.B 18, 3126 (1978).
[10] K.Held, Adv.Phys.56, 829 (2007). [43] D.M.CeperleyandB.J.Alder,Phys.Rev.Lett.45,566
[11] A. Szabo and N. S. Ostlund, Modern Quantum Chem- (1980).
istry: IntroductiontoAdvancedElectronicStructureThe- [44] D. Bohm and D. Pines, Phys.Rev.92, 609 (1953).
ory, (Macmillan, New York,1982). [45] T. Kato, Commun. Pure Appl.Math. 10, 151 (1957).
[12] G.H.Booth,A.Gru¨neis,G.Kresse,andA.Alavi,Nature [46] R. T. Pack and W. B. Brown, J. Chem. Phys. 45, 556
(London) 493, 365 (2013). (1966).
[13] W. M. C. Foulkes, L. Mitas, R. J. Needs, and G. Ra- [47] To fully satisfy the cusp condition for the singlet and
jagopal, Rev.Mod. Phys. 73, 33 (2001). triplet spin pairs, one should adopt the operator repre-
[14] J. Kolorenˇc and L. Mitas, Rep. Prog. Phys. 74, 026502 sentation fortheJastrow function[48], whichintroduces
(2011). the troublesome non-terminating series of the effective
[15] M. Marsman, A. Gru¨neis, J. Paier, G. Kresse, J. Chem. interaction in the TC Hamiltonian. Thus, we used this
Phys.130, 184103 (2009). approximated cusp condition in the same manner as in
[16] A.Gru¨neis,M.Marsman,G.Kresse,J.Chem.Phys.133, many QMC studies [13, 14].
074107 (2010). [48] S. Ten-no, J. Chem. Phys. 121, 117 (2004).
[17] J.J.Shepherd,A.Gru¨neis,Phys.Rev.Lett.110,226401 [49] The following HF and (Bi)TC calculations were per-
(2013). formed using the tc++ code [24, 25], which makes
[18] J. J. Shepherd, T. M. Henderson, G. E. Scuseria, Phys. use of the LDA [50] orbitals obtained using the tapp
Rev.Lett. 112, 133002 (2014). code [51, 52] as an initial guess of the HF and (Bi)TC
[19] J. J. Shepherd, T. M. Henderson, G. E. Scuseria, J. orbitals. Singularities of the electron-electron Coulomb
Chem. Phys.140, 124102 (2014). repulsion term and the Jastrow function in the k-space
[20] S. F. Boys and N. C. Handy, Proc. R. Soc. London Ser. were handled using a method proposed by Gygi and
A 309, 209 (1969); ibid. 310, 43 (1969); ibid. 310, 63 Baldereschi [53] using an auxiliary function of the same
(1969); ibid. 311, 309 (1969). form as Ref. [54]. For calculations of wurtzite ZnO, we
[21] N.C. Handy,Mol. Phys. 21, 817 (1971). useda3×3×3k-mesh,acutoffenergyof361Ry,and150
[22] S. Ten-no, Chem. Phys. Lett. 330, 169 (2000); ibid. 175 bands to expand the subspace for block-Davidson diag-
(2000). onalization in the HF and (Bi)TC methods (n = 150
[23] N. Umezawa and S. Tsuneyuki, J. Chem. Phys. 119, with the notation used in Ref. [27]). Lattice parame-
10015 (2003). ters were extracted from the experiment: a = 3.2496
[24] R. Sakuma and S. Tsuneyuki, J. Phys. Soc. Jpn. 75, ˚A, c = 5.2042 ˚A, and u = 0.3819 [55]. The same non-
103705 (2006). local norm-conserving pseudopotentials [56], where the
[25] M. Ochi, K. Sodeyama, R. Sakuma, and S. Tsuneyuki, semicore 3s and 3p orbitals of zinc atoms were treated
J. Chem. Phys. 136, 094108 (2012). as valence electrons, were used for our LDA, HF, and
[26] M. Ochi, K. Sodeyama, and S. Tsuneyuki, J. Chem. (Bi)TC calculations. [57]
Phys.140, 074112 (2014). [50] J. P. Perdew and A. Zunger, Phys. Rev. B 23, 5048
[27] M. Ochi, Y. Yamamoto, R. Arita, and S. Tsuneyuki, J. (1981).
Chem. Phys.144, 104109 (2016). [51] J. Yamauchi,M. Tsukada, S.Watanabe, and O. Sugino,
[28] M. Ochi and S. Tsuneyuki, J. Chem. Theory Comput. Phys. Rev.B 54, 5586 (1996).
10, 4098 (2014). [52] O. Sugino and A. Oshiyama, Phys. Rev. Lett. 68, 1858
[29] M. Ochi and S. Tsuneyuki, Chem. Phys. Lett. 621, 177 (1992).
(2015). [53] F. Gygi and A. Baldereschi, Phys. Rev. B 34, 4405
[30] W. Kutzelnigg, Theor. Chim. Acta 68, 445 (1985). (1986).
[31] H. Yokoyamaand H. Shiba, J. Phys. Soc. Jpn. 56, 1490 [54] S. Massidda, M. Posternak, and A. Baldereschi, Phys.
(1987). Rev. B 48, 5058 (1993).
[32] M. Capello, F. Becca, M. Fabrizio, S. Sorella, and E. [55] H. Karzel, W. Potzel, M. K¨offerlein, W. Schiessl, M.
Tosatti, Phys.Rev.Lett. 94, 026406 (2005). Steiner, U. Hiller, G. M. Kalvius, D. W. Mitchell, T. P.
[33] H.Yokoyama,T.Miyagawa,andM.Ogata,J.Phys.Soc. Das, P. Blaha, K. Schwarz, and M. P. Pasternak, Phys.
Jpn.80, 084607 (2011). Rev. B 53, 11425 (1996).
[34] T. Miyagawa and H. Yokoyama, J. Phys. Soc. Jpn. 80, [56] N. Troullier and J. L. Martins, Phys. Rev. B 43, 1993
084705 (2011). (1991).
[35] S.Tsuneyuki,Prog.Theor.Phys.Suppl.176,134(2008). [57] See SupplementalMaterial for more information.
[36] J. M. Wahlen-Strothman, C. A. Jim´enez-Hoyos, T. [58] O. K. Andersen,Phys.Rev. B 12, 3060 (1975).
M. Henderson, and G. E. Scuseria, Phys. Rev. B 91, [59] S. V. Faleev, M. van Schilfgaarde, and T. Kotani, Phys.
041114(R) (2015). Rev. Lett.93, 126406 (2004).
[37] H.Luo,W. Hackbusch,and H.-J. Flad, Mol. Phys.108, [60] M. Usuda, N. Hamada, T. Kotani, and M. van Schilf-
6
gaarde, Phys.Rev.B 66, 125101 (2002).
[61] G.ZwickerandK.Jacobi, Solid StateCommun.54,701
(1985).
[62] Numerical Data and Functional Relationships in Sci-
ence and Technology, edited by K. H. Hellwege and
O. Madelung, Landolt-Bo¨rnstein, New Series, Group III
(Springer,Berlin, 1982), Vols. 17a and 22a.
[63] L. Y. Lim, S. Lany, Y. J. Chang, E. Rotenberg, A.
Zunger, and M. F. Toney, Phys. Rev. B 86, 235113
(2012).
[64] A.R.H.Preston,B.J.Ruck,L.F.J.Piper,A.DeMasi,
K.E.Smith,A.Schleife,F.Fuchs,F.Bechstedt,J.Chai,
and S. M. Durbin,Phys.Rev. B 78, 155114 (2008).
[65] T. Kotani, M. van Schilfgaarde, and S. V. Faleev, Phys.
Rev.B 76, 165106 (2007).
[66] M. Shishkin, M. Marsman, and G. Kresse, Phys. Rev.
Lett.99, 246403 (2007).
[67] F. Ma, S. Zhang, and H. Krakauer, New J. Phys. 15,
093017 (2013).
[68] J. Yu, L. K. Wagner, and E. Ertekin, J. Chem. Phys.
143, 224707 (2015).
[69] L. Ley, R. A. Pollak, F. R. McFeely, S. P. Kowalczyk,
and D. A.Shirley,Phys. Rev.B 9, 600 (1974).
[70] R.A.Powell,W.E.Spicer,andJ.C.McMenamin,Phys.
Rev.Lett. 27, 97 (1971).
[71] C. J. Vesely, R. L. Hengehold, and D. W. Langer, Phys.
Rev.B 5, 2296 (1972).
[72] Also for weakly correlated systems, we have observed
that localized states (e.g., the 3s states of bulk silicon)
tend tobecome deeperin their band structures[25–27].
[73] In our observation, the similar trend is also realized for
other semiconductors and insulators with sp electrons,
which include both the ionic and covalent crystals.
[74] Refs.[60]and[63]usedthelatticeparametersthatarea
bit differentfrom ours. However,such differenceyieldsa
differenceinthebandgapandthed-bandpositionforless
than0.1eVinLDAcalculations. Therefore, wecan con-
cludethatthedifferenceofthelatticeparametersamong
several references has a minor effect on the comparison
presented in Table 1.
[75] K.Momma andF. Izumi,J.Appl.Crystallogr. 44, 1272
(2011).
[76] This scaling is valid only when a single Slater determi-
nant is used (i.e. multi-determinant methodologies are
not used) and the Jastrow function is represented as
u(x,x′)=Ppwp(x)vp;σ,σ′(|r−r′|)wp(x′)(i.e.notagen-
eral two-body function).
[77] TheratioofcomputationaltimerequiredfortheTCand
HF calculations is about 15 to 20 [25]. As for the BiTC
method, this ratio becomes about 20 to 40 owing to the
nonequivalencyof thebra and ket orbitals.
[78] P.E. Bl¨ochl, Phys.Rev.B 50, 17953 (1994).
[79] G.KresseandD.Joubert,Phys.Rev.B59,1758(1999).
[80] Our present implementation requires about 20,000 core
hours for solving the BiTC-SCF equation, which is de-
termined by a very large cutoff energy of plane waves.
If we can make the cutoff energy an order of magnitude
lower, which is a usual demand for other methods, by
development of efficient core treatment, this cost will be
reduced also byan order of magnitude.
7
SUPPLEMENTAL MATERIAL
CONVERGENCE WITH RESPECT TO THE NUMBER OF k-POINTS
We haveperformedband-structurecalculationsusing a 2 2 2 k-mesh(32 atomsin the unit cell) for the HFand
BiTC methods to confirm the convergence with respect to×the×number of k-points. For the BiTC calculations, the
bandgapandtherelativeenergyoftheZn-d-bandbottomtothevalence-bandtopare2.92and 9.73eV,respectively,
for a 2 2 2 k-mesh, whereas they are 3.09 and 9.68 eV, respectively, for a 3 3 3 k-m−esh (108 atoms in the
× × − × ×
unit cell). For the HF calculations, the corresponding values are 11.53 and 10.07 eV, respectively, for a 2 2 2
k-mesh, whereasthey are 11.42and 9.90eV,respectively,for a 3 3 3k-mesh. The difference betweenthe×val×ues
calculated using the two k-meshes is−less than 0.2 eV for all of thes×e va×lues.
In our previous study on the TC method [25], we saw that the calculated band gap E (N ) with N k-points
gap k k
obeys an approximate relation:
E (N ) E (N = )+C/N , (6)
gap k gap k k
≃ ∞
where C is a constant. Using this approximate relation for the band energies, the residual finite-size error in our
calculations using a 3 3 3 k-mesh is estimated to be less than 0.1 eV. Therefore,we can say that the convergence
with respect to the nu×mbe×r of k-points is sufficiently achieved. Although it is possible that the approximate relation
that is used above does not hold for this case [68], we can expect that the convergence will be sufficient for our
discussion because the finite-size erroris expected to be significantly reduced in a 3 3 3 k-mesh from its value for
× ×
a 2 2 2 one.
× ×
IMAGINARY PART OF THE TC EIGENVALUES WITHOUT THE BI-ORTHOGONAL FORMULATION
We do not show the band structure calculated using the TC method without the bi-orthogonal extension in the
main text. This is because we found that the imaginary part of the eigenvalues of the Zn-3d bands obtained by
solving the one-body SCF equation became as large as 0.1 eV in the TC method, which is more than 1,000 times
larger than that in the BiTC method. Because the non-Hermiticity in the (Bi)TC formulation originates from the
similarity-transformation, the presumed equivalency between the TC solution and the true many-body eigenstate is
suggested to break when the imaginary parts of the eigenvalues become too large. The difference between the band
structuresobtainedwiththeTCandBiTCmethodswassmallforweaklycorrelatedsystemsexceptLiFthatisrather
localized-electron system [29].
ELECTRON-DENSITY DIFFERENCES AMONG THE LDA, HF, AND BITC METHODS
The differences in the electron densities calculated with the LDA, HF, and BiTC methods are shown in Fig. S1.
In Fig. S1, we can verify the tendencies explained in the main text, i.e., that the BiTC electron density tends to be
localizedneartheatomicsites. Inaddition,wecanseeanincreaseoftheBiTCelectrondensityalongtheZn-Obonds
near the oxygen atoms.
8
FIG. 1: The differences between the HF and LDA electron densities and between the BiTC and LDA electron densities are
presented in panels (b) and (c), respectively. The electron densities on the plane shown in panel (a) are depicted here.