Next: Conclusions Up: First-principles density-functional calculations using Previous: Density-functional calculations

# Results of the calculations

In this section we present and discuss the results obtained from density-functional calculations using the spherical-wave basis set. We study the convergence of the total energy as a function of the cutoff energy, ; the radii of the basis spheres, ; the maximum angular momentum component, ; and the number of basis spheres, . Physical properties are deduced from total energy calculations on the systems. For molecules, we calculate the equilibrium bond lengths and force constants. For bulk crystalline silicon, we calculate the equilibrium lattice parameter and bulk modulus. These results are compared with those obtained using a plane-wave code[18], and from experiment[26].

In Fig. 1 we plot the total energy of the chlorine molecule with a bond length of 1.6 Å as a function of cutoff energy and basis sphere radius . The figure shows that the total energy decreases rapidly as the cutoff energy and the basis sphere radius are increased, which is to be expected from the additional variational freedom that is introduced. Convergence in the total energy is achieved for cutoff energies above 800 eV.

Fig. 2 shows that the rate of convergence of the total energy with respect to the cutoff energy is the same for both the localized spherical-wave and extended plane-wave basis sets. This confirms that the energy cutoff concept can be equally applied in the spherical-wave basis set.

Using an energy cutoff above 900 eV, we calculate the total energy of the chlorine molecule for a variety of bond lengths as a function of the basis sphere radius . Fig. 3 shows that the total energy converges exponentially with respect to . We also note that the total energy converges slightly faster with respect to for molecules with smaller bond lengths. This reflects the fact that for a given , the basis set is more complete for a smaller molecule than a larger one because the basis spheres are closer to one another in the smaller molecule.

Since the total energy also depends on other parameters such as and the number of basis spheres , we have performed calculations on the chlorine molecule with a bond length of 2.4 Å. The results in Fig. 4 show the convergence of the total energy of the system as a function of for different basis-sphere radii and numbers of basis spheres . The rapid convergence of the total energy with respect to is evident from the figure. We note that the best'' result obtained from the spherical-wave calculation with ,  4.50 Å, and gives a total energy of  eV, which lies 0.023 eV above the plane-wave total energy of  eV. This difference, which is due to the incompleteness of the spherical-wave basis set, could be reduced further by increasing the basis-sphere radius and . However, we are content with this accuracy because the error due to the incompleteness of the spherical-wave basis set is only about of the total energy obtained from the plane-wave calculation. The number of spherical-wave basis functions in this case is only 672, which is a small fraction (0.6%) of 112452, the number of plane waves.

To study the effect of , and on the calculated physical properties such as the equilibrium bond length and force constant , we perform a series of calculations on the chlorine molecule for a range of bond lengths from 1.70 Å to 2.25 Å. A typical result is shown in Fig. 5. The results of the calculations of and with 1, 2, and 3 are displayed in Tables 1, 2, and 3, respectively. The errors in and displayed in the columns headed under and are deduced from the results of the plane-wave calculations.

Table 1: Results for the equilibrium bond length and force constant of the chlorine molecule with = 1. When , two basis spheres of the same radius centered on the atoms are used. When , three basis spheres of the same radius are used, the third basis sphere being centered between the atoms. The experimental values for and are 1.9878 Å and 3.23 N/cm, respectively. The plane-wave calculations give values of 1.9661 Å and 3.790 N/cm, respectively, where we have used the same pseudopotentials, Brillouin zone sampling, and cutoff energy.
 = 2 = 3 (Å) (Å) (Å) 2.00 2.0214 2.81 4.880 28.8 1.9565 0.49 4.492 18.5 2.50 2.1240 8.03 2.542 32.9 1.9988 1.66 3.747 1.1 3.00 2.1536 9.54 1.874 50.6 2.0098 2.22 3.557 6.1 3.50 2.1595 9.84 1.772 53.2 2.0123 2.35 3.496 7.8 4.00 2.1623 9.98 1.703 55.1 2.0127 2.37 3.484 8.1 4.50 2.1625 9.99 1.701 55.1 2.0128 2.38 3.481 8.2

In Table 1, the values of and converge rapidly with respect to . However, the results with two basis spheres and shows that the converged results contain unacceptably large systematic errors. The inclusion of a third sphere reduces the errors significantly because the bonding region between the atoms is described better by the third sphere. The results show it is impossible to improve the results simply by enlarging when is used.

Table 2: Results for the equilibrium bond length and force constant of the chlorine molecule with = 2. The meanings of and are explained in the caption of Table 1.
 = 2 = 3 (Å) (Å) (Å) 2.00 1.8941 3.66 5.276 39.2 1.9158 2.56 4.591 21.1 2.50 1.9609 0.26 3.959 4.5 1.9537 0.63 4.017 6.0 3.00 1.9733 0.37 3.754 0.9 1.9636 0.13 3.849 1.6 3.50 1.9769 0.55 3.694 2.5 1.9658 0.02 3.802 0.3 4.00 1.9777 0.59 3.673 3.1 1.9663 0.01 3.789 0.0 4.50 1.9778 0.60 3.670 3.2 1.9663 0.01 3.786 0.1

We repeat the calculations for and with , for which the results are presented in Table 2. The converged results with and are better than the converged results with and , which indicates the importance of over for the minimal basis set'' calculations. With and , the error of the converged results for and are 0.50% and 13.6% compared to the experimental values, respectively. These accuracies are acceptable within the LDA.

Table 3: Results for the equilibrium bond length and force constant of the chlorine molecule with = 3.
 (Å) (Å) 2.00 1.8833 4.21 5.212 37.5 2.50 1.9481 0.92 4.128 8.9 3.00 1.9636 0.13 3.860 1.8 3.50 1.9668 0.04 3.792 0.1 4.00 1.9674 0.07 3.777 0.3 4.50 1.9675 0.07 3.773 0.4

Finally in Table 3, we present the values of and using two basis spheres centered on the atoms . As expected, the converged values of and agree very well with the plane-wave results. We note that calculations with are expensive, since the number of basis functions is almost double that for .

Next we calculate the total energy of hydrogen molecule with a bond length of 1.0 Å as a function of the cutoff energy , and the basis sphere radius , for which the results are displayed in Fig. 6. The total energy converges rather slowly with respect to the cutoff energy because a bare Coulomb potential due to the hydrogen atom is used. Such behavior is also observed in the plane-wave calculations. However, the convergence of energy differences is achieved when the cutoff energy exceeds 800 eV.

We perform a series of total energy calculations on the hydrogen molecule for a range of bond lengths to determine the values of and . The results are tabulated in Table 4 and show that we can use a value of as small as 3.00 Å to obtain an accuracy of less than 1% in and with only two basis spheres. This should be contrasted with the case of the chlorine molecule where with , , and  Å, the values of and agree only fortuitously with the plane-wave results.

Table 4: Results for the equilibrium bond length and force constant of the hydrogen molecule with = 2. The meanings of and are explained in the caption of Table 1. Cut-off energies above 1000 eV are used. The experimental values for and are 0.7414 Å and 5.75 N/cm respectively. The equivalent plane-wave calculations give values of 0.7711 Å and 5.197 N/cm, respectively.
 = 2 = 3 (Å) (Å) (Å) 2.00 0.7476 3.05 5.998 15.4 0.7503 2.70 5.865 12.9 2.50 0.7643 0.88 5.420 4.3 0.7668 0.56 5.232 0.7 3.00 0.7695 0.21 5.237 0.8 0.7712 0.01 5.198 0.0 3.50 0.7709 0.03 5.193 0.1 4.00 0.7710 0.01 5.178 0.4

We can explain why, to obtain the same accuracy, the chlorine molecule requires a larger than the hydrogen molecule. The equilibrium bond length of the hydrogen molecule (which is about 0.74 Å) is smaller than the equilibrium bond length of the chlorine molecule (which is about 1.99 Å). The bonding region between the hydrogen atoms is thus described better by the basis functions because the basis spheres are closer to one another. The hydrogen molecule is also smaller'' (in the sense of the extent of the charge distribution) than the chlorine molecule.

In Fig. 7 we show the total energy of the hydrogen chloride molecule with a bond length of 1.60 Å as a function of the cutoff energy and the radius of the basis sphere. The energy differences converge when the cutoff energy exceeds 800 eV. Calculations are performed to obtain and , and the results are tabulated in Table 5.

Table 5: Results for the equilibrium bond length and force constant of the hydrogen chloride molecule with = 2. The meanings of and are explained in the caption of Table 1. For the hydrogen chloride molecule, the experimental values for and are 1.2746 Å and 5.16 N/cm, respectively. The equivalent plane-wave calculations give values of 1.2948 Å and 5.458 N/cm.
 = 2 = 3 (Å) (Å) (Å) 2.00 1.2601 2.68 6.351 16.4 1.2712 1.82 6.030 10.5 2.50 1.2862 0.66 5.683 4.1 1.2885 0.49 5.605 2.7 3.00 1.2926 0.17 5.518 1.1 1.2933 0.12 5.489 0.6 3.50 1.2936 0.09 5.475 0.3 1.2941 0.05 5.461 0.1 4.00 1.2942 0.05 5.466 0.1 1.2946 0.02 5.452 0.1 4.50 1.2945 0.02 5.464 0.1 1.2948 0.00 5.450 0.1

We repeat the and calculations for the hydrogen chloride molecule, where the radius of the basis sphere centered on the chlorine atom is fixed at 4.00 Å but the radius of the basis sphere centered on the hydrogen atom is varied. The results are presented in Table 6, which shows that we can use a smaller basis sphere of a radius of 2.0 Å centered on the hydrogen atom to obtain an accuracy of less than 1%. It is thus possible to use different basis spheres depending on the atomic species, which is important because this can reduce the computation time significantly.

Table 6: Results for and of the hydrogen chloride molecule. Two basis spheres are used. The radius of the basis sphere centered on the chlorine atom is fixed at 4.00 Å but the radius of the basis sphere centered on the hydrogen atom is varied.
 (Å) (Å) 2.00 1.2928 0.15 5.508 0.9 2.50 1.2937 0.08 5.488 0.5 3.00 1.2943 0.04 5.472 0.3 3.50 1.2939 0.07 5.468 0.2 4.00 1.2942 0.05 5.466 0.1

Fig. 8 shows the total energy of the silane molecule with a Si-H bond length of 1.83 Å, as a function of and . Total energy differences converge for cutoff energies above 800 eV. The results of the calculations of and (for the breathing mode) are summarized in Table 7. We find that the accuracy is acceptable when 3.00 Å.

Table 7: Results for and of the silane molecule. Five basis spheres of the same radius centered on the atoms are used. The cubic simulation cell has sides of length 12 Å. . The experimental value of is 1.4798 Å, while the equivalent plane-wave calculations give the values of 1.4910 Å and 11.38 N/cm for and , respectively.
 (Å) (Å) 2.00 1.4402 3.41 14.640 28.6 2.50 1.4811 0.66 12.010 5.5 3.00 1.4893 0.11 11.523 1.3 3.50 1.4906 0.03 11.415 0.3 4.00 1.4912 0.01 11.382 0.0 4.50 1.4914 0.03 11.379 0.0

We repeat the and calculations on the silane molecule with the radius of the basis sphere centered on the silicon atom fixed at 4.00 Å, but with the radius of the basis spheres centered on the hydrogen atoms varied. The results in Table 8 show that an accuracy of 1% can be achieved by using  Å, which is 1 Å larger than the basis spheres centered on the hydrogen atom in the hydrogen chloride molecule calculation (c.f. Table 6).

Table 8: Results for and of the silane molecule. Five basis spheres centered on the atoms are used. The radius of the basis sphere centered on the silicon atom is fixed at 4.00 Å but the radius of basis spheres centered on the hydrogen atoms is varied.
 (Å) (Å) 2.00 1.4851 0.40 11.847 4.1 2.50 1.4856 0.36 11.948 5.0 3.00 1.4902 0.05 11.438 0.5 3.50 1.4906 0.03 11.419 0.3 4.00 1.4912 0.01 11.382 0.0

From the pseudo-charge density of the hydrogen chloride molecule (Fig. 9), we observe that the valence electrons are concentrated towards the chlorine atom, as expected from the relative electronegativites of hydrogen and chlorine. This enables us to use a smaller basis sphere centered on the hydrogen atom to obtain accurate results. However, from the pseudo-charge density of the silane (Fig. 10), we observe substantial charge density around the hydrogen atoms, reflecting the fact that hydrogen is more electronegative than silicon. Hence for the silane molecule calculations, the radius of the basis spheres centered on the hydrogen atoms need to be larger than that for hydrogen chloride. These observations lead to the conclusion that the relative electronegativities of neighboring atoms in a calculation should be taken into account when choosing basis sphere radii.

We have chosen bulk crystalline silicon to test the performance of the basis set on an extended system. Fig. 11 shows the total energy per atom for a 64-atom cell of silicon with a lattice parameter of 5.43 Å as a function of the cutoff energy and . The total energy converges at a cutoff energy of about 250 eV. The rapid convergence of the total energy with respect to is evident from the figure.

To determine the equilibrium lattice parameter, , and the bulk modulus, , we perform a series of calculations on the bulk silicon system for a range of lattice parameters from 5.31 Å to 5.51 Å. The results of the calculations with 1 and 2 are tabulated in Tables 9 and 10, respectively. It is found that even with , the results with  Å agree quite well with the plane-wave and experiment results. The calculations with improve the results slightly. The reason why calculations give rather good results is because silicon atoms mix the and states to form four orbitals which are obviously well-described by a basis set with .

Table 9: Results for the equilibrium lattice parameter and bulk modulus of the 64-atom bulk crystalline silicon, with . 64 basis spheres of the same radius centered on the atoms are used. The experimental values for and are 5.43 Å and 100.0 GPa, respectively. The plane-wave calculations, with a cutoff energy of 250 eV, give the results of 5.395 Å and 92.3 GPa, respectively.
 (Å) (Å) (GPa) 2.60 5.353 0.78 139.1 50.7 2.80 5.413 0.33 104.5 13.2 3.00 5.445 0.93 97.0 5.1 3.20 5.453 1.08 89.8 2.7

Table 10: Results for the equilibrium lattice parameter and bulk modulus of the 64-atom bulk crystalline silicon, with . 64 basis spheres of the same radius centered on the atoms are used.
 (Å) (Å) (GPa) 2.60 5.310 1.58 238.0 157.9 2.80 5.353 0.78 129.8 40.6 3.00 5.377 0.33 111.4 20.7 3.20 5.385 0.19 96.2 4.2

Finally we present Table 11 which shows the numbers of basis functions for the spherical-wave and plane-wave basis set calculations. Since in general the number of spherical-wave basis functions is very small for molecules compared to that of plane-wave basis functions, we conclude that spherical-wave basis sets can be used to study molecules and possibly clusters with high efficiency.

Table 11: Numbers of the basis functions required so that the agreement between the results for and (for molecules); and and (for bulk silicon) from the spherical-wave and plane-wave basis set calculations is less than 1%. The numbers of the spherical-wave and plane-wave basis functions are denoted by and , respectively. For the spherical-wave basis set calculations, the choice of atom-centered and is used. is the cutoff energy.
 System (eV) Chlorine molecule 800 234 88663 Hydrogen molecule 1000 270 124097 HCl molecule 800 288 88663 Silane molecule 800 720 88663 Bulk silicon 250 4608 10827

Next: Conclusions Up: First-principles density-functional calculations using Previous: Density-functional calculations
Peter D. Haynes 2002-10-31