next up previous
Next: 6. Non-local pseudopotential Up: Localised spherical-wave basis set Previous: 4. Overlap matrix elements

5. Kinetic energy matrix elements

The kinetic energy matrix elements for any two basis functions $\chi_{n \ell
m}^{\alpha}$ and $\chi_{n' \ell' m'}^{\beta}$ centred at ${\bf R}_{\alpha}$ and ${\bf R}_{\beta}$ respectively are defined by
$\displaystyle T_{\alpha \beta}$ $\textstyle =$ $\displaystyle -{\textstyle{1 \over 2}}\int_{\mathrm{all~space}} \mathrm{d}^3 r ...
...\bf R}_{\alpha}) \nabla^2 \chi_{n' \ell'
m'}^{\beta}({\bf r} - {\bf R}_{\beta})$  
  $\textstyle =$ $\displaystyle {1 \over 2(2
\pi)^3} \int \mathrm{d}^3 k~k^2 \mathrm{e}^{-{\mathr...
...hi}_{n \ell m}^{\alpha}({\bf k}) {\tilde
\chi}_{n' \ell' m'}^{\beta}({\bf k}) .$ (26)

Because of the discontinuity in the first derivatives of the basis functions at the sphere boundaries, a delta-function arises when the Laplacian operates on a basis function. This is integrated out when the matrix element is calculated and this contribution is included when transforming the real-space integral to reciprocal-space in equation (26).

The second line of equation (26) is identical to equation (11) apart from a factor of ${\textstyle{1 \over 2}}k^2$. The same separation into individually regular terms can be applied here, and the result is that we need to calculate the contour integral (17) as before, except that the integer $p$ must be replaced by $(p-2)$ and a numerical factor of ${\textstyle{1 \over 2}}$ is introduced. The calculation of the residues is identical to that presented in the previous section, except that the integrand no longer always has a pole at $z = 0$ in every term.

The results for $T_{\alpha \beta}$ when ${\bf R}_{\alpha \beta} = 0$ are

    $\displaystyle \frac{ {\textstyle{1 \over 2}}\delta_{\ell \ell'} \delta_{m m'}}{...
... r_{\alpha} \geq r_{\beta}
\end{array} \right\} q_{n \ell} \not= q_{n' \ell'} ,$  
      (27)
    $\displaystyle \textstyle{1 \over 4} \delta_{\ell \ell'} \delta_{m m'} q_{n
\ell...
...quad r_{\alpha} \geq r_{\beta}
\end{array} \right\} q_{n \ell} = q_{n' \ell'} .$  

The calculation of the kinetic energy has been checked by projecting a set of wave functions expanded in the spherical-wave basis onto the plane-wave basis using equation (9a). As the kinetic energy cut-off for the plane-wave basis is increased, so the description of the wave functions becomes more accurate. The kinetic energy calculated using the results above can then be compared against the kinetic energy calculated by a plane-wave $O(N^3)$ code [2].

From the asymptotic behaviour of the spherical Bessel functions, the Fourier transform (9a) for large $k$ is

\begin{displaymath}
{\tilde \chi}_{n \ell m}^{\alpha}({\bf k}) \sim \frac{\sin(k...
...- {\ell \pi \over 2})}{k^3}~{\bar Y}_{\ell m} (\Omega_{\bf k})
\end{displaymath} (28)

and so the error in the kinetic energy due to truncating the plane-wave basis with cut-off $E_{\mathrm{cut}} = {\textstyle{1 \over 2}}
k_{\mathrm{cut}}^2$ is
\begin{displaymath}
\Delta T \sim \int_{k_{\mathrm{cut}}}^{\infty} \mathrm{d}k~k...
...c{1}{k_{\mathrm{cut}}} \sim \frac{1}{\surd E_{\mathrm{cut}}} .
\end{displaymath} (29)

In figure 1 the kinetic energy as calculated by the plane-wave code has been plotted against $1 / \surd E_{\mathrm{cut}}$ and yields a straight line as expected, which can then be extrapolated to obtain an estimate of the kinetic energy calculated for infinite cut-off: $60.66 \pm 0.01$ eV. This is in agreement with the value calculated analytically of $60.65$ eV.

Figure 1: Plot of asymptotic fit to kinetic energy data.
\includegraphics [height=100mm,angle=-90]{keplot.ps}


next up previous
Next: 6. Non-local pseudopotential Up: Localised spherical-wave basis set Previous: 4. Overlap matrix elements
Peter Haynes