17.7 Random-phase approximation

The random-phase approximation program (`rpatddft`) can be used to calculate RPA
correlation energies after a SCF calculation.
Additionnally, it can be used to calculate dynamic dipole polarizabilities,
C dispersion coefficients, and excitation energies.
The program currently works without point-group symmetry.

List of the main keywords:

`ECORR, <list of methods>`-

Calculation of RPA correlation energies [1] (see options below) `PROPERTIES, <list of methods>`-

Calculation of dynamic dipole polarizabilities and C dispersion coefficients [10] (see options below) `EXCIT, <list of methods>`-

Calculation of excitation energies [11] (see options below)

as well as:

`ORB,<orbrec>`- Record for input orbitals (required).
`STAB`- Check matrices stability conditions in RPA calculations.
When used without an
`ECORR`,`EXCIT`or`PROPERTIES`keyword, check the Hessian and RPA matrices eigenvalues and do nothing more. `INTAC,NLAMBDA=<n>[,LAMBDA=<lambda>,WEIGHT=<weight>]`-

Number of quadrature points for the Gauss-Legendre numerical integration along the adiabatic connection for RPA calculations (default is 7). If`LAMBDA`and`WEIGHT`are given, assumes a one-point quadrature with given abscissa and weight. `INTFREQ,METHFREQ=<methfrq>,NFREQ=<nfreq>`-

Options for the numerical integration over the frequency variable of RPA calculations.`METHFREQ`governs the type of quadrature used (0(default) is Gauss-Chebyshev, 1 and 2 are Gauss-Legendre, 3 is Clenshaw-Curtis) and`NFREQ`governs the number of quadrature points (default is 16). `INTEGRAL,<nbr>`- Specify the two-electron integral transformation routine: 0 (still the default for spin-unrestricted and gradient calculations) is the `old' one, 1 is the `old' one that has been cleaned up, and 2 (default otherwise) is a much more efficient transformation using Molpro's “transform” routine.
`DFTKERNEL,<funcx>[,<funcc>]`-

Specify the exchange and correlation kernel for`EXCIT`(if only one argument is given, it is understood to be the exchange-correlation kernel). `NOMP2`- The MP2 energy is calculated in certain situations where it is
available almost for free, provided that some matrices are allocated.
This behavior can be switched off by this
`NOMP2`keyword. `DIELSOLID`- Use the solid-state variant when performing DIEL calculations.
`NOSPINBLOCK`- For spin-unrestricted calculations, use a formalism where matrices are of dimensions (the default is to use a formalism with a nospinflip/splinflip block structure)
`NOSPINFLIP`- Exclude spin-flip dimensions of unrestricted RPA calculations
that use the
`NOSPINBLOCK`formalism (not suitable for all RPA variants). `C6`- Computes C6 coefficients from last two saved polarizabilities.
`TDA`- Tamm-Dancoff approximation for
`EXCIT`and`PROPERTIES`. `OCC,<nocc>`- Explicitly specify the number of occupied orbitals (useful for fake pseudopotential calculations).
`CORE,<core>`- Specify core orbitals (default: last specified core orbitals or, if none, atomic inner shells)
`PRINT,<nbr>`- Level of print expected from the output (from 0(default) to 3).

Calculation of RPA correlation energies
`ECORR, <list of methods>`

If no method is given, a `SO2-RCCD` calculation will be done (see below).

There are two main RPA variants [1]: dRPA (direct RPA, without the inclusion of an Hartree-Fock exchange kernel in the response function) and RPAx (with the Hartree-Fock exchange kernel included in the response function).

There are four main formulations in which the RPA equations can be derived. The adiabatic-connection fluctuation-dissipation theorem (ACFDT) equation involves integrations both over frequency and coupling constant: an analytical integration along the frequency variable followed by a quadrature on the coupling constant yields the adiabatic connection formulation (AC) [1], while an analytic integration on the coupling constant followed by a numerical integration over the frequency yields the dielectric formulation (DIEL) [2]. Two other formulations are obtained when the two integrations are carried out analytically: the plasmon formula (PLASMON) [1] and the ring coupled cluster doubles formulation (RCCD) [3]. These four formulations are not in general equivalent.

Most variants+formulations can readily by used in a spin-unrestricted context [6]. This is implemented in the code and does not need any further input from the user: the RPA program recognizes the spin-unrestricted character of a SCF calculation that was done beforehand and acts accordingly.

Gradients of most of the RCCD-formulation RPA energies are available, both without range-separation
with RHF orbitals and with range-separation with RSH orbitals [9].
The calculations are triggered by the presence of the keyword `FORCE` or `OPTG` after
the energy-related section (see examples at the end of the section).

The user can test the RPA program using `make rpatddfttest`,
which proposes a variety of tests for RPA correlation energy calculations.

The keywords for the methods are constructed on the model:

For the AC formulation, the available methods are:

`DRPAI-AC`- dRPA calculation (see Refs. [1] and [4]).
`DRPAII-AC`- dRPA calculation, using antisymmetrized two-electron integrals (see Refs. [1] and [4]).
`RPAXI-AC`- RPAx calculation (see Refs. [1] and [5]).
`RPAXII-AC`- RPAx calculation, using antisymmetrized two-electron integrals (see Refs. [1] and [5]).
`DRPAI-AC-ALT`- dRPA calculation using an alternative derivation (see Ref. [7]).
`RPAXI-AC-ALT`- RPAx calculation using an alternative derivation (see Ref. [7]).
`DRPAI-AC-NOINT`- dRPA calculation without integration along the adiabatic connection (using the “kinetic” and “potential” contributions, sometimes called the “alternative plasmon formula”, see Ref. [1]).
`RPAXII-AC-NOINT`- RPAx calculation without integration along the adiabatic connection (using the “kinetic” and “potential” contributions, sometimes called the “alternative plasmon formula”, see Ref. [1]).

For the DIEL formulation, the available methods are:

`DRPAI-DIEL`- dRPA calculation (see Ref. [2]).
`DRPAIIA-DIEL`- dRPA-IIa approximation (see Ref. [2]).
`RPAXIA-DIEL`- RPAx-Ia approximation (see Ref. [2]).

For the RCCD formulation, the available methods are:

`DRPAI-RCCD`- dRPA-I calculation (see Ref. [3]).
`RPAXII-RCCD`- RPAx-II calculation (Szabo-Ostlund variant 1 is calculated too, see Ref. [3]).
`SO2-RCCD`- Szabo-Ostlund variant 2 (see Ref. [3]).
`SOSEX-RCCD`- dRPA-I+SOSEX correction (see Ref. [3]).
`RPAX2-RCCD`- RPAX2 approximation (see Ref. [8]).

For the PLASMON formulation, the available methods are:

`DRPAI-PLASMON`- dRPA-I calculation (see Ref. [1]).
`RPAXII-PLASMON`- RPAx-II calculation (see Ref. [1]).

Note that to all these keywords are associated energy variables defined as :

Example of a dRPA-I calculation using the PBE functional:

{rks,pbe,orbital,2101.2} {rhf,start=2101.2,maxit=0} {rpatddft; orb,2101.2; ecorr,DRPAI-AC } e=ECORR_DRPAI_AC

Example of a range-separated RPAx-I calculation using the short-range PBE exchange-correlation functional and the range-separated parameter mu=0.5:

{int;erf,0.5} {rks,exerfpbe,ecerfpbe;rangehybrid;orbital,2101.2} {rpatddft; orb,2101.2; ecorr,RPAXI-AC }

Example of several RPA calculations in the same run:

{rhf,orbital,2101.2} {rpatddft; orb,2101.2; ecorr,DRPAI-AC,RPAXII-RCCD,DRPAI-DIEL } e1=ECORR_DRPAI_AC e2=ECORR_RPAXII_RCCD e3=ECORR_DRPAI_DIEL(the calculations are done with the same transformed integrals, i.e. without redoing the integral transformation).

Example of a dRPA-I gradient calculation:

{rks,pbe,orbital,2101.2} {rpatddft; orb,2101.2; ecorr,DRPAI-RCCD } force

Example of a geometry optimization at the LDA+dRPA-I level:

{int;erf,0.5} {rks,exerf,ecerf;rangehybrid;orbital,2101.2} {rpatddft; orb,2101.2; ecorr,DRPAI-RCCD } optg

Calculation of properties, excitation energies and oscillator strengths
`PROPERTIES, METHOD=<method>`

`EXCIT, METHOD=<method>`
The `EXCIT` calculations output shows the excitation energies in ua, eV and nm, the oscillator strengths in length and velocity gauge, as well as the major excitations involved in each mode.
The methods available are:

`DRPA`- Direct random-phase approximation (or time-dependent Hartree).
`TDHF`- Time-dependent Hartree-Fock.
`TDDFT`- Time-dependent density-functional theory.
`RS-TDDFT`- Range-separated time-dependent density-functional theory [11].

The exchange density functionals (FUNCX) available are:

`LDAXERF`- (short-range LDA exchange density functional for the erf interaction [12]).

The correlation density functionals (FUNCC) available are:

`LDAC`- (Perdew-Wang-92 LDA correlation density functional)
`LDACERF`- (short-range LDA correlation density functional for the erf interaction [13]).

Example of a range-separated time-dependent density-functional theory calculation using the short-range LDA exchange-correlation functional and the range-separated parameter mu=0.5:

mu=0.5 {int;erf,mu;save} {rks,exerf,ecerf;rangehybrid;orbital,2101.2} {int} {setmu,mu} {rpatddft; orb,2101.2; excit,method=rs-tddft; dftkernel,funcx=ldaxerf,funcc=ldacerf }

References

J. G. Ángyán, R.-F. Liu, J. Toulouse, and G. Jansen, J. Chem. Theory Comput. **7**, 3116 (2011).

G. Jansen, B. Mussard, D. Rocca, J. G. Ángyán (in prep).

J. Toulouse, W. Zhu, A. Savin, G. Jansen, and J. G. Ángyán, J. Chem. Phys. **135**, 084119 (2011).

F. Furche, Phys. Rev. B **64**, 195120 (2001).

J. Toulouse, I. C. Gerber, G. Jansen, A. Savin, and J. G. Ángyán, Phys. Rev. Lett. **102**, 096404 (2009).

B. Mussard, P. Reinhardt, J. G. Ángyán, and J. Toulouse, J. Chem. Phys. (submitted).

Heßelmann, A., Görling, A., Phys. Rev. Lett. **106**, 093001 (2011).

Heßelmann, A., Phys. Rev. A **85**, 012517 (2012).

B. Mussard, P. G. Szalay, J. G. Ángyán, J. Chem. Theory Comput. **10**, 1968 (2014).

J. Toulouse, E. Rebolini, T. Gould, J. F. Dobson, P. Seal, J. G. Ángyán, J. Chem. Phys. **138**, 194106 (2013).

E. Rebolini, A. Savin, J. Toulouse, Mol. Phys. **111**, 1219 (2013).

J. Toulouse, A. Savin, and H.-J. Flad, Int. J. Quantum Chem. **100**, 1047 (2004).

S. Paziani, S. Moroni, P. Gori-Giorgi, and G. B. Bachelet, Phys. Rev. B **73**, 155111 (2006).

molpro@molpro.net 2019-03-24