Results and discussion

We now present some results of calculations we have performed using the methods described above. We first demonstrate the accuracy of the FFT box technique as compared to using the entire simulation cell as the FFT grid. Figure 4 shows the quantity $\Delta E_{1}$, defined as

$$\Delta E_{1} \equiv E^{\mathrm{box}}[n] - E[n],$$ (46)

for the molecules silane ( $\mathrm{SiH_{4}}$) and cyclohexane ( $\mathrm{C_{6}H_{12}}$), for different FFT box sizes. For these tests we used a cubic simulation cell of side length $50 a_{0}$, although any other lattice symmetry could easily have been employed. The PAOs are confined within spherical regions of radius $6.0 a_{0}$, and vanish exactly at the region boundary [20]. They were generated using norm-conserving pseudopotentials [26], the local density approximation (LDA), and an energy cut-off of 538 eV (corresponding to our grid spacing of $0.5 a_{0}$). In the case of silane the silicon atom had one 3s and three 3p PAOs, and for cyclohexane the carbon atoms had one 2s and three 2p PAOs. The hydrogen atoms all had a single 1s PAO.

Figure 4: Top panel: $\Delta E_{1}$ plotted for a silane molecule as a function of FFT box size. Bottom panel: the same quantity plotted for cyclohexane. In both cases all PAOs were confined to atom centered localisation regions with radius $6.0 a_{0}$, and the grid spacing was $0.5 a_{0}$.

We see that the error associated with the FFT box technique are of the order of a few $\mu\mathrm{Hartree}$ per atom, which is comparable to the inherent errors present in plane-wave DFT. We note that the convergence of the total energy with FFT box size is not expected to be variational, and indeed it is not: the FFT box technique should be viewed as a good approximation to the `correct' result that would be obtained by using the entire simulation cell to perform FFTs for calculating the total energy. For a given FFT box size, however, the kinetic energy cut-off of our basis functions (and hence the grid spacing) is a variational parameter, just as in traditional plane-wave DFT.

In addition, we look at a feature common to all real space grid methods, namely the variation of the energy as the system is translated by fractions of a grid spacing with respect to the simulation cell. We define the quantity $\Delta E_{2}$ to be

$$\Delta E_{2}(d - d_{0}) \equiv E^{\mathrm{box}}(d) - E^{\mathrm{box}}(d_{0}),$$ (47)

where $d_{0}$ is a reference position of the system, and $d$ is some rigid displacement with respect to the grid. The top panel in Figure 5 shows this variation, in steps of $0.1a_{0}$, for a lithium hydride (LiH) molecule on a cubic grid with grid spacing $0.5 a_{0}$. The PAOs of both lithium and hydrogen had localisation radii of $10.0a_{0}$, and were generated in the same way as described above. $\Delta E_{2}$ can be seen to be of the order of $10^{-6} \mathrm{Hartree}$. The bottom panel shows the same quantity, but this time with optimisation of the NGWF coefficients, $\{ C_{KLM,\alpha} \}$, in addition to the elements of the density kernel [27]. We see that the variation in total energy is reduced to the $10^{-7} \mathrm{Hartree}$ level. We expect the optimisation of the NGWF coefficients to be important for performing accurate geometry optimisations.

Figure 5: Top panel: Variation in total energy, $\Delta E_{2}$, as the LiH molecule is translated by fractions of a grid spacing. The PAOs on both the lithium and hydrogen atoms had a radius of $10.0a_{0}$. Bottom panel: The same quantity, but with optimisation of the PAO coefficients [27].