13.3 Generalized Kohn Sham scheme for RIRPA

The RIRPA energy, discussed in Sec. 13.1, is non-self-consistently computed using canononical KS orbitals obtained from a prior KS-DFT based calculation. If the KS-DFT densities are poor, the errors can carry over to post-KS RIRPA energetics resulting in poor interaction energies in some cases. Such density-driven errors can be alleviated by an orbital self-consistent approach called the generalized Kohn–Sham semi-canonical projected RPA (GKS-spRPA). In the GKS-spRPA method, the spRPA energy functional \[\begin{equation} E^{\text{spRPA}}[\mathbf D, \tilde{\mathbf H}_0^{\text{KS}}[\mathbf D]] = E^{\text{HF}}[\mathbf D] + E^{\text{C spRPA}}[\mathbf D, \tilde{\mathbf H}_0^{ \text{KS}}[\mathbf D]] \end{equation}\](13.17) is minimized with respect to the non-interacting GKS density matrix \[\begin{equation} \mathbf D = \sum_{\lambda}\mathbf{P}_{\lambda}\mathbf{n}_{\lambda\lambda'} \mathbf{P}_{\lambda}\text{.} \end{equation}\](13.18) \(\mathbf{D}\) is constrained to be normalized to \(N\) electrons and have eigenvalues between 0 and 1. \(\mathbf P_\lambda\) denotes orthogonal projectors belonging to blocks of KS orbitals with degenerate occupation numbers, and \(\mathbf n_{\lambda\lambda'}=\mathbf n_{\lambda} \delta_{\lambda\lambda'}\) is diagonal, with \(\mathbf n_\lambda\) denoting occupation number matrices. For integer KS occupations \(\mathbf n_\lambda\) has eigenvalues \(n_\lambda = 1,0\). The semicanonical projected (sp) KS Hamiltonian \[\begin{equation} \tilde{\mathbf{H}}_{0}^{\text{KS}} = \sum_{\lambda}\mathbf{P}_\lambda \mathbf{H}_{0}^{\text{KS}}\mathbf{P}_\lambda\text{,} \end{equation}\](13.19) contains only the diagonal (\(\lambda=\lambda'\)) blocks of the KS Hamiltonian \[\begin{equation} \mathrm H_{0,ij}^{\text{KS}} = \mathrm h_{ij} + \sum_{pq}\mathrm V_{ipjq}\mathrm D_{pq} + \mathrm V^{\text{XC}}_{ij} [\mathrm D]\;\text{.} \end{equation}\](13.20) \(\mathbf h\) is the one-electron Hamiltonian, the second term denotes the Hartree or Coulomb potential and \(\mathbf V^{\text{XC}}\) is the exchange-correlation potential. The eigenvalues of sp KS Hamiltonian are denoted by \(\tilde \varepsilon^{KS}_i\). \(\mathbf V\) is the matrix of two-electron integrals \[\begin{equation} \mathrm V_{pqrs} = \int \int d^3r_1 d^3r_2 \frac{\phi_p^{*}(\mathbf r_1)\phi_q^{*}(\mathbf r_2)\phi_r(\mathbf r_1)\phi_s (\mathbf r_2)}{\vert \mathbf r_1 - \mathbf r_2 \vert}\;\text{.} \end{equation}\](13.21) The subscripts \(i,j,..\) denote orbital indices. The spRPA energy functional, eq. 13.17, thus generalizes the post-KS RPA energy functional for arbitrary density matrices. A GKS energy minimization in the space of density matrices leads to a set of one-particle equations, \[\begin{equation} {} \mathbf H^{\text{spRPA}} [\mathbf D] \phi_p = \varepsilon_p^{ \text{GKS-spRPA}}\phi_p \;\text{,} \end{equation}\](13.22) which are solved self-consistently. The one-particle spRPA Hamiltonian is the functional derivative of the RPA energy functional, and contains contributions from the HF potential and RPA correlation potential, \[\begin{equation} \mathbf H^{\text{spRPA}} [\mathbf D] = \frac{\delta E^{\text{spRPA}}[\mathbf D]}{\delta \mathbf D} = \mathbf H^{\text{HF}}[\mathbf D] + \mathbf V^{\text{C\;spRPA}}[\mathbf D] \;\text{.} \end{equation}\](13.23) The self-consistent procedure used in eq. 13.22 is similar to that used in Hartree-Fock and KS-DFT. At the stationary point of the minimization procedure, the GKS-spRPA procedure yields a total many-body energy and one-particle orbital energies, \(\varepsilon_p^{ \text{GKS-spRPA}}\). The latter are related to approximate ionization potentials and electron affinities. \(\varepsilon_p^{ \text{GKS-spRPA}}\) provide consistently accurate estimates of IPs and EAs due to the inclusion of static Hartree-exchange effects, via \(\mathbf H^{\text{HF}}\) contribution, and orbital-correlation and orbital-relaxation effects, via \(\mathbf V^{\text{C\;spRPA}}[\mathbf D]\) contribution. \(\varepsilon_p^{ \text{GKS-spRPA}}\) provide accurate estimates of valence ionization potentials for neutral and anionic systems [306]. The balanced inclusion of orbital-correlation and -relaxation effects leads to accurate estimates of absolute and relative core-ionization energies [313].

The computation of the complete one-particle spectrum currently has a scaling of \(\mathcal O(N^6)\). To reduce the computational costs, we carry out the minimization in two steps: (i) We approximate the occupied-occupied and virtual-virtual blocks of \(\mathbf H^{\text{spRPA}} \approx \mathbf H^{HF}\) for the self-consistent procedure and obtain the converged total energy, and then (ii) evaluate the one-particle eigenvalues at the final converged stationary point only. This two step procedure effects the rate of energy-convergence only and not the final total energy while maintaining a computational scaling of \(\log(N)\mathcal O(N^4)\) for energy evaluation in each iteration of step (i). Thus the cost of GKS-spRPA total energy is number of iterations times the cost for the computation of right-hand side vector, which has a scaling \(\log(N)\mathcal O(N^4)\). (Note: The computation of the one-particle eigenspectrum, i.e. the step (ii), will be available in a future release.)

We note that the GKS-spRPA energy functional has a parametric dependence on the KS potential. However, the variational GKS minimization of the spRPA energy functional reduces the dependency; the quality of results, for energy differences and IPs/EAs, is similar for different KS potentials. We also note that GKS energy minimization of small-gap systems may have poor or no convergence. For such cases, the convergence can be improved by a (small) artificial-shift of the KS orbital energy differences.

Prerequisites

GKS-spRPA total energy calculation requires

• a converged KS-DFT calculation which provides the starting guess MOs for GKS-spRPA. scheme.

• relevant rirpa-options:

  • npoints \(\langle \text{integer} \rangle\) - Number of frequency integration points (default is 60).

  • iter \(\langle \text{integer} \rangle\) - Turns on GKS-spRPA self-consistent iterations. \(\langle \text{integer} \rangle\) is the number of GKS-spRPA iterations (default is 0).

  • ldiis - Turns on DIIS algorithm which speeds up energy convergence.

  • eigshift \(\langle \text{real}\rangle\) - (optional) Introduces a shift to the non-interacting gap to aid in the convergence of GKS iterations. This may be required for small-gap systems. (Note: Carefully check your results if eigshift is being used. After the final iteration, rerun a rirpa energy calculation without the iter a nd eigshift keywords. The resulting energy should be used as the final GKS-spRPA energy.)

  • output \(\langle \text{string}\rangle\) - (optional) A condensed version of the output relevant to GKS-spRPA energetics is written to the file \(\langle \text{string}\rangle\) (default filename is gksrpa.dat).

• the convergence criterion for total energy is controlled by the $scfconv keyword.