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