Basis sets for relativistic calculations

Gaussian Type Orbitals (GTOs)

Cartesian Gaussians are defined as

$$G^{\alpha}_{ijk}=G^{\alpha}_i(x)G^{\alpha}_j(y)G^{\alpha}_k(z)$$

with

$$G^{\alpha}(x)=\left(\frac{2\alpha}{\pi}\right)^{1/4}\sqrt{\frac{\left(4\alpha\right)^i}{\left(2i-1\right)!!}} x^i\exp\left[-\alpha x^2\right]$$

Alternatively one may express a Cartesian Gaussian as

$$G^{\alpha}_{ijk}=N^{\alpha}_{ijk}x^iy^jz^k\exp\left[-\alpha r^2\right];\quad N^{\alpha}_{ijk}=\left(\frac{2\alpha}{\pi}\right)^{3/4}\sqrt{\frac{2^l}{F_{ijk}}}\left(\sqrt{2\alpha}\right)^l$$

where the sum $i+j+k$ is associated with orbital angular momentum $l$. The factor

$$F_{ijk}=\left(2i-1\right)!!\left(2j-1\right)!!\left(2k-1\right)!!$$

shows that the Cartesian components of a given shell may have different normalization constants. In the HERMIT integral code a single normalization is chosen for each shell by ignoring $F_{ijk}$ meaning that Cartesian Gaussians are normalized to

$$\begin{split}\left<G^{\alpha}_{ijk}|G^{\alpha}_{ijk}\right>=F_{ijk}\end{split}$$

In practice this means that s- and p-functions are normalized to one. So are $d110$, $d101$ and $d011$, whereas $d200$, $d020$ and $d002$ are normalized to three. $f111$ is normalized to one, $f210$, $f201$, $f120$, $f102$ and $f012$ are normalized to three, and $f300$, $f030$ and $f003$ are normalized to 15.

Spherical Gaussians are defined by

$$G^{\alpha}_{lm} = R^{\alpha}_l\left(r\right)Y_{lm}\left(\theta,\phi\right);\quad$$

where the angular part is given by spherical harmonics $Y_{lm}$ and the radial part by

$$R^{\alpha}_l=N^{\alpha}_lr^l\exp\left[-\alpha r^2\right];\quad N^{\alpha}_l = \frac{2\left(2\alpha\right)^{3/4}}{\pi^{1/4}}\sqrt{\frac{2^l}{\left(2l+1\right)!!}}\left(\sqrt{2\alpha}\right)$$

In passing we note that

$$N^{\alpha}_l=2\sqrt{\frac{\alpha}{2l+1}}N^{\alpha}_{l-1}$$

For given $l$ there are $\frac{1}{2}\left(l+2\right)\left(l+1\right)$ $\left(2l+1\right)$ spherical Gaussians. The latter basis functions therefore provide more compact basis set expansions. However, in 4-component relativistic calculations the use of spherical Gaussians is somewhat more complicated since the coupling of large and small component basis functions needs to be taken into account.

Kinetic balance

Kinetic balance corresponds to the non-relativistic limit of the exact coupling of large and small component basis functions for positive-energy orbitals. Starting from the radial function of a large component spherical Gaussian basis function one obtains

$$R^S\propto -\left[\frac{\partial}{\partial r}+\frac{\left(1+\kappa^L\right)}{r}\right]R^{\alpha}_l = \sqrt{\alpha\left(2l+3\right)}R^{\alpha}_{l+1}-2\sqrt{\frac{\alpha}{2l+1}}R^{\alpha}_{l-1}$$

We can now distinguish two cases

$$\begin{split}R^S \propto \left\{ \begin{array}{ll} \sqrt{\left(2l+3\right)}R^{\alpha}_{l+1}-2\sqrt{\left(2l+1\right)}R^{\alpha}_{l-1} & \mbox{if $\kappa^L = l$};\\ R^{\alpha}_{l+1} & \mbox{if $\kappa^L = -\left(l+1\right)$}.\end{array} \right.\end{split}$$

The case $\kappa^L < 0$ is straightforward, but a bit more care is needed for the implementation for the case $\kappa^L>0$. The modified spherical Gaussian to be constructed is

$$G^*_{lm} = N\left\{\sqrt{\left(2l+1\right)\left(2l-1\right)}R^{\alpha}_l-2\left(2l-1\right)R^{\alpha}_{l-2}\right\}Y_{l-2,m};\quad N = \frac{1}{\sqrt{\left(2l+1\right)\left(2l-1\right)}}$$

Constructing modified spherical harmonics for kinetic balance; the gritty details

We write the spherical harmonic $G^{\alpha}_{l-2,m}$ as a linear combination of Cartesian Gaussians

$$G^{\alpha}_{l-2,m} = \sum_{i+j+k=l-2}c^{l-2,m}_{ijk}G^{\alpha}_{ijk}$$

In the HERMIT integral code we have selected a transformation such that the solid harmonics are normalized to unity. From this we obtain

$$\begin{split}\begin{array}{lcl} G^*_{lm}&=&N\left[4\alpha r^2-2\left(2l-1\right)\right]\sum_{i+j+k=l-2}c^{l-2,m}_{ijk}G^{\alpha}_{ijk}\\ &=&N\sum_{i+j+k=l-2}c^{l-2,m}_{ijk} \left[4\alpha \left( \frac{N^{\alpha}_{ijk}}{N^{\alpha}_{i+2,j,k}}G^{\alpha}_{i+2,j,k} +\frac{N^{\alpha}_{ijk}}{N^{\alpha}_{i,j+2,k}}G^{\alpha}_{i,j+2,k} +\frac{N^{\alpha}_{ijk}}{N^{\alpha}_{i,j,k+2}}G^{\alpha}_{i,j,k+2} \right) -2\left(2l-1\right)G^{\alpha}_{ijk} \right] \end{array}\end{split}$$

The ration of normalization constants in the above relation is given by

$$\frac{N^{\alpha}_{ijk}}{N^{\alpha}_{i+2,j,k}}=\frac{1}{4\alpha}\sqrt{\frac{F_{i+2,j,k}}{F_{ijk}}}$$

However, the expression is much simplified by the fact that factors $F_{ijk}$ area set to one in HERMIT, such that

$$G^*_{lm}=N\sum_{i+j+k=l-2}c^{l-2,m}_{ijk} \left[ \left(G^{\alpha}_{i+2,j,k}+G^{\alpha}_{i,j+2,k}+G^{\alpha}_{i,j,k+2}\right) -2\left(2l-1\right)G^{\alpha}_{ijk} \right]$$