text:
# Generalized Relativistic Effective Core Potential Method: Theory and calculations

## I Introduction

The Dirac-Coulomb (DC) Hamiltonian which is used for calculation of heavy atoms and heavy-atom molecules has the form ($`e=m_e=\mathrm{}=1`$):

$$𝐇^{DC}=\underset{i}{}𝐡^D(i)+\underset{i>j}{}\frac{1}{r_{ij}},$$
(1)
where the one-electron Dirac operator $`𝐡^D`$ is

$`𝐡^D=ıc(\stackrel{}{\alpha }\stackrel{}{})+c^2(\beta 1)+V_{nuc},`$

$`\stackrel{}{\alpha },\beta `$ are the $`4\times 4`$ Dirac matrices, $`c`$ is the speed of light and $`V_{nuc}`$ is the nuclear potential which can take account of effect of finite nuclear size etc.

Higher approximation levels involve modifying the instantaneous Coulomb interaction between electrons. These modifications are derived from the Quantum ElectroDynamic (QED) theory. The simplest approximation is obtained by including the one-photon exchange. In the Coulomb gauge it leads to so-called Dirac-Coulomb-Breit (DCB) Hamiltonian (see, e.g., )

$$𝐇^{DCB}=𝐇^{DC}+B_{ij},$$
(2)
where

$`B_{ij}(\omega _{ij})={\displaystyle \underset{i>j}{}}[(\stackrel{}{\alpha }_i\stackrel{}{\alpha }_j){\displaystyle \frac{\mathrm{cos}(\omega _{ij}r_{ij})}{r_{ij}}}+(\stackrel{}{\alpha }_i\stackrel{}{}_i)(\stackrel{}{\alpha }_j\stackrel{}{}_j){\displaystyle \frac{\mathrm{cos}(\omega _{ij}r_{ij})1}{\omega _{ij}^2r_{ij}}}],`$

and $`\omega _{ij}`$ denotes the transition frequency between the electrons $`i`$ and $`j`$.

A low-frequency expansion of the cosines yields the incomplete Breit interaction $`B_{ij}(0)`$:

$`B_{ij}(0)={\displaystyle \frac{1}{2}}[\stackrel{}{\alpha }_i\stackrel{}{\alpha }_j+(\stackrel{}{\alpha }_i\stackrel{}{r}_{ij})(\stackrel{}{\alpha }_j\stackrel{}{r}_{ij})/r_{ij}^2]/r_{ij}.`$

These terms describe the instantaneous magnetostatic interaction and retardation in the electric interaction between electrons.

In calculations on heavy-atom molecules, the DC and DCB Hamiltonians are usually replaced by an effective Hamiltonian

$$𝐇^{Ef}=\underset{i_w}{}[𝐡^{Schr}(i_w)+𝐔^{Ef}(i_w)]+\underset{i_w>j_w}{}\frac{1}{r_{i_wj_w}},$$
(3)
written only for valence and some outer core electrons denoted by indices $`i_w`$ and $`j_w`$; $`𝐔^{Ef}`$ is an Relativistic Effective Core Potential (RECP) operator simulating interactions of the explicitly treated (“pseudo-valence”) electrons with those which are excluded from the RECP calculation. In Eq. (3), $`𝐡^{Schr}`$ is the nonrelativistic one-electron Schrödinger operator

$`𝐡^{Schr}={\displaystyle \frac{1}{2}}\stackrel{}{}^2+V_{nuc},`$

Contrary to the four-component (relativistic) wave function used in DC(B) calculations, the pseudo-wave function in the RECP case can be both two- and one-component.

## II Features of shape-consistent RECPs

The RECP approximations with the radially local operator for “shape-consistent” (or “norm-conserving”) pseudoorbitals (pseudospinors)  are the most widely used in modern calculations of molecules with heavy elements. Below we discuss mainly these RECP versions and their developments because they provide the largest computational savings for the same level of accuracy as compared with other RECP versions.

### A Advantages

* Chemically inactive electrons are excluded from RECP calculations.
* The valence orbitals (spinors) are smoothed in the core regions of heavy atoms to generate nodeless pseudoorbitals (pseudospinors). Therefore, the number of the one-electron basis functions may be minimized, thus reducing dramatically both the number of two-electron integrals and the computational time.
* The small components of four-component spinors are eliminated and the nonrelativistic kinetic energy operator is used. The RECP method allows one to use a well-developed nonrelativistic technique of calculation and relativistic effects are taken into account with the help of spin-dependent semilocal potentials. As a result, the most part of difficulties inherent for the DC(B) molecular calculations can be avoided.
* In principle, correlations of the explicitly treated electrons with those which are excluded from the RECP calculation can be considered within “correlated” RECP versions. Reducing the number of explicitly correlated electrons with the help of the correlated RECPs is a very promising way to minimize efforts when performing high-precision molecular calculations.

### B Disadvantages

* By now, different schemes of the RECP generation are suggested which use the radially local operator for the effective potential. However, as it was demonstrated in many test calculations, they provide the comparable level of accuracy for the same number of the explicitly treated electrons.

It is clear that the explicit inclusion of the outer core electrons into the RECP calculation is the way to increase the computational accuracy. However, the extension of the space of these electrons more than some limit does not improve the accuracy as is obtained in all our atomic calculation with RECPs and PPs. Although the errors caused by the nodeless RECP approximations can be small enough in competent treatments, they still range up to 2000 $`cm^1`$ and more for the dissociation and transition energies even for lowest-lying excitations that can be unsatisfactory in many applications.

One can conclude that the scheme of smoothing and the number of the explicitly treated electrons are not only responsible for refining the accuracy of the RECP calculations.
* The reliability of the radially local RECP versions is not high for transitions with excitations from the outer core shells, for the chemistry of transition, rare-earth elements, etc.
* Moreover, the direct calculation of such properties as electronic densities near heavy nuclei, hyperfine structure, and matrix elements of other operators singular on heavy nuclei is impossible as a result of the smoothing of the orbitals in the core regions of heavy elements.

To overcome these disadvantages of the shape-consistent RECP method, a few developments were done by the authors concerning the RECP operator first of all . A detailed theoretical analysis of the RECP method is presented in .

### C Historical background of RECP method

Some of the most important achievements in the development of the effective core potential (or pseudopotential) theory are listed below:

* The first ideas and papers devoted to the PseudoPotential (PP) theory of Hellmann (1935)  and Gombas (1935)  were published more than 60 years ago.
* The procedure of treatment of the valence electrons in the frozen sea of core electrons was suggested by Fock et al.  in 1940.
* The smoothing of orbitals in atomic cores was suggested and the orbital angular momentum ($`l`$) dependence of the PP operator for the case of pseudoorbitals smoothed in the cores was first emphasized in the pioneering paper of Phillips & Kleinman (1959) .
* The PP operator with the angular projectors was used first by Heine and Abarenkov (1964, 1965)  to generate one-electron model potentials mainly for the solid-state calculations.
* Improvement of the PP accuracy when treating several valence electrons explicitly was suggested by Weeks & Rice  in 1968.
* Goddard III suggested a scheme (1968)  for the generation of the radially local effective potentials on the base of inverting the atomic Hartree-Fock (HF) equations.
* The ab initio radially local effective potentials were efficiently applied by Kahn et al. (1976)  to molecular calculations.
* The nonrelativistic shape-consistent ECP was suggested by Durand & Barthelat (1975) . Christiansen et al. (1977)  have modified their smoothing scheme.
* The Relativistic ECP method for two-component atomic pseudospinors was proposed by Lee et al. (1977) .
* Hafner & Schwarz (1979)  have split the RECP operator on the spin-independent and spin-orbit parts that allowed one to take into account the spin-orbit interaction only at the final stage of a calculation thus reducing computational efforts.

### D Notations

Below we shall use the following notations for the principal quantum numbers $`n`$ of orbitals: $`n_f`$ for the inner core (IC) orbitals, $`n_c`$ for the outer core (OC) orbitals, $`n_v`$ for the valence (V) orbitals, and $`n_a`$ for the virtual (A) orbitals. The quantum numbers of the orbital momentum, total momentum, and its projection will be designated by the $`l`$, $`j`$, and $`m`$ indices, respectively. We shall write one-electron states as

$`\stackrel{~}{\psi }_{nljm}(\tau )=\stackrel{~}{\phi }_{nlj}(r)\chi _{ljm}(\mathrm{\Omega },\sigma )`$

for the pseudospinors and as

$`\psi _{nljm}(\tau )=\left(\begin{array}{cc}P_{nlj}(r)\hfill & \chi _{ljm}(\mathrm{\Omega },\sigma )\hfill \\ Q_{nlj}(r)\hfill & i\chi _{l^{}jm}(\mathrm{\Omega },\sigma )\hfill \end{array}\right)`$

for the four-component spinors where

$`l^{}=2jl,\tau (\stackrel{}{r},\sigma )(r,\mathrm{\Omega },\sigma ),`$

$`\stackrel{~}{\phi }_{nlj}`$ are the radial parts of the pseudospinors, $`P_{nlj}`$ and $`Q_{nlj}`$ are the radial parts of the large and small components of the Dirac spinors, and $`\chi _{ljm}`$ is the two-component spin-angular function.

We shall assign an orbital to the OC or V subspaces using the average radius criterion rather than the orbital energy criterion. The difference in such a definition is essential mainly for transition metals, lanthanides and actinides in which orbitals with different principal quantum numbers may have close orbital energies. For the problems of the RECP generation and application, which we are going to discuss, the space criterion is more appropriate.

## III GRECP operator

Following basically the scheme developed by K. Pitzer’s group (Lee et al. , Christiansen et al. ), the numerical pseudospinors $`\stackrel{~}{\phi }_{nlj}(r)`$ are constructed of the large components $`P_{nlj}(r)`$ of the outer core and valence Dirac-Fock (DF) spinors so that the innermost pseudospinors of them (for each $`l`$ and $`j`$) are nodeless, the next pseudospinors have one node, and so forth:

$$\stackrel{~}{\phi }_{nlj}(r)=\{\begin{array}{cc}P_{nlj}(r),\hfill & rR_c,\hfill \\ f(r)=r^\gamma _{i=0}^5a_ir^i,\hfill & r<R_c,\hfill \end{array}$$
(4)
$`\begin{array}{ccccc}l=0,1,\mathrm{},L,& & j=|l\pm \frac{1}{2}|,& & n=n_c,n_c^{},\mathrm{},n_v,\end{array}`$

where $`L`$ is one more than the highest orbital angular momentum of the IC spinors. The leading power $`\gamma `$ in the polynomial is typically chosen to be close to $`L`$ in order to ensure a sufficient ejection of the pseudovalence electrons from the IC region.

To derive the GRECP components $`U_{nlj}(r)`$, the HF equations are inverted for the V and OC pseudospinors (for each $`l=0,\mathrm{},L`$ and $`j=|l\pm \frac{1}{2}|`$) so that $`\stackrel{~}{\phi }_{nlj}`$ are solutions of the nonrelativistic-type HF equations in the jj-coupling scheme for a “pseudoatom” with the removed IC electrons (Goddard III, 1968)

$`U_{nlj}(r)`$ $`=`$ $`\stackrel{~}{\phi }_{nlj}^1(r)[({\displaystyle \frac{1}{2}}{\displaystyle \frac{d^2}{dr^2}}+{\displaystyle \frac{1}{r}}{\displaystyle \frac{d}{dr}}{\displaystyle \frac{l(l+1)}{2r^2}}+{\displaystyle \frac{Z}{r}}\stackrel{~}{𝐉}(r)+\stackrel{~}{𝐊}(r)`$ (5)
$`+`$ $`\epsilon _{nlj})\stackrel{~}{\phi }_{nlj}(r)+{\displaystyle \underset{n^{}n}{}}\epsilon _{n^{}nlj}\stackrel{~}{\phi }_{n^{}lj}(r)],`$ (6)
where $`Z`$ is the nuclear charge, $`\stackrel{~}{𝐉}`$ and $`\stackrel{~}{𝐊}`$ are the Coulomb and exchange operators calculated with the pseudospinors, $`\epsilon _{nlj}`$ is the one-electron energy of the corresponding spinor, and $`\epsilon _{n^{}nlj}`$ are the off-diagonal Lagrange multipliers.

In general, different potentials $`U_{nlj}`$ are obtained for different pseudospinors $`\stackrel{~}{\phi }_{nlj}`$ with the help of Eq. (6). It is necessary to construct such an RECP operator which would act on each individual pseudospinor as the corresponding potential. This operator may be presented in the non-Hermitean (right-handed) form as:

$$𝐔^{Ef}=\underset{n}{}\underset{l=0}{\overset{\mathrm{}}{}}\underset{j=|l1/2|}{\overset{l+1/2}{}}U_{nlj}(r)\stackrel{~}{𝐏}_{nlj},$$
(7)
where

$`\stackrel{~}{𝐏}_{nlj}={\displaystyle \underset{m=j}{\overset{j}{}}}\left|\stackrel{~}{nljm}\stackrel{~}{nljm}\right|,`$

$`n=n_c,n_c^{},\mathrm{},n_v,n_a,n_a^{},\mathrm{},`$

i.e., $`|\stackrel{~}{nljm}\stackrel{~}{nljm}|`$ is the projector on the $`\stackrel{~}{\psi }_{nljm}`$ pseudospinor and the $`n`$ index runs over the OC, V, and A pseudospinor subspaces. The direct use of expression (7) in calculations is impossible because it includes the infinite summations on the $`n`$ and $`l`$ indices.

The $`U_{nlj}`$ potentials contain the contribution from the Coulomb and exchange interaction with the IC electrons, the contribution from smoothing the spinors, and the contribution from relativistic effects as a result of the replacement the Dirac Hamiltonian by the nonrelativistic-type Hamiltonian.

With the help of the property

$$\stackrel{~}{\phi }_{nlj}\left|U_{nlj}U_{n^{}lj}\right|\stackrel{~}{\phi }_{n^{}lj}=0,$$
(8)
when neglecting small terms, the operator (7) (Generalized RECP or GRECP operator) can be written in the Hermitian form:

$`𝐔`$ $``$ $`U_{n_vLJ}(r)+{\displaystyle \underset{l=0}{\overset{L}{}}}{\displaystyle \underset{j=|l1/2|}{\overset{l+1/2}{}}}\left[U_{n_vlj}(r)U_{n_vLJ}(r)\right]𝐏_{lj}`$ (9)
$`+`$ $`{\displaystyle \underset{n_c}{}}{\displaystyle \underset{l=0}{\overset{L}{}}}{\displaystyle \underset{j=|l1/2|}{\overset{l+1/2}{}}}\left\{\left[U_{n_clj}(r)U_{n_vlj}(r)\right]\stackrel{~}{𝐏}_{n_clj}+\stackrel{~}{𝐏}_{n_clj}\left[U_{n_clj}(r)U_{n_vlj}(r)\right]\right\}`$ (10)
$``$ $`{\displaystyle \underset{n_c,n_c^{}}{}}{\displaystyle \underset{l=0}{\overset{L}{}}}{\displaystyle \underset{j=|l1/2|}{\overset{l+1/2}{}}}\stackrel{~}{𝐏}_{n_clj}\left[{\displaystyle \frac{U_{n_clj}(r)+U_{n_c^{}lj}(r)}{2}}U_{n_vlj}(r)\right]\stackrel{~}{𝐏}_{n_c^{}lj},`$ (11)
where

$`𝐏_{lj}={\displaystyle \underset{m=j}{\overset{j}{}}}\left|ljmljm\right|`$

is the ordinary spin-angular projector and the first line in Eq. (11) presents the standard radially local RECP operator. It can be shown  that the action of this operator on the OC, V, and A pseudospinors is equivalent (with the required accuracy) to the action of the corresponding potentials on them.

Note that the nonlocal terms with the projectors on the most important correlating functions $`\stackrel{~}{\phi }_{n_kl_kj}^{corr}(r)`$ localized mainly in the OC and V regions and with the corresponding potentials, $`U_{n_kl_kj}(r)`$, in principle (i.e. in the cases in which it is essential), can be taken into account in the expression for the GRECP operator additionally to the terms with the OC projectors $`\stackrel{~}{𝐏}_{n_clj}`$.

With the help of the following identities for the $`𝐏_{lj}`$ projectors (Hafner & Schwarz, 1979) :

$$𝐏_{l,j=l\pm 1/2}(\mathrm{\Omega },\sigma )=\frac{1}{2l+1}\left[\left(l+\frac{1}{2}\pm \frac{1}{2}\right)𝐏_l(\mathrm{\Omega })\pm 2𝐏_l(\mathrm{\Omega })\stackrel{}{𝐥}\stackrel{}{𝐬}𝐏_l(\mathrm{\Omega })\right],$$
(12)
the GRECP operator can be easily written in the spin-orbit representation.

### A Calculations of Hg and Pb

The accuracy of the GRECP approximation was tested in correlation structure calculations of the Hg and Pb atoms in which we were able to use very flexible basis sets. The all-electron Fock-space Relativistic Coupled Cluster (RCC) method has been described in review  and its references. RECP errors in reproducing the all-electron RCC transition energies in Hg for the case of 20 correlated electrons and the equivalent correlation $`[7,9,8,6,7,7]`$ basis sets (for all-electron and RECP calculations, see ) are presented in Table I (compiled from ). The highest absolute error in reproducing the transition energies with excitation or ionization of single electron is 94 cm-1 for the GRECP , 729 cm-1 for the RECP of Ross et al.  and 1747 cm-1 for the energy-adjusted PP . The same number of electrons, twenty, is explicitly treated in all these RECP versions. The larger errors for RECPs  and  are mainly due to the neglect of the difference between the outer core and valence potentials in these RECP versions (see  for more details).

Calculations of transition energies in the Pb atom are performed for five lowest-lying states in  (see Table II here) by the four-component  and two-component (with GRECP)  versions of the combined method of the Configuration Interaction (CI) and the second-order Many-Body Perturbation Theory (MBPT2) . The 22e-GRECP version for Pb  was tested and very good agreement of the GRECP/MBPT2/CI calculations with the four-component DC/MBPT2/CI results is obtained.

## IV Frozen core approximation for outer core shells

To perform precise calculations of the chemical and spectroscopic properties, the correlations should be taken into account not only within the valence regions of heavy atoms and heavy-atom molecules but within the core regions and between the valence and core electrons as well. In practice, the goal is to achieve a given level of accuracy correlating as small a number of electrons as possible, thus reducing the computational effort. However, the performance of the RECPs generated for a given number of explicitly treated electrons can not always satisfy the accuracy requirements of the correlation structure calculation.

It can be illustrated on example of Tl . To attain a good accuracy, e.g., within 400 $`cm^1`$ for a group of low-lying states of molecules containing Tl, one should correlate at least 13 electrons of Tl (i.e., including the $`5d`$ shell) and it would be optimal to use the RECPs with this number of electrons treated explicitly (i.e., 13e-RECP) like the RECP of Ross et al.  or our valence RECP version . However, all the known nodeless 13e-RECPs cannot ensure the above-mentioned accuracy and one should use the RECPs with at least 21 electrons (e.g., 21e-GRECP for Tl of our group ). The $`5s`$ and $`5p`$ pseudospinors can be treated as frozen, providing, nevertheless, the above-declared accuracy. The $`5p`$ orbitals have the orbital energies about four times higher and their average radii are 1.4 times shorter than those for the $`5d`$ orbitals. Moreover, their angular correlations are suppressed as compared with those for the $`5d`$ shell because the most “profitable” polarization states ($`5d`$ in this case) are completely occupied in the lowest-lying states. Therefore, the $`5s,5p`$ states are substantially less active in chemical processes.

One can apply the energy level-shift technique in order to “freeze” the $`5s`$ and $`5p`$ pseudospinors . Following Huzinaga et al. , one should add the matrix elements of the SCF field operators (the Coulomb and spin-dependent exchange terms) with these OC pseudospinors to the one-electron part of the Hamiltonian together with the level-shift terms

$$\underset{n_{c_f},l,j}{}B_{n_{c_f}lj}\underset{m}{}\stackrel{~}{|n_{c_f}l_{c_f}jm}\stackrel{~}{n_{c_f}l_{c_f}jm|},$$
(13)
where $`B_{n_{c_f}lj}`$ is at least of order $`|2\epsilon _{n_{c_f}lj}|`$ of magnitude and $`\epsilon _{n_{c_f}lj}`$ are the orbital energies of the OC pseudospinors $`\stackrel{~}{\varphi }_{n_{c_f}lj}(r)`$ which are frozen. These nonlocal terms are needed in order to exclude a collapse of the molecular orbitals to the “frozen” states (the $`5s_{1/2},5p_{1/2,3/2}`$ pseudospinors for Tl).

All terms with the frozen core pseudospinors described here (the Coulomb and exchange interactions, the level-shift operator) can easily be presented in the spin-orbit form with the help of Eq. (12), as was said above with respect to the GRECP operator (the spin-orbit form of the GRECP operator can be found, e.g., in ). More importantly, these OC pseudospinors can be frozen (as spinors) in calculations with the spin-orbit basis sets and they can be frozen at the stage of calculation of the one-electron matrix elements of the Hamiltonian (as is implemented in the MOLGEP code ). Thus, any integrals with the indices of the frozen spinors are completely excluded already after the integral calculation step.

### A Calculation of TlH

The Multi-Reference single- and Double-excitation Configuration Interaction (MRD-CI) method  with the Spin-Orbit (SO) configuration selection scheme  was applied in the first precise GRECP calculation of molecules. In the GRECP/MRD-CI calculations of spectroscopic properties for the $`0^+`$ ground state of TlH  (see Table III here), the 21e/8fs-GRECP for Tl was used (i.e., the 21 electron GRECP  with 8 electrons occupying the frozen OC pseudospinors $`5s_{1/2}`$ and $`5p_{1/2,3/2}`$). The contracted  basis set for thallium generated in  and the  basis for hydrogen (see http://www.qchem.pnpi.spb.ru) contracted from the primitive (6,3,1) gaussian basis set of Dunning  were used.

Fourteen electrons are correlated in the calculation of spectroscopic constants in TlH and, as one can see from Table III, very good agreement is found with the experimental data contrary to other known calculations. Close results are obtained very recently in the GRECP/RCC calculations of TlH in our joint calculations with group of Prof. Kaldor (Tel-Aviv University). The GRECP/RCC and GRECP/MRD-CI calculations of HgH and HgH+ are in progress now.

## V Correlations with inner core electrons and Breit effects

On example of the mercury atom  the importance of correlation of the inner core electrons (removed from RECP calculations) was investigated with the help of the combined MBPT2/CI method  and with the RCC method . It is obtained that at least 34 external electrons (occupying the IC $`4f`$ shell, outer core $`5s,5p`$ and $`5d`$ shells, valence $`6s,6p`$ etc. shells) of Hg should be correlated and the one-electron basis set should contain up to $`h`$ angular momentum functions in order to attain a reliable agreement with the experimental data for the low excitation energies within 100–200 $`cm^1`$, whereas the errors of the gaussian approximation of the GRECP components and the effects of different nuclear models are negligible for this accuracy. Otherwise, e.g., for the 20 electron RECPs for Hg known from the literature, energies of excitations from the ground state can not be calculated with the accuracy better than 200–500 $`cm^1`$.

However, our test calculations show that the main contribution from the correlations with the $`4f`$ shell is due to the one-electron correction from the self-energy diagrams  and, therefore, this contribution can be taken into account for the 20 electron GRECP at the generation stage, i.e. the correlations with the IC $`4f`$ shell can be described by the GRECP. The technique of generating the “correlated” GRECP is discussed by the authors in . When using such GRECPs in calculations, both required basis sets for the corresponding atoms and computational efforts can be seriously reduced when performing correlation structure calculations.

Turning back to the TlH calculation, further improvement of the accuracy can be attained when correlations with the outer core $`4f,5p`$ and $`5s`$ shells of Tl and Breit effects are taken into account. We expect that this can be efficiently done in the framework of the “correlated” 21e/8fs-GRECP version in which 13 electrons are treated explicitly as in the above considered GRECP/MRD-CI calculation of TlH.

The first correlated GRECP we have generated for the Hg atom in the framework of the MBPT2 method. Despite the considerable improvement in accuracy with this GRECP version when calculating excitation energies in Hg, the reliability of the applied generation scheme was found to be not very high. A new scheme of generation of correlated GRECPs in the framework of the RCC method is in progress now.

Breit and other QED effects are relatively large for lanthanides and actinides. The incorporation of these corrections to calculations with RECPs is discussed in our theoretical paper on the GRECP method . Some recent test calculations have shown that the Breit and other QED corrections for outer core and valence electrons can be efficiently represented by a one-electron spin-dependent radially local operator and, thus, can be easily incorporated into GRECPs.

## VI RECPs for transition metals, lanthanides, and actinides

The disadvantage of the commonly used RECP versions for the transition, rare-earth, etc., elements is that they still require the explicit treatment of essentially larger number of electrons as compared with other (“normal”) elements in order to attain the same computational accuracy.

In paper , we added “self-consistent” (SfC) terms to the RECP operator which take account of the circumstance that the RECPs generated for the different occupation numbers of the outermost $`d`$ or $`f`$ shells of these elements are somewhat distinguished. Below we will refer to these shells as OuterMost Core (OMC) shells following the partitioning made in subsection II D. The significant improvement was attained in reproducing the atomic excitation energies in comparison with the conventional RECP calculations.

Whereas the form of the standard radially local RECP operator accounts for the property that the dependence of the corresponding potentials $`U_{n_vlj}(r)`$ on the relaxation of the V spinors ($`\phi _{n_vlj}`$) is weak in the bond making, the SfC correction to this operator takes into account the fact of a relatively small relaxation of the OC spinors $`\phi _{n_clj}`$ caused by the change of the OMC shell occupation number in a majority of cases.

The main features of the SfC RECP constructing for the case of the radially local RECP are:

1. The all-electron DF calculations of two generator states with different occupation numbers of the OMC $`d`$ ($`f`$) shell are carried out for an atom under consideration. These occupation numbers are designated $`N_1`$ and $`N_2`$.
2. The “generating” RECPs $`𝐔^{N_1}`$ and $`𝐔^{N_2}`$ are constructed for these generator states employing the standard scheme. The operator for a radially local RECP (incorporating the SfC corrections into the GRECP is made similarly) has the form

$$𝐔^N=E_{core}^N+U_{LJ}^N(r)+\underset{l=0}{\overset{L}{}}\underset{j=|l\frac{1}{2}|}{\overset{l+\frac{1}{2}}{}}[U_{lj}^N(r)U_{LJ}^N(r)]𝐏_{lj},$$
(14)
where $`E_{core}^N`$ is the core energy and $`U_{lj}^N(r)`$ are the radial components of the RECP for each pair of quantum numbers $`l`$ and $`j`$. They are derived for the OMC $`d`$ ($`f`$) shell occupation number $`N`$.

This operator differs from the standard (“nodeless”) RECP operator by the $`E_{core}^N`$ constant. This value corresponds to the contribution of the core electrons (excluded from calculations with the conventional RECP) to the total atomic energy and it is usually excluded from the RECP calculations because the $`E_{core}^N`$ only shifts all the energy levels by a constant and does not influence on the chemical properties. In the case of the SfC RECP, inequality $`E_{core}^{N_1}E_{core}^{N_2}`$ should be taken into account when transitions with the change of the occupation number $`N`$ are under consideration.
3. The one-electron SfC RECP operator is written as

$$𝐔^{SfC}=\frac{𝐔^{N_1}+𝐔^{N_2}}{2}+\frac{𝐔^{N_1}𝐔^{N_2}}{N_1N_2}\left(N\frac{N_1+N_2}{2}\right)+B\left(N\frac{N_1+N_2}{2}\right)^2,$$
(15)
where $`B`$ is the adjusting parameter and $`N`$ is the average value $`\mathrm{\Psi }𝐍\mathrm{\Psi }`$ of the OMC $`d`$ ($`f`$) shell occupation number operator which in the second quantization representation is

$$𝐍=\underset{jm}{}\stackrel{~}{𝐚}_{n_{omc}l_{omc}jm}^{}\stackrel{~}{𝐚}_{n_{omc}l_{omc}jm},$$
(16)
where $`n_{omc}`$ and $`l_{omc}`$ are the principal and orbital quantum numbers of the OMC shell.

At every iteration of the MC SCF, etc., procedures, the $`N`$ value can be determined for the wave function $`\mathrm{\Psi }`$ without a serious complication of the calculation.

We should add here that in a high-precision electronic structure calculation of a heavy-atom molecule, correlation with the electrons occupying the IC or/and frozen OC states (see section IV) can be substantial for calculated properties and dependent from the occupation numbers of states of V and other OC electrons as was demonstrated in  for the OC-OC correlation of $`5s`$ and $`5d`$ electrons in Hg. Therefore, in principle, the SfC terms dependent on the occupation numbers of the V states can appear in the expression for the SfC RECP operator (15). In this case, the use of the “correlated” RECPs $`𝐔^{N_1}`$ and $`𝐔^{N_2}`$ also can be required in (15).

### A SfC GRECP calculations of uranium

In Table IV (compiled from ), one-configuration numerical calculations of transition and ionization energies in uranium are presented for our GRECP versions, for the RECP of Ermler et al.  and for the energy-adjusted pseudopotential of Küchle et al.  as compared with the DF energies.

The GRECP correcting terms described in Section III allows one to describe accurately the V region of uranium when the $`5f`$ shell occupation number close to that for the generator state. The addition of the SfC correction increases the accuracy of reproducing the energies of transitions without changing the $`5f`$ shell occupation number in two to eight times and in those with changing these occupation numbers in 1.5–3 times. One can see from Table IV that the addition of the quadratic SfC (QSfC) correction (the last term in (15)) leads to a significant increasing (up to 20 times) of the SfC GRECP accuracy for the transitions with the change of the OMC shell occupation number.

## VII Term-splitting correction

As one can see from Table V, the above-discussed SfC corrections give no improvement in the term splitting of the $`5f_{5/2}^36d_{3/2}^17s_{1/2}^2`$ configuration of uranium as compared with the generating GRECPs. Analysis of these errors shows that the main contribution (about 90 %) is due to reproducing the OMC spinors with the help of the nodeless pseudospinors.

Consider the GRECPs with 24 and 32 explicitly treated electrons for which the OMC $`5f`$ pseudospinors are nodeless and the $`6d`$ ones have radial nodes (for the RECPs with the 14-electron pseudovalence space, for which these OMC shells are described by the nodeless pseudospinors, the errors are greater in general and the analysis is more complicated).

The simplest way to reduce drastically these errors is to generate and use such GRECP versions for which the $`5f`$ shell is characterized by nodal pseudospinors, whereas the $`4f`$ pseudospinors are nodeless. To reduce computational efforts, the $`4f`$ shell can be treated as frozen in the way described in section IV. However, when approximating the nodal $`5f`$ pseudospinors, some additional basis functions and computational efforts will be required.

If the small magnitude of the OMC shell ($`5f`$ here) relaxation is taken into account, there is another way out which can be applied for the cases of low-lying states. We suggest to add the terms

$`𝐔^{TSc}`$ $`=`$ $`{\displaystyle \underset{x_1,x_2,x_3,x_4}{}}\lambda _{x_1x_2,x_3x_4}\stackrel{~}{|x_1}\stackrel{~}{x_2|}\stackrel{~}{|x_3}\stackrel{~}{x_4|}`$ (17)
$``$ $`2{\displaystyle \underset{w}{}}{\displaystyle \underset{x_1,x_2,x_3}{}}(\lambda _{x_1x_2,x_3x_3}\lambda _{x_1x_3,x_3x_2})\delta _{wx_3}\stackrel{~}{|x_1}\stackrel{~}{x_2|},`$ (18)
where $`x_i=(n_{omc},l_{omc},j_i,m_i)`$ *** The $`j_i`$ and $`m_i`$ indices run over all the possible values for a given ($`n_{omc},l_{omc}`$). In the case of the partially occupied OMC $`5f`$ shell, the $`x_i`$ indices correspond to both OC and virtual pseudospinors (i.e., occupied and unoccupied one-electron states in the generator state). , to the GRECP operator. These terms are correcting the one- and two-electron integrals with the $`5f`$ pseudospinors (and only these integrals).

The results of our calculations (Table V) show that the addition of the Term-Splitting (TS) correction (called in  by the SO correction) allows one to reduce the most serious errors about 10 times for the term splittings, thus reducing the errors for the transition energies between terms to the same order of magnitude as are the errors for transitions between the states averaged over the different configurations (Table IV).

## VIII Electronic densities and singular operators in atomic cores

To evaluate the matrix elements of the operators being singular in the core regions of heavy atoms in a molecule after the RECP calculation of this molecule at the first step, the proper shapes of the valence four-component molecular spinors must be restored in the atomic core regions at the second step.

The applicability of the two-step algorithms to the electronic structure calculation of heavy-atom molecules is a consequence of the fact that the valence and core electron subspaces may be represented as the subsystems, interaction between which is described mainly by some integrated and not by detailed properties of these subsystems. As a result, the use of the GRECP allows one to reproduce accurately the original electronic densities in the valence region for excited states as it is demonstrated on U in  and on Hg, Tl, Pb, and Bi in . Another result is that the electronic structure can be restored in the core regions of heavy atoms after the GRECP calculation of a molecule with a good accuracy.

Both the nonvariational procedures (Pacios & Christiansen, 1985; Titov, 1986)  and the variational technique (Titov, 1996)  for the electronic structure restoration after the shape-consistent RECP calculation were suggested.

In , the nonvariational restoration technique was applied to evaluation of the hyperfine constants and the parity non-conservation effects in the PbF, HgF, BaF, and YbF molecules. The molecular pseudoorbitals in the form of linear combinations of atomic pseudoorbitals were evaluated at the RECP calculation stage. Then, the atomic pseudoorbitals were replaced by the unsmoothed four-component DF spinors derived for the same atomic configurations as the basis pseudoorbitals. The MO LCAO coefficients were preserved after the RECP calculation. Such constructed molecular relativistic spinors were used for calculation of the P,T-odd spin-rotational Hamiltonian parameters.

In a general case, the first-order reduced density matrix, $`\{\stackrel{~}{D_{pq}}\}`$, should be evaluated after the molecular RECP calculation on the basis set of atomic pseudospinors whereas the matrix elements $`\{A_{pq}^K\}`$ of the one-electron operators of hyperfine interaction etc., $`\{𝐀^K\}`$, should be evaluated in the equivalent basis set of atomic four-component spinors. Then the physical values are calculated as

$$𝐀^K=\underset{pq}{}\stackrel{~}{D_{pq}}A_{pq}^K.$$
(19)
The latest ab initio calculations of the hyperfine, $`P`$-odd, and $`P,T`$-odd constants for the BaF and YbF molecules are performed by our group with the help of the GRECP/RASSCF/EO scheme  (RASSCF is the Restricted Active Space SCF method  and EO is Effective Operator technique based on MBPT2 ). This technique allows one to take into account polarization and correlation in the valence and core regions. The electronic wave function of BaF and YbF for the ground ($`{}_{}{}^{2}\mathrm{\Sigma }_{1/2}^{}`$) states is calculated in the GRECP approximation. The molecular four-component spinors in the core region of heavy atoms are restored in the framework of the non-variational procedure. Core polarization/correlation effects are taken into account with the help of the atomic MBPT2 when constructing the effective valence Hamiltonians for the Ba and Yb atoms.

### A Calculations of YbF

The results for the YbF molecule by Mosyagin et al. (1998)  are presented in Table VI compared with other results. For the isotropic hyperfine constant $`A=(A_{}+2A_{})/3`$, the accuracy of our calculation is about 3% as compared to the experimental datum. The dipole constant $`A_\mathrm{d}=(A_{}A_{})/3`$ (which is much smaller in magnitude), though better than in all previous calculations known from the literature, is still underestimated by almost 23%. Being corrected within a semiempirical approach for a perturbation of $`4f`$-shell in the core of Yb due to the bond making, this error is reduced to 8%. Our value for the effective electric field on the unpaired electron is 4.9 a.u.=$`2.5\times 10^{10}`$ V cm-1.

One can see from Table VI that the values of the $`W_d`$ constant

$`W_dd_\mathrm{e}=2^2\mathrm{\Sigma }_{1/2}|H_d|^2\mathrm{\Sigma }_{1/2},`$

where $`H_d`$ describes interaction of the electron EDM (i.e. Electric Dipole Moment) $`d_\mathrm{e}`$ with the internal molecular field $`𝐄^{\mathrm{mol}}`$:

$`H_d=2d_\mathrm{e}\left(\begin{array}{cc}0& 0\\ 0& \sigma \end{array}\right)𝐄^{\mathrm{mol}},`$ (22)
from the unrestricted DF calculation of Parpia (1998) , the most recent semiempirical calculation of Kozlov (1997)  and our latest GRECP/RASSCF/EO calculation are in a very close agreement now. It should be noted that the valence electron contribution to the $`W_d`$ in  is in 7.4% agreement with the corresponding RECP-based calculation of Titov et al. (1996) . Another recent DF calculation of Quiney et al. (1998)  gives the value, which is two times smaller.

It is clear from the calculation that the $`4f`$ electrons in the core of ytterbium should be explicitly treated in the planned GRECP calculation of YbF by the relativistic coupled cluster method in order to obtain the parameters of $`P,T`$-odd spin-rotational Hamiltonian with better accuracy.

### B Variational restoration

In the variational technique of the restoration , the proper behavior of the four-component molecular spinors in the core regions of heavy atoms may be restored as an expansion on the spherical harmonics inside the sphere with an atomic core radius $`R_{rest}`$. The outer parts of spinors are frozen after the RECP calculation of the molecule considered. This method enables one to combine the advantages of two well-developed approaches: molecular RECP calculations (with gaussian basis sets) and one-center calculations (of atomic type with, e.g., numerical functions) in the most optimal way.

The most promising application of the two-step method is in a possibility to “split” the correlation structure calculation of a molecule on two consequent calculations in the valence and core regions when the OC and V electrons are explicitly involved into the GRECP calculation and the OC and IC space regions are treated at the one-center restoration stage.

Roughly speaking, the computational efforts for the correlation structure calculations in the core and valence regions are added within the two-step technique, whereas in the ordinary one-step scheme, they are multiplied by each other.

The efficiency of the two-step technique was confirmed in the BaF and YbF calculations . Even these not yet perfect calculations have provided a high accuracy (5–10 %) for the P,T-odd spin-rotational Hamiltonian parameters (which can be compared with the 50 % error for the hyperfine constants in YbF with the correlations in the V region taken into account within RASSCF), thus confirming a considerable promise of the two-step scheme for the correlation structure calculations.

###### Acknowledgements.

This work was supported by the INTAS grant N 96-1266, the DFG/RFBR grant N 96–03–00069, and the RFBR grant N 99–03–33249. The most consuming molecular calculations were carried out at the computer centers of the Bergische Universität GH Wuppertal and Tel Aviv university with the help of the JECS codes .