For Project 1, the lattice constants for Hf are calculated. Since Hf is an hcp metal with c/a=1.58, we need to calculate the energy by changing the lattice constant c while keeping c/a fixed. And the optimized a, c can be obtained by minimizing the free energy.
1. The hcp structure of Hf
Figure 1 shows the unit cell of Hf hcp structure. OA, OB, OC are three lattice vectors. The lengths and angles of the lattices vectors (OA, OB, OC) are a, b, c and α, β, γ, respectively. To keep the symmetry, the constraints for the unit cell are a=b, α=β=90 degrees, γ=120 degrees.
2. The parameters for the calculations
Some important parameters and inputs for the calculation are listed as following:
Atomic structure for Hf:1s2 2s2 2p6 3s2 3p6 3d10 4s2 4p6 4d10 4f14 5s2 5p6 5d2 6s2
The core configuration of the pseudopotential: 1s2 2s2 2p6 3s2 3p6 3d10 4s2 4p6 4d10
Functional: GGA Perdew Burke Ernzerhof (PBA) functional
Pseudopotential: OTFG ultrasoft
Smearing width: 0.1
Origin shift of K points: 0, 0, 0
The k points grid: investigated in the next section
The energy cutoff: investigated in the next section
3. The convergence of the calculations
Before starting the actual calculation for the lattice constant, it is important to explore the convergence of the calculations with respect to the number of k points and the energy cutoff. So the most accurate and efficient results can be obtained. For the following calculation, the lattice constant is fixed at c= 5.0511Å, a= 3.1946Å, while the size of k point grid and the energy cutoff are varied.
(1) The number of k-points in the irreducible part of the BZ
In this section 3.1, the energy cutoff is fixed at 400 eV, and the size of k point grid changes and Number of k-points in the irreducible part of the BZ (the number of k points) varies from 1 to 162, shown in Table 1. According to Table 1, the relations between Energy/atom, calculation time and the number of k points can be plotted, which are shown in Figure 2. As shown in Figure 2, the energy/atom becomes stable with the increase of k points. However, the calculation time also increases rapidly when the number of k points increases. By combining the results of these two graphs, 9*9*6 k points with 36 k points is chosen to obtain an accurate result with less calculation time.
A little discussion here, since my calculation involves unit cells with different volumes, the k points need to be chosen so the density of k points in reciprocal space is comparable for the different supercells. But as shown in section 4, since the volumes of unit cell change just a little during the calculation, fixing the number of k points will not make much difference. So in section 4, the number of k points is chosen as 36 for all calculations.
(2) The energy cutoff
In this section, the number of k point is fixed at 36, and the energy cutoff varies from 200 to 800, shown in Table 2. And Figure 3 shows the relation between energy/atom, calculation time and energy cutoff. Although the energy/atom keeps increasing with the increase of energy cut off, the differences between each energy/atom are becoming smaller and smaller. On the other hand, the calculation time increases rapidly with the increase of energy cut off. In order to save the time and keep the accuracy, the energy cutoff is chosen as 435.4 eV.
4 The calculation of the lattice constants
Using the results from section 3, the k points and the energy cutoff are chosen as 9*9*6 and 435.4 eV. Then lattice constants c and a are changed to obtain the minimum energy, shown in Table 3. As shown in Figure 4, the unit cell and the lattice constants can be predicted. As shown in Figure 4a, the quadratic relation between energy/atom and volume/atom only valid for a small range of lattice constants around the equilibrium value. A small range of volume/atom is chosen and the relation is re-plotted in Figure 4b. According to the fitting equation in Figure 4b, the fitting volume/atom is 22.64 Å^3 and the lattice constant are c=5.0673 Å, a=3.2072 Å. And for the experiment data, the lattice constants are a=3.1964 Å, c=5.0511 Å (1). So the calculation is accurate with an error of 0.3%.
5. Some discussion and future works
Although the calculation result is accurate with an error of 0.3%, there are still some ways to further increase the accuracy. First, the energy cut off can be increased so the energy/atom can be more precise. Second, we can calculate more points around the value c=5.05 Å and maybe a more accurate result can be obtained. Third, instead of using the quadratic fitting, we can use a higher order polynomial fitting.
For some future works: first, it is worthwhile changing the ratio c/a to determine the optimum value of c/a. Second, it should be interesting to analyze the behavior of energy at small or large lattice constants.
6. Reference:
(1) https://www.webelements.com/hafnium/crystal_structure.html