Results

We present some illustrative examples of the importance of kinetic energy preconditioning for the convergence of calculations with our method, described in more detail in Ref. [16]. In all test cases we use norm-conserving pseudopotentials in Kleinman-Bylander [38] form, the local-density approximation [39,40] for the exchange and correlation term, and the $\Gamma$-point only for the $k$-point sampling.

A silane molecule is placed in a cubic simulation cell of side length $40 \ a_{0}$, with a grid spacing $0.5 \ a_{0}$ (corresponding to a plane wave cut-off of 537 eV) in each direction. The orbitals are initialized to atom-centered fireballs [41] which are strictly localized within spheres of radius $6.0 \ a_{0}$. Each orbital is allowed to vary freely within its localization region. There is one orbital on each hydrogen atom and four on the silicon. In Fig. 2 we show the convergence of the total energy as a function of iteration number. The effect of using different fixed values of the kinetic energy preconditioning parameter $k_{0}$ may be seen. The limit $k_{0} = \infty$ corresponds to the case of no preconditioning. It can be seen that improved performance is achieved for a range of values of $k_{0}$.

Figure 2: Convergence of the total energy as a function of the iteration number for the calculation on a silane molecule. The grid spacing was $0.5 \ a_{0}$ and the localized orbitals were confined within spheres of radius $6.0 \ a_{0}$. $E_{0}$ is the converged value of the total energy for each run, and $k_{0}$ is given in units of $a_{0}^{-1}$.

In Figs. 3 and 4 we show the convergence of the total energy as a function of iteration number for different grid spacings and localization radii, respectively. For the calculations presented in these two figures we used a kinetic energy preconditioning parameter $k_{0}=4.0 \ a^{-1}_{0}$. As the grid spacing is reduced, or the localization radius increased, the size of the basis set and the number of variational parameters in the minimization increases. From Figs. 3 and 4 it is clear that the preconditioning scheme is working well as the number of iterations required to reach a given accuracy does not vary a great deal with the size of the problem. For instance, in Fig. 3, we see that the calculation with a grid spacing of $1.0 \ a_{0}$ (134 eV) reaches an energy convergence of $10^{-6}$ Ha after 11 iterations, whilst with a grid spacing of $0.4 \ a_{0}$ (839 eV), i.e., almost 16 times as many basis functions, the same level of convergence is achieved in just 14 iterations.

Figure 3: Convergence of the total energy as a function of the iteration number for the calculation on a silane molecule. The localized orbitals were confined within spheres of radius $6.0 \ a_{0}$ and kinetic energy preconditioning with $k_{0}=4.0 \ a^{-1}_{0}$ was used. $E_{0}$ is the converged value of the total energy for each run.

Figure 4: Convergence of the total energy as a function of the iteration number for the calculation on a silane molecule. The grid spacing was $0.5 \ a_{0}$ and kinetic energy preconditioning with $k_{0}=4.0 \ a^{-1}_{0}$ was used. $E_{0}$ is the converged value of the total energy for each run.

Finally, as an alternative to preconditioning the Fourier transformed gradient $\tilde{g}^{\beta}_{p}$ by multiplying it with the preconditioner $\xi_{p}$ in reciprocal space, in accord with Eq. (44), we have developed a real space implementation of the preconditioning scheme. In this we convolve the real space gradient $g^{\alpha}_{i}$ with the inverse fast Fourier transform (FFT) of $\xi_{p}$. Of course, a full convolution would be costly: if the gradient and preconditioner are both of size $N_{\mathrm{grad}}$, then the computational effort required to perform a full convolution scales as $N^{2}_{\mathrm{grad}}$. Thus, we truncate the preconditioner in real space at a radial cut-off $R_{0}$ so that it is nonzero over only a small number of points $N_{\mathrm{prec}} \ll N_{\mathrm{grad}}$. The computational cost of performing a convolution between the gradient $g^{\alpha}_{i}$ and this truncated preconditioner is much more favourable and scales as $N_{\mathrm{prec}} N_{\mathrm{grad}}$. For typical values of $N_{\mathrm{grad}}$ and $N_{\mathrm{prec}}$, this is comparable to the cost of preconditioning in reciprocal space.

Truncating the preconditioner in real space is not simply a matter of improving the computational efficiency, for it also makes physical sense: the reason behind preconditioning is to smear out large kinetic energy variations over short distances, thus a convolution that is localized in real space over just a few grid points is all that should be required. This is demonstrated by the results presented in Fig. 5, which shows the convergence of the total energy with this real space scheme for the above-introduced silane molecule. The different curves correspond to various radial cut-offs $R_{0}$ for the inverse FFT of the preconditioning function $\xi_{p}$. Comparing Figs. 2 and 5 we see that preconditioning via local convolution in real space is as successful as the conventional reciprocal space approach, and that there is little sensitivity to the choice of the cut-off radius $R_{0}$, which may be as small as $1.0 \ a_{0}$.

Figure 5: Convergence of the total energy as a function of the iteration number for the calculation on a silane molecule. The grid spacing was $0.5 \ a_{0}$ and the localized orbitals were confined within spheres of radius $6.0 \ a_{0}$. The top curve is for the case of no preconditioning ( $k_{0} = \infty$), while for the others $k_{0} = 3.0 \ a_{0}^{-1}$. $R_{0}$ is the convolution radius in real space. $E_{0}$ is the converged value of the total energy for each run.