A simple derivation of the exact quasiparticle theory and its extension to arbitrary initial excited eigenstates

TL;DR

This paper presents a straightforward derivation of the exact quasiparticle theory, extending its applicability to arbitrary initial excited eigenstates, and validates it using the GW approximation on atoms and molecules.

Contribution

It provides a simple derivation of the quasiparticle theory and extends its validity to excited states, demonstrating its applicability beyond ground states.

Findings

QP energies and wave functions can be self-consistently determined.

The theory applies to arbitrary initial excited eigenstates.

Validation with GW approximation confirms the theory's accuracy.

Abstract

The quasiparticle (QP) energies, which are minus of the energies required by removing or produced by adding one electron from/to the system, corresponding to the photoemission or inverse photoemission (PE/IPE) spectra, are determined together with the QP wave functions, which are not orthonormal and even not linearly independent but somewhat similar to the normal spin orbitals in the theory of the configuration interaction, by self-consistently solving the QP equation coupled with the equation for the self-energy. The electron density, kinetic and all interaction energies can be calculated using the QP wave functions. We prove in a simple way that the PE/IPE spectroscopy and therefore this QP theory can be applied to an arbitrary initial excited eigenstate. In this proof, we show that the energy-dependence of the self-energy is not an essential difficulty, and the QP picture holds exactly if there is no relaxation mechanism in the system. The validity of the present theory for some initial excited eigenstates is tested using the one-shot GW approximation for several atoms and molecules.

Tables

Table 1. Table 1: G 0 ​ W 0 subscript 𝐺 0 subscript 𝑊 0 G_{0}W_{0} QP energies of neutral (A), cationic (A + ), anionic (A - ), and photoabsorbed (A ∗ ) states of lithium, aluminum, and beryllium atoms in units of eV. The right atom/ion has one more electron than the left atom/ion; and the left and the right have N ​ ( − 1 ) 𝑁 1 N(-1) and N ​ ( + 1 ) 𝑁 1 N(+1) electron configurations, respectively ( L 2 ​ S + 1 superscript 𝐿 2 𝑆 1 {}^{2S+1}L states are also listed); In the same row, the QP energies of the forward, → → \rightarrow , electron absorption process and the reverse, ← ← \leftarrow , electron release process should coincide with each other. The values inside the parentheses indicate the usual G 0 ​ W 0 subscript 𝐺 0 subscript 𝑊 0 G_{0}W_{0} QP energy calculated for the neutral atom at the ground state. EC in the superscript in the Reference values means the value calculated by the EOM coupled cluster singles and doubles (EOM-CCSD) with aug-cc-pV5Z basis set using Gaussian 09 Gaussian . All the other reference values a-g are experimental data Ref5 ; Ref6 ; Ref7 ; Ref8 ; Ref9 ; Ref11 .

Atom	$N\left(-1\right)$ Electron	${}_{}{}^{2S+1}L$		Atom	$N\left(+1\right)$ Electron	${}_{}{}^{2S+1}L$	QP energies (eV)
/ Ion	Configuration	State		/ Ion	Configuration	State	$\to$	$\leftarrow$	Reference
Li+	${\left(1s\right)}^2$	1S	$\leftrightarrow$	Li	${\left(1s\right)}^2\left(2s\right)$	2S	5.1	(5.6)	5.39a
Li+	${\left(1s\right)}^2$	1S	$\leftrightarrow$	Li∗	${\left(1s\right)}^2\left(2p\right)$	2P	3.0	3.2	3.54b
Li+∗	$\left(1s\right)\left(2s\right)$	1S	$\leftrightarrow$	Li	${\left(1s\right)}^2\left(2s\right)$	2S	63.9	(63.4)	64.46b
Li+∗	$\left(1s\right)\left(2p\right)$	1P	$\leftrightarrow$	Li∗	${\left(1s\right)}^2\left(2p\right)$	2P	65.0	66.0	65.76b
Li	${\left(1s\right)}^2\left(2s\right)$	2S	$\leftrightarrow$	Li-	${\left(1s\right)}^2{\left(2s\right)}^2$	1S	(0.6)	0.3	0.61c
Al+	${\left(3s\right)}^2$	1S	$\leftrightarrow$	Al	${\left(3s\right)}^2\left(3p\right)$	2S	5.7	(5.6)	5.98d
Al+	${\left(3s\right)}^2$	1S	$\leftrightarrow$	Al∗	${\left(3s\right)}^2\left(4s\right)$	2P	2.6	3.1	2.84e
Al+∗	$\left(3s\right)\left(3p\right)$	1S	$\leftrightarrow$	Al	${\left(3s\right)}^2\left(3p\right)$	2P	12.6	(10.2)	13.40e
Al	${\left(3s\right)}^2\left(3p\right)$	2P	$\leftrightarrow$	Al-	${\left(3s\right)}^2{\left(3p\right)}^2$	3S	(0.3)	-	0.43f
Al	${\left(3s\right)}^2\left(3p\right)$	2P	$\leftrightarrow$	Al-∗	${\left(3s\right)}^2\left(3p\right)\left(4s\right)$	3P	(-0.1)	0.5	0.2EC
Al∗	${\left(3s\right)}^2\left(4s\right)$	2S	$\leftrightarrow$	Al-∗	${\left(3s\right)}^2\left(3p\right)\left(4s\right)$	3P	3.5	2.7	3.3EC
Be+	$\left(2s\right)$	2S	$\leftrightarrow$	Be	${\left(2s\right)}^2$	1S	9.3	(9.3)	9.32f
Be+	$\left(2s\right)$	2S	$\leftrightarrow$	Be∗	$\left(2s\right)\left(2p\right)$	1P	4.2	5.4	4.40g
Be+	$\left(2s\right)$	2S	$\leftrightarrow$	Be∗	$\left(2s\right)\left(2p\right)$	3P	6.4	6.2	6.59g

Table 2. Table 2: G 0 ​ W 0 subscript 𝐺 0 subscript 𝑊 0 G_{0}W_{0} QP energies of neutral (M), cationic (M + ), and photoabsorbed (M ∗ ) states of nitrogen, oxygen, and lithium molecules in units of eV. The right molecule has one more electron than the left molecule; and the left and the right have N − 1 𝑁 1 N-1 and N 𝑁 N electron configurations, respectively (electronic states are also listed); In the same row, the QP energies of the forward, → → \rightarrow , electron absorption process and the reverse, ← ← \leftarrow , electron release process should coincide with each other and also with the reference values Molecule1 ; Molecule2 ; Molecule3 ; Molecule4 . The values inside the parentheses indicate the usual G 0 ​ W 0 subscript 𝐺 0 subscript 𝑊 0 G_{0}W_{0} QP energy calculated for the neutral atom at the ground state.

Molecule	State		Molecule	State	QP energies (eV)
					$\to$	$\leftarrow$	Reference
N${}_{}{}^{+}{}_2{}^{}$	${}_{}{}^2\Sigma _g^{+}$	$\leftrightarrow$	N2	${}_{}{}^1\Sigma _g^{+}$	15.69	(15.40)	15.58a
N${}_{}{}^{+}{}_2{}^{}$	${}_{}{}^2\Sigma _g^{+}$	$\leftrightarrow$	N${}_{}{}^{\ast }{}_2{}^{}$	${}_{}{}^1\Pi _g^{}$	7.28	7.21	6.99b
O${}_{}{}^{+}{}_2{}^{}$	${}_{}{}^2\Pi _g^{}$	$\leftrightarrow$	O2	${}_{}{}^3\Sigma _g^{-}$	12.16	(12.35)	12.30c
O${}_{}{}^{+}{}_2{}^{}$	${}_{}{}^2\Pi _g^{}$	$\leftrightarrow$	O${}_{}{}^{\ast }{}_2{}^{}$	${}_{}{}^3\Pi _g^{}$	3.80	3.73	4.17d
Li${}_{}{}^{+}{}_2{}^{}$	${}_{}{}^2\Sigma _g^{+}$	$\leftrightarrow$	Li2	${}_{}{}^1\Sigma _g^{-}$	5.30	(5.32)	5.11e
Li${}_{}{}^{+}{}_2{}^{}$	${}_{}{}^2\Sigma _g^{+}$	$\leftrightarrow$	Li${}_{}{}^{\ast }{}_2{}^{}$	${}_{}{}^1\Sigma _u^{-}$	3.74	3.95	3.35f

Equations

$$E_{\mu }^{M-1}-E_{\gamma }^M=h\nu -K$$

$$E_{\gamma }^M-E_{\nu }^{M+1}=h\nu -K.$$

$${\epsilon }_{\mu }=E_{\gamma }^M-E_{\mu }^{M-1},$$

$${\epsilon }_{\nu }=E_{\nu }^{M+1}-E_{\gamma }^M.$$

$${\varphi }_{\mu }(\mathbf{r},s)$$

$${\varphi }_{\nu }^{\ast }(\mathbf{r},s)$$

$${\rho }_s(\mathbf{r})$$

$$=\sum _{\mu }^{\mathrm{occ}}⟨{\Psi }_{\gamma }^M|{\psi }_s^{†}(\mathbf{r})|{\Psi }_{\mu }^{M-1}⟩⟨{\Psi }_{\mu }^{M-1}|{\psi }_s(\mathbf{r})|{\Psi }_{\gamma }^M⟩$$

$$=\sum _{\mu }^{\mathrm{occ}}{\left|{\varphi }_{\mu }(\mathbf{r},s)\right|}^2,$$

$$\sum _{\mu }^{\mathrm{occ}}|{\Psi }_{\mu }^{M-1}⟩⟨{\Psi }_{\mu }^{M-1}|=1$$

$$H^{(1)}=\hat{T}+\hat{v}=\sum _s\int {\psi }_s^{†}(\mathbf{r})h_s^{(1)}(\mathbf{r}){\psi }_s(\mathbf{r})d\mathbf{r},h_s^{(1)}(\mathbf{r})=-\frac{{\hslash }^2}{2m}{\nabla }^2+v(\mathbf{r})$$

$$⟨\hat{T}⟩$$

$$=-\frac{{\hslash }^2}{2m}\sum _{\mu }^{\mathrm{occ}}\sum _s\int {\varphi }_{\mu }^{\ast }(\mathbf{r},s){\nabla }^2{\varphi }_{\mu }(\mathbf{r},s)d\mathbf{r},$$

$$⟨\hat{v}⟩$$

$$=\sum _{\mu }^{\mathrm{occ}}\sum _s\int v(\mathbf{r}){\left|{\varphi }_{\mu }(\mathbf{r},s)\right|}^{2}d\mathbf{r}.$$

$${\rho }_s({\mathbf{r}}^{\prime },\mathbf{r})=⟨{\Psi }_{\gamma }^M|{\psi }_s^{†}({\mathbf{r}}^{\prime }){\psi }_s(\mathbf{r})|{\Psi }_{\gamma }^M⟩=\sum _{\mu }^{\mathrm{occ}}{\varphi }_{\mu }^{\ast }({\mathbf{r}}^{\prime },s){\varphi }_{\mu }(\mathbf{r},s).$$

$$⟨\hat{T}⟩$$

$$⟨\hat{v}⟩$$

$$\sum _{\mu }^{\mathrm{occ}}{\varphi }_{\mu }(\mathbf{r},s){\varphi }_{\mu }^{\ast }({\mathbf{r}}^{\prime },s^{\prime })+\sum _{\nu }^{\mathrm{emp}}{\varphi }_{\nu }(\mathbf{r},s){\varphi }_{\nu }^{\ast }({\mathbf{r}}^{\prime },s^{\prime })$$

$$=\sum _{\lambda }^{\mathrm{all}}{\varphi }_{\lambda }(\mathbf{r},s){\varphi }_{\lambda }^{\ast }({\mathbf{r}}^{\prime },s^{\prime })=\delta (\mathbf{r}-{\mathbf{r}}^{\prime }){\delta }_{s{s}^{\prime }}.$$

$$\left\{{\psi }_s(\mathbf{r}),{\psi }_{s^{\prime }}^{†}({\mathbf{r}}^{\prime })\right\}=\delta (\mathbf{r}-{\mathbf{r}}^{\prime }){\delta }_{s{s}^{\prime }},$$

$$\sum _{\nu }^{\mathrm{emp}}|{\Psi }_{\nu }^{M+1}⟩⟨{\Psi }_{\nu }^{M+1}|=1$$

$$S_{\lambda {\lambda }^{\prime }}$$

$$=\sum _{s{s}^{\prime }}\int {\varphi }_{\lambda }^{\ast }(\mathbf{r},s)\delta (\mathbf{r}-{\mathbf{r}}^{\prime }){\delta }_{s{s}^{\prime }}{\varphi }_{{\lambda }^{\prime }}({\mathbf{r}}^{\prime },s^{\prime })d\mathbf{r}d{\mathbf{r}}^{\prime }$$

$$=\sum _{s{s}^{\prime }}\int {\varphi }_{\lambda }^{\ast }(\mathbf{r},s)\sum _{{\lambda }^{\prime \prime }}^{\mathrm{all}}{\varphi }_{{\lambda }^{\prime \prime }}(\mathbf{r},s){\varphi }_{{\lambda }^{\prime \prime }}^{\ast }({\mathbf{r}}^{\prime },s^{\prime }){\varphi }_{{\lambda }^{\prime }}({\mathbf{r}}^{\prime },s^{\prime })d\mathbf{r}d{\mathbf{r}}^{\prime }$$

$$=\sum _{{\lambda }^{\prime \prime }}^{\mathrm{all}}\sum _s\int {\varphi }_{\lambda }^{\ast }(\mathbf{r},s){\varphi }_{{\lambda }^{\prime \prime }}(\mathbf{r},s)d\mathbf{r}\sum _{s^{\prime }}\int {\varphi }_{{\lambda }^{\prime \prime }}^{\ast }({\mathbf{r}}^{\prime },s^{\prime }){\varphi }_{{\lambda }^{\prime }}({\mathbf{r}}^{\prime },s^{\prime })d{\mathbf{r}}^{\prime }$$

$$=\sum _{{\lambda }^{\prime \prime }}^{\mathrm{all}}{S}_{\lambda {\lambda }^{\prime \prime }}{S}_{{\lambda }^{\prime \prime }\lambda },$$

$${\tilde{\varphi }}_{\alpha }(\mathbf{r},s)=\sum _{\lambda }^{\mathrm{all}}{\varphi }_{\lambda }(\mathbf{r},s)U_{\lambda \alpha },$$

$${\tilde{S}}_{\alpha {\alpha }^{\prime }}=\int {\tilde{\varphi }}_{\alpha }^{\ast }(\mathbf{r},s){\tilde{\varphi }}_{{\alpha }^{\prime }}(\mathbf{r},s)d\mathbf{r}={\delta }_{\alpha {\alpha }^{\prime }}.$$

$$\sum _{{\alpha }^{\prime \prime }}^{\mathrm{restrict}}{\tilde{S}}_{\alpha {\alpha }^{\prime \prime }}{\tilde{S}}_{{\alpha }^{\prime \prime }\alpha }$$

$$=\sum _{s{s}^{\prime }}\int {\tilde{\varphi }}_{\alpha }^{\ast }(\mathbf{r},s)\sum _{{\alpha }^{\prime \prime }}^{\mathrm{restrict}}{\tilde{\varphi }}_{{\alpha }^{\prime \prime }}(\mathbf{r},s){\tilde{\varphi }}_{{\alpha }^{\prime \prime }}^{\ast }({\mathbf{r}}^{\prime },s^{\prime }){\tilde{\varphi }}_{{\alpha }^{\prime }}({\mathbf{r}}^{\prime },s^{\prime })d\mathbf{r}d{\mathbf{r}}^{\prime },$$

$$\sum _{\lambda }{\tilde{\varphi }}_{\lambda }(\mathbf{r},s){\tilde{\varphi }}_{\lambda }^{\ast }({\mathbf{r}}^{\prime },s^{\prime })=\delta (\mathbf{r}-{\mathbf{r}}^{\prime }){\delta }_{s{s}^{\prime }}$$

$$\displaystyle S^{\prime}=U^{\prime\dagger}SU^{\prime}=\left(\begin{array}[]{cc}\mbox{\boldmath$\mu$}&\;0\\ 0&\;0\end{array}\right),$$

$$\displaystyle\widetilde{\phi}_{\alpha}(\bm{r},s)=\frac{1}{\sqrt{\mbox{\boldmath$\mu$}_{\alpha\alpha}}}\sum_{\lambda}^{\rm all}\phi_{\lambda}(\bm{r},s)U^{\prime}_{\lambda\alpha}$$

$$\displaystyle\phi_{\lambda}(\bm{r},s)=\sqrt{S_{\lambda\lambda}}\sum_{\alpha}^{\rm restrict}\sqrt{\mbox{\boldmath$\mu$}_{\alpha\alpha}}\widetilde{\phi}_{\alpha}(\bm{r},s)U^{\prime\dagger}_{\alpha\lambda}.$$

$${\varphi }_{\lambda }(\mathbf{r},s)=\sum _{\alpha }^{\mathrm{restrict}}{\tilde{\varphi }}_{\alpha }(\mathbf{r},s)U_{\alpha \lambda }^{†}.$$

$$\displaystyle\sqrt{\mbox{\boldmath$\mu$}_{\alpha\alpha}}U^{\prime\dagger}_{\alpha\lambda}\sqrt{S_{\lambda\lambda}}=U^{\dagger}_{\alpha\lambda}.$$

$$\displaystyle\sum_{\alpha}^{\rm restrict}U^{\prime}_{\lambda\alpha}\,\mbox{\boldmath$\mu$}_{\alpha\alpha}U^{\prime\dagger}_{\alpha\lambda^{\prime}}=\left(1/S_{\lambda\lambda}\,\right)\delta_{\lambda\lambda^{\prime}}$$

Full text

I Introduction

There are variety of excited states of atoms, molecules and any kind of materials, for example, irradiated by laser-, electron-, and ion-beams, and they show variety of fascinating behaviors. However, lack of a first-principles theory, which can treat arbitrary electronic excited states of materials, has strongly prevented the progress of theoretical studies of such excited states of materials. In 1964, Hohenberg and Kohn HohenbergKohn invented density functional theory (DFT) on the basis of the variational principle in quantum mechanics. Next year, adopting DFT in the local density approximation (LDA), Kohn and Sham KohnSham proposed to use the so-called Kohn–Sham (KS) equation, which has a simple form of a single particle Schrödinger equation with the Hartree and exchange-correlation potentials. If one could solve the KS equation for continuous fractional occupation numbers $f_{\mu}$ ($f_{\nu}$) for certain occupied (empty) KS levels $\mu$ ($\nu$), one might obtain the quasiparticle (QP) energies, $E^{N}_{f_{\mu}=1}-E^{N-1}_{f_{\mu}=0}=\int_{0}^{1}\varepsilon_{\mu}(f_{\mu})df_{\mu}$ and $E^{N+1}_{f_{\nu}=1}-E^{N}_{f_{\nu}=0}=\int_{0}^{1}\varepsilon_{\nu}(f_{\nu})df_{\nu}$, which should rigorously correspond to the experimental photoemission and inverse photoemission (PE/IPE) spectra, thanks to Janak’s theorem, $\varepsilon_{\lambda}(f_{\lambda})=\partial E/\partial f_{\lambda}$ for $\lambda=\mu$ or $\nu$, by integration Janak . Moreover, in 1985, Almbladh and von Barth Almbladh demonstrated that the correct treatment of the exchange-correlation energy produces the KS eigenvalue for the HOMO (highest-occupied molecular orbital) level that is identical to the minus of the experimental, i.e., exact, ionization potential (IP). However, even if one could perform such calculations, the resulting information is just the usual PE/IPE energy or spectra of the system at the ground state only.

In 1984, in order to treat the non-steady states and the electronic excited states of materials, the time-dependent density functional theory (TDDFT) was invented by Runge and Gross RungeGross . Although this is a very important theory that may treat arbitrary excited states of materials, there is a serious problem that no systematic way is given to determine the exchange-correlation kernel $\hat{\mu}_{xc}$, which is a functional of a whole temporal history of the electron density $\{\rho(\bm{r},t)\}$ and the initial wave function $\Psi(\{\bm{r}_{i}\},t_{0})$ at $t=t_{0}$. So, one usually adopts the adiabatic LDA (A-LDA), in which the whole temporal history is completely ignored in the time domain and the simplest LDA is used in the space domain in the functional. The reliability of the A-LDA has been discussed by many authors Dobson ; Vignale ; Kohn ; Hessler ; D'Amico ; LevyNagy ; Ohno . Recently, it has been also emphasized to use accurate exchange-correlation functional in the study of excited states CrawfordUrangaRubio .

Of course some of the electronic excited states may be treated with highly accurate quantum chemistry (QC) methods Sabo ; CC such as the configuration interaction (CI), complete active space (CAS), Møller–Plesset (MP), and coupled cluster (CC) methods. However, these wave function-based methods beyond the Hartree–Fock approximation (HFA) are computationally very heavy and again restricted to the singly, doubly, or triply (perhaps at most quadruply) excited neutral/charged states only.

On the other hand, the Green’s function method on the basis of the many-body perturbation theory (MBPT) can directly describe the PE/IPE energy spectra, which are usually called the quasiparticle (QP) energies. The number of the papers using this method such as the $GW$ approximation has increased recent years. One difficulty of using this method is that the MBPT and the Green’s function are not easy to understand for the researchers who are not familiar to the field theory. Therefore, the Green’s function is not at all popular in the fields of quantum chemistry and materials science. Moreover, there is a serious problem related to the energy-dependence of the self-energy, which has been a long standing issue and the QP picture has sometimes been considered as an approximation. And again, the main issue here is that the Green’s function has been always defined on the neutral ground state only, while the initial electronic state in the PE/IPE spectroscopy can be any of the neutral or charged excited states.

The purpose of this paper is to present a simple derivation of the QP theory, which is applicable to any of the initial electronic excited eigenstates $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ of the $M$-electron system ($M$ can be different from the total number of protons $N$ in the system) without invoking the complicated field theoretical MBPT. In the course of this derivation, we will show that the energy-dependence of the self-energy is not an essential difficulty, and the QP picture holds rigorously if there is no relaxation mechanism such as electron-photon and electron-phonon interactions or cascade damping in the continuum energy bands. (Even when there is such a relaxation mechanism in the system, the QP picture should still hold within the lifetime Abrikosov ; Lifshitz , or equivalently within the mean free path.)

Moreover, it becomes clear that the QP wave functions $\phi_{\lambda}(\bm{r},s)$, which are defined as the overlap between the initial state $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ operated by a creation or annihilation operator ($\psi^{\dagger}_{s}(\bm{r})$ or $\psi_{s}(\bm{r})$) and the final state $\langle\,{\Psi_{\lambda}^{M\pm 1}}\,|\,$ in a PE/IPE process, satisfy the completeness condition but are generally non-orthonormal and linearly dependent. (A relation to the normal spin orbitals in the CI theory Sabo ; Lowdin1 will be mentioned in Section II.) Although the QP theory looks like a one-electron approximation, it is not an approximation at all. It should be emphasized that the QP energies and the QP wave functions rigorously include all the necessary many-body information and formally exact in a sense that these quantities yield exact electron density as well as the exact kinetic, local potential, and electron-electron interaction energies. The present theory is exact for any initial excited eigenstate.

The rest of this paper is organized as follows. In Section II, the QP energies and the QP wave functions are defined in connection to the PE/IPE spectroscopy, and the equations satisfied by them are derived. The exact expressions for the electron density and the expectation values for the kinetic, local potential, and interaction energies are derived. In Section III, the same equations are rederived using the Green’s functions. In Section IV, a concrete method to calculate the self-energy is derived by invoking Brillouin–Winger’s perturbation theory. In Section V, several test calculations of isolated atoms, ions and molecules using the one-shot $GW$ approximation demonstrate the validity of the present theory. Section VI presents several important remarks: The virial theorem, the Ward–Takahashi identity, which is equivalent to the local charge conservation (the continuity equation), and the macroscopic conservation laws hold, and the Luttinger–Ward functional exists for the total energy, for the excited eigenstate. Finally the DFT universal functional of the density is extended to the excited eigenstate, for which the occupation numbers of all KS levels become its additional dependence, and a variational principle still holds in a restricted Fock space that are orthogonal to all lower lying excited eigenstates. Finally, Section VII is devoted to summarizing this paper.

II Basic theory

The present formulation can be simply extended to infinite crystal systems, although the summation over $\bm{k}$ is not written explicitly here. The applicability to crystals is obvious in the case of the ground state because there have been so many papers using GW approximation for crystals. How about the case of the excited states? If only some finite number of electrons are excited from the ground state, nothing changes in the infinite system. Therefore, the calculations are only meaningful when the finite density of electrons are excited from the ground state. Such calculations would correspond to strong excitations induced by e.g. laser irradiations.

We start from the definition of the QPs. The initial state can be a ground or any excited state of the $N$-electron neutral system or any excited state of the $M(\neq N)$-electron charged system. We write this state as $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ and its total energy as $E_{\gamma}^{M}$. They satisfy the eigenvalue equation $H\,|\,{\Psi_{\gamma}^{M}}\,\rangle=E_{\gamma}^{M}\,|\,{\Psi_{\gamma}^{M}}\,\rangle$. The electrons are interacting each other by the Coulomb interaction, and therefore cannot be treated independently. However, if we introduce the QP picture, we have a simple way to treat $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ in terms of the QPs. The idea is based on Einstein’s photoelectronic effect, which is called in modern terminology the potoelectron or PE/IPE spectroscopy. After the photoemission process, the system becomes a $(M-1)$-electron excited state $\,|\,{\Psi_{\mu}^{M-1}}\,\rangle$, which satisfies the eigenvalue equation $H\,|\,{\Psi_{\mu}^{M-1}}\,\rangle=E_{\mu}^{M-1}\,|\,{\Psi_{\mu}^{M-1}}\,\rangle$. The corresponding total energy $E_{\mu}^{M-1}$ satisfies

[TABLE]

Therefore, the total energy differences are the measurable quantities. We call these energies the QP energies and define them as

[TABLE]

The QP energy $\varepsilon_{\mu}$ corresponds to the minus of the energy required to remove one electron from the $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ state and add it at at the vacuum level; see Fig. 1. Similarly, the QP energy $\varepsilon_{\mu}$ corresponds to the minus of the energy gain when one electron is added to the $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ state from the vacuum level. Note that the final states includes many ($M\pm 1$)-electron states, which can be created by adding or removing one electron from the initial state. This does not necessarily mean that the final states include all possible ($M\pm 1$)-electron states. For example, multiple scattered eigenstates where many electron configurations are different from the initial eigenstate are excluded or their amplitudes are quite small even if they are included.

The QP wave functions can be defined as an overlap between the initial state $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ with one-electron removed or added at a point $(\bm{r},s)$ ( $\bm{r}$ and s are spacial and spin coordinates of the electron) and the final state $\,|\,{\Psi^{M-1}_{\mu}}\,\rangle$ or $\,|\,{\Psi^{M+1}_{\nu}}\,\rangle$ (therefore, they are sometimes called the overlap amplitudes):

[TABLE]

Here $\psi^{\dagger}_{s}(\bm{r})$ and $\psi_{s}(\bm{r})$ are the (creation and annihilation) field operators. These QP wave functions represent the amplitude of the one-electron difference between the initial and the final states (see Fig. 2). The total electron spin density of the $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ state is given by

[TABLE]

where the completeness condition for the eigenstates $\,|\,{\Psi_{\mu}^{M-1}}\,\rangle$ of the ($M-1$)-electron Hamiltonian

[TABLE]

was used. Similarly, the expectation value of the one-body part of the Hamiltonian

[TABLE]

($v(\bm{r})$ is the local potential such as the Coulomb potential caused by nuclei) are obtained by the QP wave functions as follows:

[TABLE]

In fact, the QP wave functions diagonalize the density matrix

[TABLE]

The spin density $\rho_{s}(\bm{r})$ is given by $\rho_{s}(\bm{r},\bm{r})$, and the expectation values of the kinetic and local potential energies in the $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ state are given by

[TABLE]

In this sense, the QP wave functions $\phi_{\mu}(\bm{r},s)$ have a property similar to the natural spin orbitals introduced by Löwdin more than 60 years ago Lowdin1 ; Lowdin2 in connection with the CI theory Sabo .

Here, it is beneficial to note that the completeness condition

[TABLE]

holds for the QP states. This can be easily seen by sandwiching the anticommutation relation of the field operators

[TABLE]

with $\langle\,{\Psi_{\gamma}^{M}}\,|\,$ and $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ and using the other completeness condition for the eigenstates $\,|\,{\Psi_{\nu}^{M+1}}\,\rangle$ of the ($M+1$)-electron Hamiltonian

[TABLE]

together with (5). On the other hand, the QP wave functions are not orthonormal, i.e., not normalized to unity, and even not linearly independent HedinLundqvist ; Goscinski , and therefore they are different from the normal spin orbitals mentioned above in an exact sense. However, one should notice that the overlap matrix $S$ between the QP wave functions satisfies the idempotency relation $S=S^{2}$, i.e.,

[TABLE]

where the completeness condition (12) was used in the third equality. So, if one diagonalizes the overlap matrix $S$ by a unitary transformation as $\widetilde{S}=U^{\dagger}SU$ according to the canonical orthonormalization procedure introduced by Löwdin Lowdin1 ; Lowdin2 ; Goscinski , it is realized that the resulting diagonal matrix $\widetilde{S}$ satisfies $\widetilde{S}(\widetilde{S}-1)=0$ because the original matrix $S$ satisfies $S(S-1)=0$. Therefore the diagonal element of $\widetilde{S}$, i.e., the eigenvalues of $S$ must be either zero or unity. Thus the restricted set of the transformed functions

[TABLE]

which correspond to the eigenvalue unity only, becomes linearly independent and satisfies the orthonormality condition

[TABLE]

This restricted set of the transformed functions satisfies the idempotency relation $\widetilde{S}=\widetilde{S}^{2}$. Since $\widetilde{S}^{2}$ is given by

[TABLE]

which must be equal to (17), the completeness condition

[TABLE]

is fulfilled also for the transformed functions Goscinski . In practical calculations, however, one usually assumes that the diagonal elements $S_{\lambda\lambda}$ are unity. In such cases, one has to redetermine these norms afterwards, because they are generally not unity at all. In this procedure, the unitary transformed $S^{\prime}$ is diagonal but its diagonal elements are no more zero or unity. It has a form Goscinski

[TABLE]

where the submatrix $\mu$ is a diagonal matrix, whose elements are not necessarily unity. Then, the orthonormal set must be chosen as

[TABLE]

instead of (16). From this equation, the QP wave functions with the correct norm $S_{\lambda\lambda}$ are found to be

[TABLE]

From Eq. (16), they should be equal to

[TABLE]

Comparing these two equations, we have

[TABLE]

Because $U_{\alpha\lambda}$ is unitary,

[TABLE]

is obtained. Therefore, one can determine the norm of each QP wave function $S_{\lambda\lambda}$ using this procedure. Apart from this complexity, the non-orthogonality and linear dependence of the QP wave functions do not bring any essential difficulty in the present formulation.

What we want to emphasize here is that there are vast amount of irrelevant ($M\pm 1$)-electron states, which can be ignored from the summation over all empty states or all occupied states in Eqs. (4)-(9) and (12). For example, one may consider a lot of independent excitations that are not directly associated with the photoemission process under consideration. Even such excited states are the ($M-1$)-electron states, which can be included in the summation over all occupied states. Of course, such states do not contribute to the intermediate states in Eqs. (4)-(9) and (12), because such multiple excitations cannot occur with the single creation or annihilation operator. This is a clear example of the non-orthonormality: The amplitude for such QP states is zero or almost zero. The simplest way is to simply remove such ($M\pm 1$)-electron states from the QP states. (We will use this approximation in our test calculations given in Section V.)

Now, let us derive the equation satisfied by the QP wave functions and the QP energies. Since the Hamiltonian $H$ is given by $H=H^{(1)}+H^{(2)}$ with the one-body part (6) and the two-body part, i.e., the electron-electron Coulomb interaction,

[TABLE]

the commutation relation between $\psi_{s}(\bm{r})$ and $H$ is given by

[TABLE]

Then, by sandwiching this with $\langle\,{\Psi^{M-1}_{\mu}}\,|\,$ and $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ or with $\langle\,{\Psi_{\gamma}^{M}}\,|\,$ and $\,|\,{\Psi^{M+1}_{\nu}}\,\rangle$, and using Eqs. (2) and (3), we derive

[TABLE]

Here, the second term in the right-hand side (r.h.s.) of these equations are rewritten as

[TABLE]

by introducing the energy-dependent self-energy $\Sigma_{s}(\bm{r},\bm{r}^{\prime};\varepsilon_{\lambda})$ with $\varepsilon_{\lambda}=\varepsilon_{\mu}$ or $\varepsilon_{\nu}$ given by Eq. (2). (Here and hereafter, the self-energy includes the Hartree term, which describes the classical Coulomb interaction.) It is clear that both identities (31a) and (31b), which define the self-energy, are consistent with each other, in particular, in the excited-state formulation. In fact, since $\langle\,{\Psi_{\mu}^{M-1}}\,|\,$ and $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ are arbitrary excited eigenstates in Eq. (31a), they can be changed to $\langle\,{\Psi_{\gamma}^{M}}\,|\,$ and $\,|\,{\Psi_{\nu}^{M+1}}\,\rangle$, and then Eq. (31b) is readily obtained. In other words, in order to have the two identities (31a) and (31b) simultaneously, one must recognize that the ordinary formulation for the self-energy for the ground state should hold for arbitrary excited eigenstates also.

From (31a) and (31b), the QP equations

[TABLE]

are derived. Since there is the energy-dependence in the self-energy, we have to solve these equations level to level. Note that the self-energy has a state dependence only through the dependence on a particular QP energy. In the case of the ground state, it has been demonstrated that the energy dependence can be approximately treated by linearization Kresse ; Kuwahara ; Kuwahara2 . However, the energy-dependence is very important in the QP theory. Indeed, multiplying both sides of Eq. (31a) by $\frac{1}{2}\langle\,{\Psi_{\gamma}^{M}}\,|\,\psi^{\dagger}_{s}(\bm{r})\,|\,{\Psi^{M-1}_{\mu}}\,\rangle=\frac{1}{2}\phi^{\ast}_{\mu}(\bm{r},s)$, summing up with respect to $\mu$ and $s$, integrating with respect to $\bm{r}$, and using the completeness condition (5) again, we obtain

[TABLE]

for the expectation value of the electron-electron Coulomb interaction for the $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ state, which is known as the Galitskii and Migdal’s formula Galitskii

[TABLE]

(Notations appearing here will be explained in Section VI.) in the QP representation Sanchez-Friera for the case of the ground state. The crucial points here are that the energy dependence of the self-energy appears in the same way in both Eqs. (32a) and (33) in a general scheme of the non-orthonormal QP wave functions and that Eq. (33) holds not only for the ground state but also for arbitrary excited eigenstates. Equation (33) implies that the information of the self-energy is enough to determine the total energy $E_{\gamma}^{M}$ of the $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ state of the system, because the QP wave functions are uniquely obtained by solving the QP equation (32). In other words, if the self-energy is known, the QP energies and QP wave functions can be calculated by solving the QP equation (32) one by one, and then the electron spin density $\rho_{s}(\bm{r})$, the kinetic energy $\langle\hat{T}\rangle$, the local potential energy $\langle\hat{v}\rangle$, and the interaction energy $E_{\rm int}$ are obtained via Eqs. (4), (7), (8), and (33), respectively. Consequently the total energy $E_{\gamma}^{M}$ is obtained from the knowledge of the self-energy only. How to determine the self-energy will be described in Section IV. In next section (Sec. III), we will show that the equations derived in Section II can be rederived by using the Green’s function.

III Relation to the Green’s Function Formalism

The equations (31) and (32) can be also derived from the Dyson equation for the Green’s function as follows. The Green’s function $G$ on the eigenstate $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ is defined as

[TABLE]

where the convention to write $x=(\bm{r},t)$ is used; $T$ denotes the time-ordered product, and $\psi^{\dagger}_{s}(x)=e^{iHt}\psi^{\dagger}_{s}(\bm{r})e^{-iHt}$ and $\psi_{s}(x)=e^{iHt}\psi_{s}(\bm{r})e^{-iHt}$ are, respectively, the creation and annihilation operators in the Heisenberg representation. By means of the completeness conditions (5) and (14),

[TABLE]

is readily derived. Here we put

[TABLE]

and $\phi_{\lambda}(\bm{r},s)$ and $\varepsilon_{\lambda}$ ($\lambda=\mu$ or $\nu$) are the QP wave functions and the QP energies defined in Eqs. (2) and (3).

The Heisenberg equation of motion for the annihilation operator is given by

[TABLE]

Operating both sides $\psi_{s}^{\dagger}(x^{\prime})$ from the right and the time-ordered operator $T$ from the left, and sandwiching them with $\langle\,{\Psi_{\gamma}^{M}}\,|\,$ and $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$, we have

[TABLE]

where $x^{\prime\prime+}$ with $t^{\prime\prime}=t$ means $(\bm{r}^{\prime\prime},t^{+})$ with $t^{+}=t+0^{+}$. The self-energy is introduced as

[TABLE]

and then the Dyson equation

[TABLE]

holds for any eigenstate. Equation (41) is rewritten in the Fourier space as

[TABLE]

The Fourier transform of the Green’s function (36) with (37) is given by

[TABLE]

where $\delta_{\lambda}=0^{+}$ for occupied states and $\delta_{\lambda}=-0^{+}$ for empty states, i.e., $\delta_{\mu}=0^{+}$ and $\delta_{\nu}=-0^{+}$. Inserting this into Eq. (42) and approaching $\omega$ to one of the $\varepsilon_{\lambda}$’s ($\lambda$ can be either $\nu$ (emp) or $\mu$ (occ)), we obtain the QP equation

[TABLE]

which is equivalent to Eqs. (32a) and (32b). If the eigenvalues are $m$-fold degenerate, one can choose the corresponding $m$ QP wave functions to be orthogonal to each other in this degenerate subspace. Then, substituting Eq. (43) with $\omega=\varepsilon_{\lambda}$ into Eq. (42), multiplying the resulting equation by the complex conjugate of one of the $m$ QP wave functions $\phi^{\ast}_{\lambda}(\bm{r}^{\prime},s)$, and integrating with respect $\bm{r}^{\prime}$, we derive Eq. (44).

With the aid of Eqs. (36) and (37), the left-hand side (l.h.s.) of Eq. (40) can be written as

[TABLE]

whereas the r.h.s. of Eq. (40) can be written as

[TABLE]

Comparing these two equations, we find

[TABLE]

for $t>t^{\prime}$. Here we put

[TABLE]

for $t>t^{\prime}$. Note that $f_{\mu}(x,s;t^{\prime})$ and $f_{\nu}(x,s;t^{\prime})$ in Eq. (47) are implicitly defined by Eq. (48). Since the r.h.s. of Eqs. (47a) and (47b) does not depend on $t^{\prime}$, the l.h.s. should not depend on $t^{\prime}$ also. Therefore, we can put $t^{\prime}=\infty$ in Eqs. (47a) and (48a), and $t^{\prime}=-\infty$ in Eqs. (47b) and (48b). Since $f_{\mu}(x,s;\infty)=f_{\nu}(x,s;-\infty)=0$, we finally obtain

[TABLE]

These equations are definitely identical to Eqs. (31a) and (31b).

IV Perturbation Series for the Self-Energy

Next, we derive a concrete method to calculate the self-energy $\Sigma_{s}$. Because here we discuss the energies at an instant time only, even when the system is slowly varying in time, we simply drop all the time dependence in the following expressions.

First, by the eqautions

[TABLE]

and $a_{\nu}\,|\,{\Psi_{\gamma}^{M}}\,\rangle=a^{\dagger}_{\mu}\,|\,{\Psi_{\gamma}^{M}}\,\rangle=0$ for an empty QP state $\nu$ and an occupied QP state $\mu$, we introduce the annihilation and creation operators, $a_{\lambda}$ and $a^{\dagger}_{\lambda}$, as

[TABLE]

where $\lambda$, $\lambda^{\prime}$ includes both $\mu$ and $\nu$. Then, the anticommutation relations of the field operators are written as

[TABLE]

From these equations, we can assume that the QP annihilation and creation operators satisfy the anticommutation relations

[TABLE]

Note that these relations hold even when the QP states $\lambda,\lambda^{\prime}$ are not orthogonal to each other, because Eq. (52) with Eq. (53) satisfies the completeness condition (12) only.

Second, we introduce a hypothetical non-interacting system, where the QPs are not interacting each other, although in real systems the QPs generally interact with each other at least at short distances. This hypothetical system is described by the Hamiltonian

[TABLE]

where $H^{(2)}_{\rm self}$ is given by (28) whose operation is restricted to construct the self-energy $\Sigma_{s}$ only and ineffective for any other interaction. Its eigenvalue equation is given by

[TABLE]

where $M+n$ indicates the total number of electrons in this system. When $n=0$ and $\alpha=\gamma$, we assume that $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$ has the same energy eigenvalue ${E^{\prime}}^{M}_{\gamma}=E_{\gamma}^{M}$ as the true system. All the other states of the hypothetical system are constructed by operating several $a_{\lambda}$’s and $a^{\dagger}_{\lambda}$’s on $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$. We assume $H^{\prime}$ as the unperturbed Hamiltonian and

[TABLE]

as a perturbation, and adopt Brillouin–Wigner’s perturbation theory, which is well known to be applicable not only to the ground state but also to any excited eigenstate of the full Hamiltonian. Multiplying both sides of

[TABLE]

by the projection operator

[TABLE]

and noting that $H^{\prime}$ commutes with $P$, we have

[TABLE]

Equation (59) can be rewritten as

[TABLE]

Then, by iteration, we have

[TABLE]

So, the expectation value of $A(\bm{r},\bm{r}^{\prime})=\psi_{s}^{\dagger}(\bm{r})\psi_{s^{\prime}}^{\dagger}(\bm{r}^{\prime})\psi_{s^{\prime}}(\bm{r}^{\prime})\psi_{s}(\bm{r})$ can be evaluated as

[TABLE]

where we replaced $H^{\prime\prime}$ with $H^{(2)}$ without restriction because it is obvious that these interactions are not related to the construction of the self-energy but simply related to the construction of the expectation value of $A(\bm{r},\bm{r}^{\prime})$. Then, multiplying the both sides of Eq. (62) by $\frac{1}{2}h_{ss^{\prime}}^{(2)}(\bm{r},\bm{r}^{\prime})$, summing over $s,s^{\prime}$, and integrating with respect to $\bm{r}$ and $\bm{r}^{\prime}$, we have

[TABLE]

In the $k$-th order term in this perturbation series, there are $4(k+1)$ field operators sandwiched by $\langle\,{\Phi_{\gamma}^{M}}\,|\,$ and $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$, which can be expanded as (51a) and (51b). For example, the first term in the r.h.s. of Eq. (63) is expressed as

[TABLE]

In this matrix element, the QP states $\mu_{4}$ and $\mu_{3}$ are removed from and $\mu_{2}$ and $\mu_{1}$ are added to the $M$-particle state $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$, and after these operations, the $M$-particle state must be restored to the original state $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$; otherwise this matrix element becomes zero. That is, either ($\mu_{1}=\mu_{4}$ and $\mu_{2}=\mu_{3}$) or ($\mu_{1}=\mu_{3}$ and $\mu_{2}=\mu_{4}$) should be satisfied. Of course, this first order approximation rigorously corresponds to the Hartree–Fock energy in the extended HFA Morokuma . In fact, Eq. (64) becomes

[TABLE]

Its analogy holds for all orders. That is, all the intermediate QP states $\lambda_{1},\lambda_{2},\lambda_{3},\lambda_{4},\cdots$ must be pairwise ($\lambda_{i}$ can be either occupied states $\mu_{i}$ or empty states $\nu_{i}$); each pair consists of one annihilation operator $a_{\lambda_{i}}$ and one creation operator $a^{\dagger}_{\lambda_{i}}$ of the same QP state $\lambda_{i}$. All possible combinations of the product of the pairs have to be summed up to evaluate the total contribution at this order. Because of the existence of the projection operator $P$ in (63), the intermediate states cannot be the initial state $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$. Moreover, the energy denominator in (63) is just a simple addition or subtraction of the QP energies $\varepsilon_{\lambda}$’s. If we draw a pair of the QP state by a thick sold line and $H^{(2)}$ by a dotted line, (63) is expressed by the diagrams shown in Fig. 3. This expression is exactly the same as the interaction energy expressed by the skeleton diagrams, i.e., the diagrams expressed by the clothed Green’s function (a thick solid line in Fig. 3), in the usual MBPT Feynman ; Luttinger ; Pines ; Nozieres . They also correspond to the same order contribution in the Møller–Plesset theory if the initial excited-state configuration and the full self-consistency are imposed in the latter. Thus we derived exactly the same formula only by using the true eigenstates of the full Hamiltonian and their derivatives such as the QP energies and QP wave functions. An important point here is that we have never assumed the orthogonality of the QP wave functions in these equations.

Third, in order to obtain the expression for the self-energy $\Sigma_{s}$, we compare the expression (63) with the Galitskii and Migdal’s formula of the interaction energy (33). Looking at Eq. (33), we find that if we drop $\phi_{\mu}(\bm{r},s)\phi^{\ast}_{\mu}(\bm{r}^{\prime},s)$ from (63), i.e., a thick solid line from Fig. 3, we should have the self-energy $\Sigma_{s}(\bm{r},\bm{r}^{\prime};\varepsilon_{\mu})$ equivalent to Fig. 4. In Fig. 4, we drew two external thin solid lines to indicate the dependence of the self-energy $\Sigma_{s}$ on the two positions $(\bm{r},s)$ and $(\bm{r}^{\prime},s)$. This is again exactly the same skeleton diagrams for the self-energy in the MBPT Pines ; Nozieres . When we remove a pair of $\phi_{\mu}(\bm{r},s)\phi^{*}_{\mu}(\bm{r}^{\prime},s)$ from the product of $H^{(2)}$’s in the expression of $E_{\rm int}$ (63), we also have to remove the summation with respect to this special QP state, $\mu$. Therefore, the corresponding QP energy $\varepsilon_{\mu}$ appearing in some parts of the energy denominators $(E^{M}_{\gamma}-H^{\prime})^{-1}$ remains unsummed with the special index $\mu$ after taking the summations with respect to all the intermediate QP states, resulting in the $\varepsilon_{\mu}$ dependence in the final expression. Consequently, the resulting self-energy must have the explicit $\varepsilon_{\mu}$ dependence, where $\mu$ corresponds to the QP state of the external thin line in Fig. 4.

Last, we make a brief comment on the relationship between the present formulation using the resolvent operator $(E^{M}_{\gamma}-H^{\prime})^{-1}$ and the formulation using the clothed Green’s function (43). To see this relationship, we insert the completeness condition

[TABLE]

in between each two $H^{(2)}$’s in Eq. (63). Then $H^{(2)}$ is replaced by the sum over all possible intermediate QP states $\lambda_{1},\lambda_{2},\lambda_{3},\lambda_{4}$ of the product of four QP wave functions and the matrix element $\langle\,{\Phi_{\lambda}^{M}}\,|\,\cdots\,|\,{\Phi_{\lambda}^{M}}\,\rangle$ replaced with $\langle\,{\Phi_{\alpha}^{M}}\,|\,\cdots\,|\,{\Phi_{\beta}^{M}}\,\rangle$ in a form (64). (Note that the number of electrons remains $M$ in the intermediates states $\alpha$ and $\beta$, because $H^{(2)}$ has two annihilation and two creation operators.) Then, all the QP states $\lambda_{i}$ are made pairwise, each of which has a form $\phi_{\lambda}(\bm{r},s)\phi^{*}_{\lambda}(\bm{r}^{\prime},s)$ identical to the numerator of the Green’s function (43). Simultaneously, each resolvent operator is replaced with its expectation value $(E^{M}_{\gamma}-E^{\prime M}_{\alpha})^{-1}$ sandwiched by an intermediate state $\,|\,{\Phi^{M}_{\alpha}}\,\rangle$, which corresponds to the denominator of the result of the $\omega$-integration of the even-number product of the Green’s functions. As an example, let us consider the third term of Fig. 3. It consists of two bubble diagrams, each of which corresponds to the polarization function $P(\bm{r}_{i},\bm{r}^{\prime}_{i};\omega_{i})$ ($i=1,2$) in the random-phase approximation (RPA), connected by two dotted lines, $h_{s_{1}s_{2}}^{(2)}(\bm{r}_{1},\bm{r}_{2})$ and $h_{s_{1}s_{2}}^{(2)}(\bm{r}^{\prime}_{1},\bm{r}^{\prime}_{2})$. The polarization function is evaluated as

[TABLE]

as usual, where $f(\varepsilon)$ is the Fermi distribution function. Multiplying $P(\bm{r}_{1},\bm{r}^{\prime}_{1};\omega)$ by $P(\bm{r}^{\prime}_{2},\bm{r}_{2};\omega)=P(\bm{r}_{2},\bm{r}^{\prime}_{2};-\omega)$ and $i/2\pi$, and integrating it with respect to $\omega$ from $-\infty$ to $\infty$, we obtain a possible term,

[TABLE]

whose denominator is identical to the expectation values of the resolvent operator $(E^{M}_{\gamma}-H^{\prime})^{-1}$ in the second term of Eq. (63). Final result can be obtained by operating

[TABLE]

The eight product of the QP wave functions comes from the expansions (51a) and (51b). All these operations correspond to taking the trace, $\mathop{\rm Tr}\nolimits[Ph^{(2)}Ph^{(2)}]$, in the Green’s function method Kuwahara .

V Test Calculations

In order to verify the applicability of the QP theory for some initial excited eigenstates, let us first consider a hydrogen atom with $N=1$ and for the case $M=0$ as the simplest case. Its ($M=0$)-electron excited state is none other than a proton. Then, the lowest-energy (empty) QP state having $M+1=1$ electron is the hydrogen 1s orbital. The self-consistent calculation surely gives this result exactly, because there is no occupied electron. Note that there is no Hartree term, no exchange term, and no higher-order terms. This is so even in the LDA or the generalized gradient approximation (GGA) of DFT. The same thing holds not only for the 1s state but also for any excited eigenstates of the hydrogen atom. Similar situation occurs also for an alkali metal atom like lithium, sodium, and potassium, if we can ignore the core contributions. In this sense, the present excited-state formalism has a chance to serve an extremely good result compared to the ordinary ground-state formalism in some special cases.

Next we consider a general isolated atom A or molecule M. It is obvious that the ionization potential (IP) of A or M is identical to the electron affinity (EA) of its cation A+ or M+, provided that the atomic geometries of M and M+ are the same (This is not a bad approximation for typical molecules). Similarly, there are several other identities between the electron removal and absorption energies of cationic, anionic and photoabsorbed states of A or M. These identities allow us to check the validity of the present theory. Here, we perform the perturbation expansion starting from the cationic, anionic, or photoabsorbed state as well as the neutral ground state of A or M to calculate the EQP energies, which should coincide with each other.

We use the standard one-shot $GW$ approximation ($G_{0}W_{0}$) by beginning with the LDA HybertsenLouie for both the ground-state and excited-state configurations of isolated lithium, aluminum, and beryllium atoms and nitrogen, oxygen and lithium (diatomic) molecules with the bond lengths fixed at 1.106 $\mathrm{\mathring{A}}$, 1.216 $\mathrm{\mathring{A}}$, and 2.723 $\mathrm{\mathring{A}}$, respectively, for N2, O2, and Li2 (calculated with B3LYP/6-31G* NIST ). We use the all-electron mixd basis, TOMBO code TOMBO , in which both plane waves (PWs) and atomic orbitals (AOs) are used as basis functions. The face-centered cubic (fcc) unit cell with an edge length of 18 $\mathrm{\mathring{A}}$ is used for Li, Li2, and Al, and that of 12 $\mathrm{\mathring{A}}$ is used for Be, N2, and O2. Cutoff energies for PWs, the Fock exchange and the correlation part of the self-energy (as well as the polarization part) are set as 11.05 Ry, 44.22 Ry, 11.05 Ry for Li, Li2, and Al, 19.65 Ry, 78.61 Ry, 19.65 Ry for Be and N2, and 27.71 Ry, 122.83 Ry, 30.71 Ry for O2, respectively. In the summation in the correlation part, 400 levels are included for Li, Al, and Be; 3000, 6000, and 7000 levels are included, respectively, for Li2, N2, and O2. We use the generalized plasmon pole model HybertsenLouie to avoid the $\omega$ integration. The energy dependence of the self-energy is dealt by introducing the $Z$ factor (linearization around the Kohn-Sham eigenvalues) as usual.

The computational time scales as $O(N_{\rm emp}N_{\rm occ}N_{G}^{2})$ where $N_{\rm emp}$, $N_{\rm occ}$, and $N_{G}$ are the numbers of empty states, occupied states, and $G$ vectors needed for the non-local Fourier components of the polarization function and the correlation part of the self-energy.

There is no essential difference between the QP theories using the neutral ground state and the other (neutral or charged) excited states as the initial state, and the GW calculations for the latter excited-state case can be simply performed by just changing the order of levels after the diagonalization. However, in order to compare the results, we will call the former neutral ground-state case the ordinary QP (OQP) theory and the latter excited-state case the extended QP (EQP) theory.

The calculated results of the QP energies of several QP states are listed in Table 1 for Li, Al, and Be, and in Table 2 for N2, O2, and Li2, together with the reference values. EC in Table 1 means the value calculated by the EOM coupled cluster singles and doubles (EOM-CCSD) with aug-cc-pV5Z basis set using Gaussian 09 Gaussian . References Ref5 ; Ref6 ; Ref7 ; Ref8 ; Ref9 ; Ref11 ; Molecule1 ; Molecule2 ; Molecule3 ; Molecule4 are all experimental data. We find good agreements between the EQP energies of the forward ($\rightarrow$) electron absorption process and the backward ($\leftarrow$) electron release process. Maximum deviation from the reference value is 1.0 eV in the case Be+ $\leftarrow$ Be∗ for atoms, and 0.6 eV in the case Li${}_{2}^{+}$ $\leftarrow$ Li${}_{2}^{*}$ for molecules, which are smaller than the deviations in the OQP theory (values inside parentheses in Tables 1 and 2); 3.2 eV in the case Al${}^{+}*$ $\leftarrow$ Al. The average deviations are 0.4 eV and 0.3 eV, respectively, for atoms and molecules; they are the same orders for both the OQP and EQP theories (values outside parentheses in Tables 1 and 2). This means that there is no specific deviation in the EQP theory, demonstrating the validity of the present theory. To obtain better agreement, it would be necessary to use the self-consistent GW($\Gamma$) methods Kuwahara ; Kuwahara2 . The resulting QP energy spectra of Li${}^{+}*$ calculated by $G_{0}W_{0}$ are shown in Fig. 5. It is seen from this figure that the agreement between the theory and experiment is very good already in the present approximation.

VI Several Important Remarks

Before ending this paper, we will give several important remarks.

(1) First of all, it is easy to see that the virial theorem holds for any excited eigenstate of the $M$-electron Hamiltonian

[TABLE]

where $\{\bm{r}_{i}\}$ and \{\mbox{\boldmathR}_{k}\} denote the positions of electrons and nuclei. In fact, since the $M$-electron wave function $\Psi_{\gamma}^{M}\{\bm{r}_{i}\}$ of this excited eigenstate $\gamma$ is normalized as $\int\Psi^{M\ast}_{\gamma}\{\alpha\bm{r}_{i}\}\Psi_{\gamma}^{M}\{\alpha\bm{r}_{i}\}\prod_{i}^{M}d(\alpha\bm{r}_{i})=1$ when all lengths are scaled by $\alpha$, the total energy of this excited eigenstate at the equilibrium atomic positions \{\mbox{\boldmathR}_{k}\},

[TABLE]

has the $\alpha$-dependence only in the form of $\langle\hat{T}\rangle/\alpha^{2}+\langle\hat{V}_{\rm total}\rangle/\alpha$, where H^{M}\{\bm{r}_{i},\mbox{\boldmathR}_{k}\} is given by (71), $\langle\hat{T}\rangle$ and $\langle\hat{V}_{\rm total}\rangle$ denote the kinetic and potential contributions to $E^{M}_{\gamma}$. If the atomic configuration is fully relaxed at this excited state, the total energy must be stationary with respect to the change in $\alpha$. Therefore, differentiating $E^{M}_{\gamma}$ with respect to $\alpha$ and setting $\alpha$ at unity, the virial theorem $2\langle\hat{T}\rangle+\langle\hat{V}_{\rm total}\rangle=0$ is readily derived.

(2) Second, the adiabatic switch on process in the MBPT can be extended to excited eigenstates as is written in Section 5-4 in Ref. Nozieres . In fact, Wick’s theorem holds also for an excited eigenstate if we transform creation and annihilation operators from the coordinate representation to the one-particle state representation of the corresponding noninteracting ($H^{(2)}=0$) system, exchange the role of creation and annihilation operators for occupied and empty states of the excited eigenstate, and introduce the normal product as usual. When the interaction $H^{(2)}$ is gradually switched on, the QP levels may apparently change the order. However, this does not matter at all, because, in the case of an excited state, it is not necessary to require that the noninteracting ground state moves onto the interacting ground state without level crossing. Even if the order of occupied and empty levels changes during the switch on process, we can keep track of each QP state; see Fig. 6. In the noninteracting system, the QP states refer to the one-particle states. There is a one-to-one correspondence between the one-particle states of the noninteracting system and the QP states of the interacting systems, although the orders of these states may be different. There is no constraint on the orders of the one-particle states and the QP states, and they can be different to each other. Therefore, we should not use the standard words “the adiabatic switch on process” for this process because the word “adiabatic” involves the meaning of no level crossing. Instead, it would be better to use the words “the gradual switch on process” to avoid confusion. The gradual switch on of the interaction from time $t=-\infty$ to $t=0$ transforms the noninteracting excited state (at $t=-\infty$) into the interacting excited state (at $t=0$).

(3) Third, according to the Klein’s derivation (see, for example, Section 5-6 in Ref. Nozieres ), the total energy can be alternatively written as

[TABLE]

where $G^{(0)}_{s}$ is the Green’s function of the non-interacting system, i.e., the system described by the one-body part of the Hamiltonian $H^{(1)}$ only, $\mathop{\rm Tr}\nolimits=i\sum_{s}\int d\bm{r}\int_{-\infty}^{\infty}e^{i\eta\omega}d\omega/2\pi$ with $\eta=0^{+}$, and $\Phi[G]$ is the the Luttinger–Ward functional LW ; Nozieres ; Kuwahara

[TABLE]

Here $\Sigma^{(k)}_{s}(\bm{r},\bm{r}^{\prime};\omega)$ represents the $k$-th order skeleton diagram of the self-energy. In Eq. (73), $\mathop{\rm Tr}\nolimits[G_{s}/G^{(0)}_{s}-1]$ is no other than $\mathop{\rm Tr}\nolimits[\Sigma_{s}G_{s}]=-2E_{\rm int}$ (see Galitsukii–Migdal’s formula (34) and the Dyson equation (42) with $G^{(0)\,-1}_{s}=\omega-h^{(1)}_{s}$), and $-\mathop{\rm Tr}\nolimits\ln G_{s}/G^{(0)}_{s}$ can be evaluated as $\sum_{\mu}^{\rm occ}\varepsilon_{\mu}-E^{M}_{\gamma}[G^{(0)}]$ and $E^{M}_{\gamma}[G^{(0)}]=\sum_{\mu}^{\rm occ}\varepsilon^{(0)}_{\mu}$ with the eigenvalue of the non-interacting system $\varepsilon^{(0)}_{\mu}$. We have also the variational principles $\delta\Phi[G]/\delta G_{s}=-\Sigma_{s}$ and $\delta E^{M}_{\gamma}[G]/\delta G_{s}=0$.LW ; Nozieres

(4) The Ward–Takahashi identity Ward ; Takahashi ; WT ,

[TABLE]

which is equivalent to the continuity equation for the electron density (gauge invariance), holds also for the excited eigenstate. Here $\Gamma^{(\alpha)}_{s}(x,y,z)$ are the scalar ($\alpha=0$) and vector ($\alpha=1,2,3$) vertex functions defined as

[TABLE]

where the four-current density $j_{\alpha}(x)$ is defined as

[TABLE]

with the vector potential \mbox{\boldmathA}(x).

(5) For the discussion of transport properties, it is important to guarantee the macroscopic conservation laws of the total energy, the total number of electrons, the total momentum, and so on. These macroscopic conservation laws have been rigorously discussed by Baym and Kadanoff BaymKadanoff ; Baym for the ground state. Their argument may be generalized to the case of an arbitrary excited eigenstate, and the macroscopic conservation laws hold for arbitrary excited states, when the necessary symmetries suggested by them are satisfied for the excited-state Green’s functions. In general, however, we expect that the total number of electrons is conserved in our excited-state theory because the local charge conservation is guaranteed by the Ward–Takahashi identity (75). Moreover, we expect that the total energy and the total momentum of the electronically excited many-atom systems should be conserved, if the system is isolated. This is simply due to the time-reversal symmetry and the translational symmetry in space of the total Hamiltonian of the many-atom system.

(6) However, the fluctuation-discipation theorem is not guaranteed for arbitrary excited states because it is based on the canonical ensemble, i.e. the thermal equilibrium, at a given temperature. Therefore, the linear response theory is not developed here and left for the future study.

(7) Last, we discuss how the DFT is modified for the excited state. Similar to the ordinary DFT, we introduce a universal function $F$ defined as

[TABLE]

where $H_{\rm eg}$ is the Hamiltonian of the interacting electron gas. By definition, $F$ is a unique functional of the eigenstate $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$ of the full Hamiltonian $H$. This statement can be also straightforwardly deduced from the Luttinger–Ward theory, because $F$ is expressed as a unique functional of the Green’s function, which is by definition (see Eq. (35)) a unique functional of the eigenstate $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$. (Note that the rest contribution to the total energy is just the external potential $v(\bm{r})$ coupled to the electron density $\rho(\bm{r})=\sum_{s}\rho_{s}(\bm{r})$.) Below we show that this eigenstate $\,|\,{\Phi_{\gamma}^{M}}\,\rangle$ and therefore $F$ are unique functionals of the electron density $\rho(\bm{r})$ and the occupation numbers of all KS levels. For this purpose, we use the reductio ad absurdum along the same line as the Hohenberg–Kohn’s theory HohenbergKohn . Let us assume that two different local potentials $v(\bm{r})$ and $v^{\prime}(\bm{r})$ give the same electron density $\rho(\bm{r})$, and show that this assumption leads to a contradiction. The Hamiltonians with the local potentials $v(\bm{r})$ and $v^{\prime}(\bm{r})$ are written as $H$ and $H^{\prime}$, respectively. Their eigenvalues and eigenstates are written as $E^{M}_{\gamma}$, ${E^{\prime}}^{M}_{\gamma}$, $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$, and $\,|\,{{\Psi^{\prime}}_{\gamma}^{M}}\,\rangle$, which satisfy the eignvalue equations $H\,|\,{\Psi_{\gamma}^{M}}\,\rangle=E^{M}_{\gamma}\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ and $H^{\prime}\,|\,{{\Psi^{\prime}}_{\gamma}^{M}}\,\rangle={E^{\prime}}_{\gamma}^{M}\,|\,{{\Psi^{\prime}}_{\gamma}^{M}}\,\rangle$. The excited eigenstate is orthogonal to all lower excited eigenstates as well as the ground state. Therefore, we can still use the variational principle that $E^{M}_{\gamma}$ (or ${E^{\prime}}_{\gamma}^{M}$) become minimum when $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ (or $\,|\,{{\Psi^{\prime}}_{\gamma}^{M}}\,\rangle$) is the true excited eigenstate inside the Fock space orthogonal to all lower excited eigenstates and the ground state. Then, inside this Fock space, we expect Ohno

[TABLE]

Adding these two inequalities, we obtain

[TABLE]

Thus we find that the initial assumption was wrong and that $v(\bm{r})$ is a unique functional of $\rho(\bm{r})$ when the corresponding state is constructed as the $\gamma$-th excited eignstate. Since $v(\bm{r})$ determines the form of the Hamiltonian $H$ uniquely, its eigenstate $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ itself is a unique functional of the density $\rho(\bm{r})$ and the information that this is the excited level $\gamma$. Since the latter information is equivalent to the knowledge of the occupation numbers of all KS states, we thus conclude that $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ and therefore $F$ are unique functionals of the electron density $\rho(\bm{r})$ and the occupation numbers of all KS states only. The variational principles still holds within the Fock space orthogonal to all lower excited states and the ground state, which can be controlled by the choice of the occupation numbers of KS states. Therefore only additional requirement in the extension of the ordinary DFT to the electronic excited state is that we have to impose the constraint on the occupation numbers of all KS states. In other words, the ordinary DFT still holds for the excited state if we additionally give the information on the occupation numbers of all KS states. This part only influences the exchange-correlation potential in the DFT. We expect that there is a close similarity between the DFT for excited states and the QP theory for initial excited eigenstates. This argument encourages a direct use of the DFT to an electronic excited state such as in MauriCar , if a reliable exchange-correlation functional is used.

VII Summary

In this paper, we have shown that the quasiparticle (QP) picture, which correspond to the photoemission or inverse potoemission (PE/IPE) spectroscopy, holds exactly for arbitrary excited eigenstate $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ of the system. The QP energy spectra correspond to the PE/IPE energy spectra, and the QP wave functions yields the electron density and the expectation values of the kinetic and local potential energies. They are the solution of the QP equations coupled with the equation for the energy-dependent self-energy. The interaction energy $E_{\rm int}$ of $\,|\,{\Psi_{\gamma}^{M}}\,\rangle$ eigenstate is expressed by the energy-dependent self-energy and the QP wave functions in Galitski–Migdal’s form. All the formulation is given either with or without introducing the Green’s function. We have tried to drive the closed set of equations as simple as possible to be readable for general readers who are not familiar to the field theoretical MBPT. In this derivation, we showed that the energy-dependence of the self-energy is not an essential difficulty. Some simple calculations of isolated atoms and molecules using the one-shot $GW$ approximation demonstrate the validity of the present theory. We also showed that the virial theorm, the Ward–Takahashi identity, and the conservation laws hold and the Luttinger–Ward functional exists for arbitrary excited eigenstates. Last, we gave a comment on the extension of density functional theory (DFT) to arbitrary excited eigenstate. We hope that the present theory makes a breakthrough toward the first-principles calculations of arbitrary electronic excited states of materials.

Acknowledgements.

We thank Steven G. Louie for helpful discussions. This work has been supported by the grant-in-aid for Scientific Research on Innovative Areas (Grant No. 25104713) from MEXT and the grant-in-aid for Scientific Research B (Grant No. 25289218) from JSPS, and also by the high performance computer infrastructure (HPCI) (Project IDs. hp150231, hp150259, hp160234) and the post ”K” projects (priority issue 7 and germinating issue 1) both promoted by MEXT.