Author Archives: Jeff Rable

Energies and Lattice Constants of Ferromagnetic and Antiferromagnetic States in Co, MnF2, and MnBi

In this project, the energies and lattice constants of the ferromagnetic and antiferromagnetic states of Co, MnF2, and MnBi were found and compared. To perform these calculations, a suitable energy cutoff and k-space sampling were selected to obtain a 0.01eV tolerance, and research was done in order to determine a proper starting c/a value of the lattice constant. From there, the the materials were placed into a ferromagnetic or antiferromagnetic state by aligning or misaligning their atomic spin polarizations and their energies were found at varying values of lattice constant a using CASTEP plane wave calculations [1].

These calculations were performed in Materials Studio using the GGA Perdew Burke Ernzerhof (PBE) functional and the pseudopotential was solved using the Koelling-Harmon atomic solver. In Materials Studio, spin polarizations are defined by the the total number of up electron spins minus down electron spins, so a spin polarization of +3, for example, implies that the atom has 3 more spin up than spin down electrons.

Cobalt

Fig. 1. Picture of cobalt cell used

For these calculations, cobalt was imported using the default Materials Studio structure, which gives a lattice constant of 2.507Å, a c/a value of 1.623, and forms the P63/MMC symmetry. From there, the spins in the 2 atom unit cell were polarized into a ferromagnetic +3, +3 state and energy calculations with varying energy cutoffs and k-samplings were performed until the results converged to a 0.01eV tolerance. The psuedopotential was calculating using the 1s2 2s2 2p6 3s2 3p6 electrons for the inner shells and 3d7 4s2 as the valence shells. Ultimately, a 420eV cutoff with a total k space sampling of 108 points were used.

Fig. 2. Energy cutoff calculations for Co

Fig. 3. k-space Sampling calculations for Co

From here, the lattice constant a was varied in 0.1Å and 0.005Å increments to find a suitable ferromagnetic lattice constant. These same calculations were repeated for the antiferromagnetic state, which were started with opposite +3 and -3 spin polarizations.

Fig. 4. Lattice constant calculations for ferromagnetic and antiferromagnetic states of Co

In the antiferromagnetic calculations, as the atoms approached one another, the calculations ceased to reach a local minimum in the antiferromagnetic state and began converging to non-spin polarized states. Additionally, some calculations failed to converge within 200 steps. Regardless, Co favors a ferromagnetic state, with a lattice constant of approximately 2.49Å, over the antiferromagnetic state, which appears to have a lattice constant between 2.465Å and 2.45Å (energy decreases past 2.465Å, but calculations converge to a non-spin-polarized state). This is in good agreement with the experimental value of 2.507Å in the ferromagnetic state.

Table 1. Results of the cobalt ferromagnetic and antiferromagnetic calculations

MnF2

Fig. 5. Picture of cell used

For these calculations, MnF2 was constructed using the analogous rutile TiO2 structure, which has the P42/MNM symmetry. The starting value of lattice constants c/a and a, 0.68 and 4.875Å respectively, used in these calculations were found by Dormann et. al [2]. From there, the spins in the atomic unit cell were polarized into a ferromagnetic state and energy calculations with varying energy cutoffs and k-samplings were performed until the results converged to a 0.01eV tolerance. The psuedopotential was calculating using the 1s2 electrons for the inner shells of fluorine and 2s2 2p5 as the valence shells. For manganese, the 1s2 2s2 2p6 electrons were used for the inner shells and 3s2 3p6 3d5 4s2 for the valence shells.

Fig. 6. Energy cutoff calculations for MnF2

Fig. 7. k-space Sampling calculations for MnF2

Calculations were made by shifting the lattice constant a in 0.1Å increments to get a rough idea of the lattice constant and 0.01Å measurements to take a finer sweep around the minimum energies of the two states. Ultimately, MnF2 converged to a lattice constant of approximately 4.930Å in the favored antiferromagnetic state and 4.945Å in the less favorable ferromagnetic state. This is relatively close to the experimental value of 4.875Å in the antiferromagnetic state.

Fig. 8. Graph of energies vs. lattice constants for the MnF2 lattice with ferromagnetic and antiferromagnetic spin polarization.

Table 2. Table of energies vs. lattice constants for MnF2.

MnBi

Fig. 9. Picture of cell used

For these calculations, MnBi was constructed using the analogous NiAs structure, which has the P63/MMC symmetry. The starting value of lattice constants c/a and a, 1.426 and 4.270Å respectively, used in these calculations were found by Koyama et. al [3]. From there, the spins in the atomic unit cell were polarized into a ferromagnetic state and energy calculations with varying energy cutoffs and k-samplings were performed until the results converged to a 0.01eV tolerance. The psuedopotential was calculating using the 1s2 2s2 2p6 3s2 3p6 3d10 4s2 4p6 4d10 4f14 5s2 5p6 electrons for the inner shells of bismuth and 5d10 6s2 6p3 as the valence shells. For manganese, the 1s2 2s2 2p6 electrons were used for the inner shells and 3s2 3p6 3d5 4s2 for the valence shells.

Fig. 10. Energy cutoff calculations for MnBi

Fig. 11. k-space Sampling calculations for MnBi

Next, the lattice constant was varied in 0.05Å increments to find an approximate energy minimum for the ferromagnetic and antiferromagnetic states. The lattice constant was then incremented in finer 0.01Å increments around the minimum value to obtain more accurate results. Ultimately, the lattice constant for both states came out to 4.23Å, with the ferromagnetic state being favored over the antiferromagnetic state. This is extremely close to the experimental value of 4.27Å in the ferromagnetic state.

Fig. 12. Graph of energy vs. lattice constant for the two ordered magnetic states of MnBi

Table 3. Table of energy vs. lattice constant for the two ordered states of MnBi

Discussion

These results all agree with the experimentally accepted lattice constant values for the three materials to a good degree of accuracy. In addition, they correctly determined that Co and MnBi are ferromagnetic, while MnF2 is antiferromagnetic. However, the effect of spin polarization in these materials is quite small – on the order of <1eV – suggesting that the energy change is relatively small even for materials considered to be ferromagnets and antiferromagnets.This provides some justification for ignoring spin polarization in our previous calculations, such as those involving copper, which is diamagnetic.

Another interesting note is that the antiferromagnetic calculations required breaking the initial lattice symmetry. This resulted in the antiferromagnetic calculations taking approximately 3-4 times longer than the ferromagnetic calculations; for example, the a = 4.18Å MnBi calculations  took 1701s for the P1 antiferromagnetic state versus 458s for the P63/MMC ferromagnetic state.

Though all the antiferromagnetic calculations in this post were performed with only P1 symmetry, the author realized afterwards that Materials Studio enables one to find a higher symmetry state factoring in the spins. By going into build->Symmetry->Find Symmetry->Options, one can force Materials Studio to consider the formal spins of each atom when searching for symmetries. Using the symmetry finder with these settings, one can reduce the  a = 4.18Å MnBi antiferromagnetic calculation time to 295s, a nearly 6-fold decrease in calculation time.

References

Zero Point Energy of Hydrogen in Cu (111) Layers Near the Surface

In this project,  the zero point energy (ZPE) of a hydrogen atom in the 2nd and 3rd layers below a Cu (111) surface was found by calculating the Hessian matrix of the atom. To do this, the hydrogen atom was shifted by 0.1Å along the x, y, and z directions of the lattice and the energy was calculated at each location using CASTEP [1]. Using these energy values, a second order finite difference approximation was used to find the Hessian, which was then diagonalized to find the vibrational frequencies of the atom. These frequencies could then be used to calculate the ZPE of the defect by approximating its oscillations as that of a quantum harmonic oscilator [2].

These calculations were performed in Materials Studio using the GGA Perdew Burke Ernzerhof (PBE) functional. The calculations used an energy cutoff of 520eV for the b-layer substituted hydrogen and a k point sampling of 15x15x2, for a total of 225 points. For the c-layer, a 540eV cutoff and a 15x15x2 k point sampling were used. The pseudopotential was solved using the Koelling-Harmon atomic solver. For Cu, the solver used 1s2 2s2 2p6 3s2 3p6 as the interior electron shells and 3d10 4s1 as the valence shells. For H, the solver used 1s1 as the valence shell.

The cutoff value and k point sampling were determined by finding where the energy converged to a 0.01eV tolerance.

Construction of the Surface

To build the desired surface, we started with a pure metal Cu lattice. This lattice was then cleaved along the (111) axis with a thickness of 4 atoms, so as to contain both B and C layers, as well as two A layers. The bottom A layer was fixed and a vacuum slab with 10Å vacuum was constructed. Then, a hydrogen atom replaced the B layer hydrogen and a symmetric slab was constructed to balance out the dipole moment of the hydrogen in the plane wave calculations. The geometry of this surface was then optimized with a 450eV cutoff and a k-sampling of 11x11x1. This structure was then used to find the proper cutoff energy and k-space sampling to get a 0.01eV convergence. The cutoff for the B layer ended up being 520eV and the k-sampling was 15x15x2. 540eV was used for the C layer cutoff. Next, the structure was geometry optimized again using the proper cutoff and k-samplings.

Fig. 1. The appearance of the b-layer defective structure before geometry optimization was performed.

Fig. 2. The appearance of the c-layer defective structure after geometry optimization was performed.

Energy at Different Hydrogen Locations

Next, the hydrogen atoms were symmetrically shifted by 0.1Å along the x, y, and z directions as needed to construct the Hessian matrix of the atom. To get the energy of only one substitutional hydrogen, these energies were halved when calculating the Hessian matrix.

Table 1. The values of the energy of the structure with a hydrogen atom located in the B layer. Here, the atoms were shifted by 0.1Å and the new energy values of the structure were found. In order to keep the structures inverted, the atoms were moved in the opposite directions, with the chart labeling of the movement of one of the atoms.

Table 2. The values of the energy of the structure with a hydrogen atom located in the C layer.

Hessian Construction and Results

A finite difference approximation and the values from Table 1 and 2 were used to construct the Hessian matrices. The diagonal elements were found using the formula:

\begin{equation} E_{xx} = \frac{E(\Delta x, 0) – 2E(0,0) + E(-\Delta x, 0)}{ \Delta x^2} \end{equation}

and the off diagonal elements were found using the formula:

\begin{equation} E_{xy} = \frac{E(\Delta x, \Delta y) – E(-\Delta x,0) – E(0, -\Delta y) + E(-\Delta x, -\Delta y)}{4 \Delta x \Delta y} \end{equation}

These elements were then placed into a matrix in Wolfram Mathematica, which was used to find the eigenvalues of the matrix, which gives

\begin{equation} 2 \pi \nu _i = \omega _i = \sqrt{ \lambda _i} \end{equation}

Using these values and the QHO approximation, we can find the vibrational energies along each direction with:

\begin{equation} E_{ZPE} = \frac{ \hbar \omega _i}{2} \end{equation}

Yielding energies of 0.047eV, 0.047eV, and 0.033eV along the x, y, and z directions respectively for the B layer hydrogen. For the C layer, ZPE energies obtained were 0.089eV, 0.089eV, and 0.065eV respectively. This gives a total ZPE of 0.128eV for the B layer HCP hydrogen defect and 0.243eV for the C layer FCC defect. Because both structures have identical classical energies within the tolerance used for this study, these results suggest that hydrogen prefers the B layer HCP site over the C layer FCC site.

References

Determination of the Optimal Geometry of Water Molecules

In this project, the optimal geometry of individual water molecules was found using CASTEP geometry optimization calculations [1]. These calculations utilized the GGA Perdew Burke Ernzerhof (PBE) functional and the Koellig-Harmon atomic solver, which was used to find  psuedopotentials for the hydrogen 1s1 electrons and on the oxygen 1s2 2s2 2p4 electrons. In these calculations, the H atoms had the 1s1 electrons used as valence electrons and the 2s2 2p4 electrons of the O atoms were used as the valence electrons. Additionally, the geometry optimization calculations were performed with a 0.01 eV/Å force tolerance.

Fig. 1. Water Molecule

Determination of Lattice Size

Because CASTEP relies on plane waves, one must construct a lattice to find a molecule’s geometry. Additionally, this lattice must be large enough that the neighboring molecules do not interact in a way that would significantly alter their geometry.  Thus, in performing these calculations, one must ensure that their chosen lattice constant is sufficiently large that the molecule geometry converges past the desired tolerance. In this case, the atomic distances converged to 0.01Å with a 9Å lattice constant, which was used for the primary calculation.

(a)

(b)

Fig. 2. Graph of atomic distances vs. Lattice Constant for (a) OH (b) HH.

Table 1. Effects of lattice constant on bond lengths and angles.

Determination of Energy Cutoff

In order to find the desired energy cutoff for our calculations, the energy cutoff was varied until the distances between atoms converged to 0.01Å. A 435eV cutoff energy was found to be sufficient, and a cutoff of 630eV was used for the main calculation.

(a)

(b)

Fig. 3. Graph of atomic distances vs. cutoff for (a) OH (b) HH.

Table 2. Effects of energy cutoff on atomic distances

Determination of k-space Sampling

To ensure that our k-space sampling was large enough, we checked that the atomic distances converged to 0.01Å at the desired sampling. A sampling of 2x2x2 points, which, after symmetry considerations, resulted in a 4 point total sampling, was found to be sufficient.

(a)

(b)

Fig. 4. Graph of atomic distances vs. cutoff for (a) OH (b) HH.

Table 3. Comparison of total k-space sampling and atomic distances. A single k point was found to be sufficient, and a 4 point (2x2x2) sampling was used in the primary calculation.

Symmetry Considerations

When performing geometry optimization calculations, one must carefully consider the symmetry of the system – in the event that the molecule contains certain symmetries, the forces between atoms may cancel and result in a calculation-breaking local minimum. In this case, when the calculation is performed with the hydrogen atoms perfectly symmetric about an axis (see fig. 5), the calculation fails; the geometry remains symmetric, with OH bond lengths of 0.943Å. However, by starting from slightly modified, symmetry breaking starting positions, we can find values that properly converge.

Fig. 5. A crystal of symmetric H2O molecules with a 9Å lattice spacing. The symmetry about the central oxygen atom causes the calculation to converge to the wrong value.

Conclusion

Utimately, these calculations found that H2O has an OH bond length of 0.97Å with a bond angle of 104.2 degrees. These closely match previously found experimental data  of 0.9572Å and 104.52 degrees [2].

References

Determining the Structure and Lattice Constant of Platinum

In this project,  the lattice constant and structure of crystalline platinum were found by performing CASTEP energy calculations and finding where the energy was minimized [1]. These calculations were performed in Materials Studio using the GGA Perdew Burke Ernzerhof (PBA) functional. The calculations used an energy cutoff of 321.1eV for the rougher calculations and 420eV for the finer calculations. The pseudopotential was solved using the Koelling-Harmon atomic solver with the interior shells 1s2 2s2 2p6 3s2 3p6 3d10 4s2 4p6 4d10, giving an energy of -13041.2296 eV. The outer shells used in the calculations consisted of the 4f14 5s2 5p6 5d9 6s1 electrons.

1. Simple Cubic Lattice (P23 Point Group)

Fig. 1. A primitive cell of the simple cubic lattice

For the simple cubic lattice rough calculations, 6x6x6 k points were sampled with 0.1eV Gaussian smearing and no shift. After applying symmetry, this sampling reduced to 11 total k points.

Table 1. The roughly calculated energy of the simple cubic lattice for varying values of the lattice constant.

Fig. 2. The roughly calculated energy of the simple cubic lattice configuration for different lattice constants. The minimum occurs at 2.6 Angstroms, which will help determine a starting value for both the finer, more time consuming simple cubic calculations and the rough calculations for other lattice structures. The lines connecting the data points act as a guide to the eye and are not data.

Fig. 3. K space sampling vs. energy (eV) for the P23 simple cubic lattice. The energy converges to 0.01eV beyond 176 points.

Upon finding the approximate lattice constant with an 11 point k space sampling and a 321.1eV energy cutoff, we can perform the same calculations around that point with a higher k space sampling and cutoff energy to find a more accurate result. For these calculations, 176 k points were sampled (16x16x16) and a cutoff of 420eV was used. These results should converge to 0.01eV. Ultimately, these calculations resulted in a lattice constant of approximately 2.62 Angstroms.

Table 2. The finer lattice constant vs. energy calculations.

Fig 4. Finer lattice constant vs. energy graph. The lines connecting the data points act as a guide to the eye and are not data.

 

2. Hexagonal Close-Packed Lattice (D3H-3 Point Group)

Fig. 5. A primitive cell of the hexagonal close-packed lattice.

For the hexagonal close-packed lattice, 24x24x18 k points were sampled with 0.1eV Gaussian smearing and a shift of (0.021inverse Angstroms ,0.021 inverse Angstroms,0). The calculation was ultimately performed using 882 k points. For this lattice structure, only calculations with a 420eV cutoff were performed resulting in approximately a 2.60 Angstrom lattice constant. For these calculations, the ratio c/a was set at 1.633 because it gives ideal hard sphere close packing.

Table 3. Table of the HCP lattice constant (angstroms) and corresponding energy. The minimum lies at approximately 2.60 and 2.61 Angstroms.

Fig. 6. Graph of the Pt HCP lattice constant (Angstroms) vs. energy (eV). The lines connecting the data points act as a guide to the eye and are not data.

Fig. 7. K space sampling vs. energy (eV) for the HCP lattice. The energy converges to 0.01eV beyond 882 points.

3. Face Centered Cubic Lattice

Fig. 8. A conventional cell of the FCC cubic lattice.

3a. F23 Point Group

For the first F23 FCC lattice calculations, 8x8x8 k points were sampled with 0.1eV Gaussian smearing and no shift. This resulted in the rough calculation being performed on 88 total k points and the finer calculations being performed on 176 k points. The finer calculation were performed on 16x16x16 k points, resulting in a total sample of 688 points, with an energy cutoff of 420eV. Both the rough and fine calculations resulted in a lattice constant of 2.80 Angstroms.

Table 4. Shows the roughly calculated energy of the F23 FCC lattice compared to the lattice constant.

Fig. 9. Shows a rough calculation of the energy of the F23 FCC lattice vs. the lattice constant. The minimum occurs around 2.8 Angstroms, which will be used as a starting point for the finer calculations. The lines connecting the data points act as a guide to the eye and are not data.

Fig. 10. K space sampling vs. energy (eV) for the F23 FCC lattice. The energy converges to 0.01eV beyond 76 points.

Table 5. Shows the lattice constant (Angstroms) vs. energy for the F23 FCC lattice with a finer k space sampling and energy cutoff. The minimum lies at approximately 2.8 Angstroms.

Fig. 11. Graph of the lattice constant (Angstroms) vs. the finely calculated energy (eV) of the configuration. The lines connecting the data points act as a guide to the eye and are not data.

3b. FM-3M Point Group

For the second FCC lattice calculation, the FM-3M point group was used, as platinum has been shown to form an FCC lattice under this point group in nature [2]. For the rough calculation, 8x8x8 k points were sampled with a 0.1eV Gaussian smearing and no shift, for a total of 20 k points. For the finer calculations, 16x16x16 k points were sampled, for a total of 120 points.

Table 6. Shows the energy of the FM-3M FCC lattice compared to the lattice constant.

Fig. 12. Shows the lattice constant vs. the energy of the FM-3M FCC Lattice for a rough calculation. The minimum occurs at 3.95 Angstroms, which will be used as a starting point for finer calculations.The lines connecting the data points act as a guide to the eye and are not data.

Fig. 13. K space sampling vs. energy (eV) for the F23 FCC lattice. The energy converges to 0.01eV beyond 35 points.

Table 7. Results for the fine calculations of the energy of the FM-3M FCC lattice at various lattice constants.

Fig. 14. Graph of the fine calculations of the lattice constant (Angstroms) vs. energy (eV) for the FM-3M FCC Lattice. The lines connecting the data points act as a guide to the eye and are not data.

For the FM-3M FCC lattice, the fine calculations were performed with a 420eV energy cutoff and a k space sampling of 16x16x16, resulting in a total of 120 k points being used. These results gave a minimum energy of -52203.84eV at 3.96 Angstroms.

Final Results

(a)

(b)

Fig. 15. (a) Comparing the rough results from the F23 FCC lattice, the P23 simple cubic lattice, and the D3H-3 HCP lattice with c/a=1.633, the F23 lattice has the lowest energy minimum. (b) Comparing all four lattice structures, the FM-3M FCC lattice energy is approximately 4 times lower and reaches a minimum at about 1.16 Angstroms higher lattice constant. The lines connecting the data points act as a guide to the eye and are not data.

For the first three lattice configurations, the F23 FCC lattice had the lowest energy minimum at lattice constant 2.8 Angstroms. However, the FM-3M FCC lattice has a lower energy by approximately a factor of 4, with a minimum at 3.96 Angstroms. This lattice configuration closely matches with experimental data, which shows that platinum forms a lattice in the FM-3M point group with lattice constant 3.92 Angstroms [2].

Appendix: Energy Cutoff Convergence

For the finer energy vs. lattice constant calculations, a 420eV cutoff energy was chosen based on the energy’s convergence to 0.01eV beyond this point.

Fig. 16. Graph of the energy cutoff (eV) vs. energy (eV) for platinum in an F23 FCC lattice. The energy converges beyond a cutoff of approximately 390eV and 420eV was ultimately used for the finer energy calculations.

Bibliography

2. Povarennych, A. & Povarennyck, A. Crystal chemical classification of minerals. 192 (Plenum Press, 1972).