Model-mapped RPA for Determining the Effective Coulomb Interaction

Model-mapped RPA for Determining the Effective Coulomb Interaction

Hirofumi Sakakibara    Seung Woo Jang    Hiori Kino    Myung Joon Han    Kazuhiko Kuroki    Takao Kotani Department of applied mathematics and physics, Tottori university, Tottori 680-8552, Japan Computational Condensed Matter Physics Laboratory, RIKEN, Wako, Saitama 351-0198, Japan Department of Physics, Korea Advanced Institute of Science and Technology (KAIST), Daejeon 305-701, Korea Research Organization of Science and Technology, Ritsumeikan university, Kusatsu, Shiga 525-8577, Japan National Institute for Materials Science (NIMS), Sengen 1-2-1, Tsukuba, Ibaraki 305-0047, Japan Department of Physics, Osaka University, Machikaneyama-Cho, Toyonaka, Osaka 560-0043, Japan
July 2, 2019

We present a new method to obtain a model Hamiltonian from first-principles calculations. The effective interaction contained in the model is determined on the basis of random phase approximation (RPA). In contrast to previous methods such as projected RPA and constrained RPA (cRPA), the new method named “model-mapped RPA” takes into account the long-range part of the polarization effect to determine the effective interaction in the model. After discussing the problems of cRPA, we present the formulation of the model-mapped RPA, together with a numerical test for the single-band Hubbard model of HgBaCuO.

71.20.-b, 71.20.Be, 71.30.+h
preprint: APS/123-QED

I Introduction

Recently, we often treat low-temperature physical phenomena in correlated materials by a two-step procedure, that is, deriving a model Hamiltonian from a first-principles calculation and solve the model Hamiltonian Koretsune and Hotta (2014); Tsumuraya et al. (2015); Kinza and Honerkamp (2015); Tsutsui et al. (2015). Thus, the procedure contains two key points:

  • How to derive a model Hamiltonian.

  • How to solve the model Hamiltonian.

For step (B), we can use various many-body calculation techniques to solve the model Hamiltonian. The techniques for step (B) can be dynamical mean-field theory (DMFT)Georges et al. (1996), quantum Monte Carlo (QMC)Ceperley et al. (1977) methods, and so on. The two-step procedure is generally applicable to strongly correlated systems such as high- superconductors and metal–insulator transitions.

In this paper, we focus on step (A), that is, how to derive a model Hamiltonian from first-principles calculations, especially the effective interaction contained in the Hamiltonian (we neglect phonons here). If step (A) is well established and combined with a reasonable technique in step (B), we can even evaluate the transition temperature of superconductivity for given crystal structures without introducing parameters by hand Kent et al. (2008); Sakakibara et al. (2010, 2012); Onari and Kontani (2012); Watanabe et al. (2013); Sakakibara et al. (2014); Misawa and Imada (2014). This means we can use the two-step procedure for material informatics combined with databases of crystal structures. In future, we may find new high- superconductors among thousands of possible candidates Arita et al. (2007); Hansmann et al. (2009) in the two-step procedure.

Let us present an overview of step (A). We have various first-principles methods to determine one-body Hamiltonian . These method are the local density approximation (LDA), the quasiparticle self-consistent GW method (QSGW)Kotani and van Schilfgaarde (2007); Kotani (2014); Friedrich et al. (2014); Klimeš et al. (2014), and so on. The one-body Hamiltonian describes an independent particle picture. The static screened Coulomb interaction can be calculated in the random-phase approximation (RPA). From , we can construct a set of the atomic-like localized orbitals which describe the low-energy bands. The orbitals can be constructed, for example, by the method of the maximally localized Wannier functions Marzari and Vanderbilt (1997). The orbitals span a model Hilbert space . The choice of  is not unique and has ambiguity. If necessary, we should use a larger  to reduce the ambiguity. However, step (B) requires a sufficiently small  to ensure tractability within current computational resources. The one-body part of the model Hamiltonian in can be determined by the projection of into . As for the effective interaction, we cannot simply project into . In advance, we have to remove the screening effect expected within the model. This is necessary to avoid double counting of the screening effect. This idea was first introduced in the projected RPA (pRPA) method by Kotani Kotani (2000). This was followed by the constrained RPA (cRPA) by Aryasetiawan et al. Aryasetiawan et al. (2004). cRPA has been applied to several compounds to analyze strongly correlated systems Nakamura et al. (2008, 2010, 2012); Nomura et al. (2012, 2015); Tsuchiizu et al. (2015). Miyake implemented a Wannier-based modified cRPA in the ecalj package Miyake et al. (2009). Şaşıoǧlu, Friedrich, and Blügel also proposed a modified cRPA applicable to cases with entangled bands Şaşıoǧlu et al. (2011). However, cRPA contain theoretical problems as we discuss in Sec. II.

The reliability of the model Hamiltonian obtained in step (A) is determined by the reliability of the first-principles calculation. Most popular calculations are in LDA. However, LDA often gives an unreliable independent particle picture, especially for transition metal oxides and -electron materials. A well-known problem is the underestimation of the band gaps. Furthermore, there are problems with the bandwidth, the positions of the 3 bands and 4 bands relative to the oxygen bands, and so on. In such cases, we need to use advanced methods such as hybrid functional methods Kim et al. (2010) or QSGW Kotani and van Schilfgaarde (2007); Kotani (2014); Friedrich et al. (2014); Klimeš et al. (2014). One of the advantages of QSGW is that and are determined simultaneously in a self-consistent manner. QSGW has even been applied to metallic ground states. For example, Han et al. recently applied QSGW to LaNiO/LaAlOHan et al. (2014), Jang et al. applied QSGW to high- materials Jang et al. (2015), and Ryee et al. applied QSGW to SrRuO and SrRuO Ryee et al. (2016). To handle such metallic systems, QSGW is more reliable than the hybrid functional methods Heyd et al. (2003). Furthermore, Deguchi, Sato, Kino, and Kotani have recently shown that a QSGW-based hybrid method can systematically give a good description for a wide range of materials Deguchi et al. (2016). These calculations were performed by the first-principles calculation package ecalj eca (), which is based on a mixed-basis full-potential method, the linearized augmented plane wave, and muffin-tin orbital method (the PMT method) Kotani et al. (2015); Kotani and Kino (2013); Kotani and van Schilfgaarde (2010). It is freely available from github eca ().

In this paper, we propose a new method named model-mapped RPA (mRPA). This is based on an assumption of the existence of a model Hamiltonian that explains the low-energy physical properties of materials. This assumption is standard in the field of model calculations. For example, we may assume that low-energy physical properties can be quantitatively understood by a Hubbard model. Then, the role of mRPA is to determine the interaction parameters in the model.

After we review cRPA and point out its problems in Sec. II, we give a formulation of mRPA in Sec. III. Then we show how it works for a test case of single-band Hubbard model for the high- superconductor HgBaCuO in Sec. IV, followed by a summary.

Ii cRPA and its problems

In the first-principles calculations, the screened Coulomb interaction in RPA is given by


where and are the Coulomb interaction and the non-interacting proper polarization, respectively. (written as below) denotes the inverse of matrix . consists of a product of two Green functions . We can represent the quantities and expanded in an improved version of the mixed product basis (MPB). The MPB was originally introduced by Kotani in Ref. Kotani and van Schilfgaarde, 2002. Then, the MPB was improved by Friedrich, Blügel, and Schindlmayr Friedrich et al. (2010). We usually use the improved MPB.

Let us recall the idea of the so-called cRPA. We first choose a model space  spanned by a basis set of atomic-like localized orbitals, , where is the index of the primitive cell and is the index used to specify an orbital in the cell. In the following, we use the notations, and


as in the manner of Ref. Lichtenstein and Katsnelson, 1998. The eigenfunctions in  are calculated from the Hamiltonian , that is, is the same as but restricted within the space . Then, we have the non-interacting proper polarization function in . Let us consider the RPA-screened Coulomb interaction in . It can be written as


where is the (not yet determined) effective interaction between quasi-particles in . cRPA determines by assuming


with Eq.  (3). From Eqs. (1), (3), and (4), we have


that is, we can calculate from and . Then, we calculate the on-site interaction for the model as


If necessary, we can calculate any elements of the interaction given as ). However, in the usual cRPA, we only calculate the set of parameters used by the model. In summary, for the choice of a localized basis set , we determine a set of interactions of the model in cRPA, where are -dependent.

The cRPA, which appears to be reasonable, however, contains the following three problems.

  • Range truncation problem

    is satisfied only when we take all possible elements of in cRPA. However, practically adopted models consider a limited number of . Note that given in Eq.  (5) is inevitably long-range. This is because the strong screening effects such as metallic screening contained in are removed from the total polarization . This problem is well illustrated when only the on-site is used. In this case, because we use only the on-site part of evaluated from the right-hand side of Eq.  (5), given in Eq.  (3) cannot satisfy the condition .

    Schüler et al. Schüler et al. (2013) proposed a method to solve an extended Hubbard model with non local interaction. However, the method is not applicable to the long-range interaction without modification. Hansmann et al. Hansmann et al. (2013) calculated the long-range behavior of the effective interaction by cRPA, and presented a method to solve a model Hamiltonian taking into account the long-range interaction. In contrast, we can handle the same problem using a model Hamiltonian with the short-range interaction given by mRPA. This is because mRPA downfolds the long-range interaction into the short-range interaction as described in Sec. III.

  • Positive definiteness and causality problem

    in the denominator of Eq.  (5) should be positive definite at . If this is not satisfied, we obtain unphysical results for having eigenvalues larger than the bare interactions . Furthermore, the imaginary (anti-hermitian) part of ) at any should be positive definite so as to satisfy causality. In the original idea of cRPA, does not necessarily satisfy this condition for the case of entangled bands. For this reason, Kotani avoided the idea of cRPA and proposed pRPA, which satisfies the above conditions Kotani (2000). Recently, two other procedures satisfying the conditions have been proposed with a modification of the definition of in the cRPA. One given by Miyake, Aryasetiawan, and Imada, neglects the off-diagonal elements between  and residual space in the one-body Hamiltonian [see Fig. 1 and Eq. (8) in Ref. Miyake et al., 2009]. Thus, the condition is clearly satisfied. The other is given by Şaşıoǧlu, Friedrich, and Blügel, where a projection procedure of eigenfunctions to  is used to satisfy the condition Şaşıoǧlu et al. (2011). Strictly speaking, neither procedure should be identified as cRPA, since the key advantage of cRPA, Eq.  (4), is no longer satisfied.

    Note the generality of the causality problem. For example, in the GW+DMFT formulation Biermann et al. (2003); Sun and Kotliar (2004) as an extension of the LDA+DMFT Lechermann et al. (2006); Pavarini et al. (2004); Lichtenstein et al. (2001), the on-site part of the GW self-energy is simply substituted with the DMFT self-energy. Then, we may have a causality problem if we have a GW self-energy whose imaginary part is larger than that of the DMFT self-energy.

  • energy window problem
    In Table 1, we have calculated the static (at ) part of and for the paramagnetic Ni, where we use the cRPA method given by Şaşıoǧlu, Friedrich, and Blügel Şaşıoǧlu (2014). We considered two cases for the energy window; the narrower one is eV and the wider is eV. In contrast to the small difference in for the different energy window, we see a large difference in . The value of 3.56 eV is in good agreement with that in Ref. Şaşıoǧlu et al., 2011. As shown in Fig. 5 in Ref. Şaşıoǧlu, 2014, they used such a wide energy window. The value of 2.25 eV for the narrower window is significantly different from 3.56 eV. This difference is because of the difference in , which describes the polarization of the 3 electrons. In the case of wider windows, we remove more polarization, resulting in larger values of . This results in an inevitable ambiguity in the cRPA because we have almost the same energy bands (the same eigenvalue dispersions in the Brillouin zone) for both windows. In addition, we have no definite criteria for choosing a certain energy window.
    We expect that a similar ambiguity also exists in other versions of cRPA. Miyake, Aryasetiawan, and Imada, successfully obtained flat low-energy behaviors, as shown in Fig. 3 of Ref. Miyake et al., 2009, similar to that in Fig. 1 of Ref. Kotani, 2000. However, the procedure of neglecting the off-diagonal elements (equivalent to how to choose the ) is ambiguous. From Fig. 3 in Ref. Miyake et al., 2009, we guess that the ambiguity of in their method can be 1 eV [from the degree of freedom in the choice of the energy window and , we may have various possible extrapolations of to ].

-8 1 eV -10 10 eV
[eV] 1.19 1.40
[eV] 2.25 3.56
Table 1: Static screened Coulomb interactions and in cRPA method Şaşıoǧlu et al. (2011) for two different outer energy windows, -81 eV, and -10 10 eV for paramagnetic Ni. For both energy windows, we have almost the same 3 bands for the model space . We use 12 12 12 points in the Brillouin zone in the tetrahedron method Kotani and van Schilfgaarde (2007); Kotani (2014). There is a large difference between the two values. See text.

Although the new method, the mRPA, formulated in Sec. III can remedy these problems, we need to pay attention to the inevitable limitations of model Hamiltonians including no long-range interactions. Recall that plasmons (charge fluctuations) do not satisfy the Goldstone’s theorem because of the behavior of the Coulomb interaction. Such model Hamiltonians cannot describe this correctly. The long-range limit of longitudinal spin fluctuations, as well. Model Hamiltonians can only be justified when these problems are irrelevant.

Iii Formulation of the mRPA

Let us assume that a model Hamiltonian in the model space  can describe low-energy excitations very well. Here we formulate mRPA, which determines the parameters included in . is given as


where is the one-body Hamiltonian, obtained from a first-principles method such as QSGW. is the spin-independent effective interaction specified by a set of parameters (here we do not consider the -dependence of these parameters). The terms are those in Eq. (1) in Ref. Lichtenstein and Katsnelson, 1998. is the one-body counter term so that the effect of is canceled out when we apply the first-principles method to the model described by Ikeda et al. (2010). In the , used elements are given by a set of a finite number of parameters . We usually allow only the short-range terms; for example, we only allow the on-site terms in the case of the Hubbard model.

Let us explain how to determine (or equivalently) in mRPA. If we apply RPA to the model Hamiltonian , we have the screened Coulomb interaction of the model as


where we use the proper non-interacting polarization calculated from . In mRPA, we only consider the case at in Eqs. (8)–(11).

Note the difference between Eq.  (3)(cRPA) and Eq.  (8)(mRPA). In Eq.  (3), inevitably become long-range as , while is short-range such as on-site only in Hubbard model. That is, and are non-zero just on the limited number of discrete index set of  in Eq.  (8).

For the theoretical correspondence, we require to satisfy


in mRPA in order to determine . Here, is the quantity calculated from in the first-principles method using Eq.  (1). It is not possible to satisfy Eq.  (9) for all the matrix elements of ; we satisfy a subset of Eq.  (9) corresponding to the degree of freedom of used in of Eq.  (7). Thus, the subset of Eq.  (9) can determine uniquely. Then, we can determine from Eq.  (8) so as to satisfy





By definition, mRPA satisfies Eq.  (9) exactly, where is expressed in terms of and in Eq.  (8). This is in contrast to the case of cRPA, which can not usually satisfy Eq.  (4) since cRPA usually discards the off-site part of . Thus, we are free from the problem (i) in Sec. II in the case of mRPA.

We may have cases that satisfying Eq.  (9) cannot be found. This is because has the upper limit ( the bandwidth of ); for , we have can be seen from Eq.  (8). Thus, we cannot determine for very large . This can be clearly seen in Fig. 1 as explained later. In such cases, we need to use a larger . This is not an intrinsic problem of mRPA but a problem associated with choosing of too small .

The causality problem (ii) in Sec. II does not arise, since -dependence of is meaningless in mRPA; we use the condition given by Eq.  (9) only at . In our opinion, we rather have to use a larger  for better results, instead of taking the -dependence into account in theoretical treatments. If we take the -dependence of the effective interaction correctly, we inevitably have to treat a quantum Langevin equation with electron thermal bath. Such a treatment is far beyond our current numerical techniques because it requires an enormous computational effort.

Let us consider problem (iii) in Sec. II in the case of mRPA. In mRPA, is only determined from the energy bands of . The choice of the Wannier functions (the choice of ) can slightly change . This yields the slight ambiguity of via Eq.  (10). This is inevitable as long as we derive a model from first-principles calculations. In contrast, cRPA has further ambiguity in the polarization of in Eq.  (5) owing to the ambiguity of the choice of Wannier functions as we have shown in Table 1.

Iv Numerical Test for a Single Band Hubbard Model of HgBaCuO

Here we present a test calculation to see how mRPA works in comparison with cRPA. We take a single-band Hubbard model for stoichiometric HgBaCuO. We treat two cases where is determined by LDA or by QSGW. The space  is chosen by a procedure based on maximally localized Wannier functions Souza et al. (2001). The term in Eq.  (7) is irrelevant in the single-band case since it gives a constant potential shift. As we use the single-band Hubbard model, we obtain and as scalars.

The two curves in Fig. 1 show as functions of given by Eq.  (8), where we use calculated from by QSGW or by LDA. As a function of , these curves are initially linear near and saturate toward . The difference between the two curves is due to the size of corresponding to the size of the bandwidth Jang et al. (2015). The two horizontal lines show the values of (0.67 and 0.85 eV) calculated by the first-principles RPA method as shown in Table 2. Using the condition Eq.  (9), we can determine for QSGW and for LDA as illustrated in Fig. 1.

The obtained values of are shown in Table 2, together with the cRPA values obtained by the method in Ref. Şaşıoǧlu et al., 2011. The values of are eV larger than those of . This is because we use the on-site interaction only in the present model. If we take into account off-site interactions, will be reduced. In other words, mRPA downfolds the off-site interactions into the on-site interaction.

[eV] [eV] [eV]
QSGW 0.67 5.2 3.7
LDA 0.85 3.9 2.0
Table 2: Calculated values of (mRPA) and (cRPA) for a single-band model of HgBaCuO , together with . Note that we only consider the values at . We use the tetrahedron method in Ref. Kotani and van Schilfgaarde, 2007 for the evaluation of and , where we use 8 8 4 points in the Brillouin zone, as was used in Ref. Werner et al., 2015. The values of in cRPA Şaşıoǧlu et al. (2011) are the same as those presented in our previous paper Jang et al. (2016). is determined by mRPA as illustrated in Fig. 1.
Figure 1: (Color online) Calculated as a function of in Eq.  (8) for a single-band model of HgBaCuO. The Green line represents QSGW, the read line represents LDA. The values of in mRPA are obtained at the intersections between the curves of and the horizontal lines of .

In Fig. 1, we see that the determined values of are sensitive to the values of and . This is because the calculated values of are close to the upper limit of RPA, at . The derivatives at are rather small, 0.060 for LDA and 0.019 for QSGW. This sensitivity may indicate the lack of suitability (or limitation) of the single-band Hubbard models for HgBaCuO. If we use a larger , we will be able to avoid such cases of .

V Summary

We have presented mRPA to determine model Hamiltonians based on first-principles calculations. mRPA is formulated starting from the assumption of the existence of a model Hamiltonian that explains the low-energy physical properties of materials. Then we determine the effective interactions contained in the Hamiltonian by matching the first-principles RPA calculations and the RPA calculations using the model Hamiltonian. mRPA is free from the theoretical problems in cRPA, which are discussed in Sec. II. Thus, mRPA is less ambiguous and logically clearer than cRPA. Through the model Hamiltonian obtained by mRPA, we will be able to predict the critical temperatures of superconductors based on first-principles calculations.

We appreciate discussions with Drs. Friedlich, Şaşıoǧlu, Miyake, and Arita. H.S. appreciates fruitful discussions with Yunoki, Shirakawa, Seki, and Shinaoka. This work was supported by JSPS KAKENHI (Grant-in-Aid for Young Scientists B, Grant No. 16J21175), and was partly supported by the Advanced Low Carbon Technology Research and Development Program (ALCA) of Japan Science and Technology Agency (JST). S.W.J. and M.J.H were supported by the Basic Science Research Program through NRF (2014R1A1A2057202). The computing resource were supported by KISTI (KSC-2015-C3-042) , the supercomputing system Great-Wave (HOKUSAI) of RIKEN, the supercomputing system of the ISSP, and the Computing System for Research in Kyushu University.


Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description