Author Archives: Sharad Maheshwari

Evaluating activation barrier of Ag hopping on Cu(111)

Leave a reply

Author: Sharad Maheshwari

Introduction

Surface diffusion can play an important role in creating a topological modification on the surface as in the case of nanoparticle growth. In the following work we aim to identify the kinetic barrier for the diffusion of a Silver (Ag) atom as it hops from a three-fold fcc site to the adjacent three-fold hcp site on (111) surface of Copper (Cu). We will use density functional theory-based methods to evaluate the activation barrier for the hopping process. We will also calculate the frequency of this hopping process.

Calculation Details

Electronic Calculations

Electronic structure calculations were performed using the Vienna Ab initio Simulation Package (VASP)[1, 2] a plane wave basis set pseudo-potential code. We used the projector augmented wave (PAW) [3] method for core-valence treatment. The exchange and correlation energies were calculated using the Perdew, Burke, and Ernzerhof (PBE) [4] functional described within the generalized gradient approximation (GGA) [5]. A plane-wave basis set cutoff energy of 450 eV was used. The ionic convergence limit was set to 0.01 eV Å-1 while the electronic convergence limit was set to 10-5 eV. The Fermi level was smeared with the Methfessel-Paxton [6] scheme using a smearing width (σ) of 0.2 eV.

Surface Model

A 3 x 3 surface slab model was used to construct periodic surfaces of Cu (111) using an experimental bulk lattice constant of 2.556 Å. The slab models were comprised of 4 layers of metal atoms.  The top two layers of the slab along with the adsorbate were allowed to relax until the force on atoms in these layers was within 0.01 eV Å-1 while the bottom two layers were kept fixed to imitate their bulk arrangement.  A vacuum region of 10 Å was inserted in the models to exclude periodic interaction between the slabs. Dipole corrections were also added in the direction normal to the surface.

The sampling of the Brillouin zone for all 3 x 3 surface cells was conducted with a k-point mesh of 5 x5 x 1 generated automatically using the Monkhorst Pack method[7]

Transition States

To locate the transition state, a series of 5 linearly interpolated images were formed between the reactant and the product state. The nudged elastic band method (NEB) [8] was used to locate the transition state. Once an approximate transition state was obtained, it was further refined using the climbing image nudged elastic band method (CI-NEB) [9]. Vibrational frequency calculations were then performed using IBRION = 7 in VASP to confirm the transition state found the first order saddle point. It was further ensured that the single imaginary frequency obtained was along the reaction pathway.

Results

Figure 1 illustrates the reactant, transition and product state for the hopping of Ag atom from the fcc to hcp site of Cu(111) surface. The activation barrier for the hopping process is 0.035 eV.

Figure 1. Top and Side view of the a) reactant b) transition c) product states for the Ag hopping from fcc site to the hcp site on Cu (111)

For the barrier evaluated using the CI-NEB method, the reactant and the product states were used to create a series of images. The energy profile obtained from the calculation for the images is shown below in Figure 2. (All energies are referenced to the reactant state and does not include zero- point vibrational energy corrections).

Figure 2. Reaction energy diagram for the hopping of Ag atom from fcc to hcp site on Cu(111)

One important thing that was observed in a trial calculation is that the NEB method missed the transition state when 4 linearly interpolated intermediate images were used as shown in Figure 3. Using small even number of images can result in missing the saddle point and thus if using the NEB method, a higher density of images are needed near the saddle point. This highlights the need for caution when evaluating the activation barrier using the NEB method and using linear interpolation to obtain the images.

Figure 3. Reaction energy diagram obtained from the calculations using NEB method with 4 images obtained through linear interpolation

On an absolute scale, the zero-point vibrational energy (ZPVE) corrections were found to be very small. However, as the activation barrier itself is very small, we corrected the activation barrier by including the ZPVE corrections to the initial and the transition state. The corrected activation barrier after including the ZPVE corrections is 0.027 eV.

To evaluate the frequency of this hopping, we use the following equation:

  

where υi and υi†  are the frequencies evaluated for reactant and transition state respectively, kB is the Boltzmann’s constant, ΔEact is the evaluated activation barrier and T is the temperature in Kelvin.  As the transition state is first-order saddle point, it is energy minimum in all direction except for the reaction coordinate, along which it is maximum. Figure 4 shows the hopping frequency of Ag atom from fcc to hcp site on Cu(111) surface at different temperatures.

Figure 4. Frequency of hopping of Ag atom from fcc to hcp site on Cu(111) surface as a function of temperature

Conclusion

We used DFT methods to estimate the kinetic barrier for the hopping of Ag atom from the fcc to hcp site on Cu(111) surface. The activation barrier for this process after including the ZPVE corrections is found to be 0.027 eV. The value found is in reasonable agreement with the reported values [10]. This work illustrates the DFT methods can be used to evaluate kinetic parameters of the process. It also highlights that while using methods to evaluate the activation barrier, one should be careful and must ensure that the transition state found is the correct transition state. Vibrational calculations must be performed and a single imaginary frequency along the reaction coordinate must be confirmed to ensure that the transition state found is indeed the transition state we were looking for.

 

References

[1] G. Kresse, J. Furthmuller, Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set, Comput. Mater. Sci., 6 (1996) 15-50.

[2] G. Kresse, J. Furthmuller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, Phys. Rev. B, 54 (1996) 11169-11186.

[3] G. Kresse, D. Joubert, From ultrasoft pseudopotentials to the projector augmented-wave method, Phys. Rev. B, 59 (1999) 1758-1775.

[4] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 77 (1996) 3865-3868.

[5] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Atoms, Molecules, Solids, And Surfaces – Applications of the Generalized Gradient Approximation for Exchange and Correlation, Phys. Rev. B, 46 (1992) 6671-6687.

[6] M. Methfessel, A.T. Paxton, High-precision sampling for Brillouin-zone integration in metals, Physical Review B, 40 (1989) 3616-3621.

[7] H.J. Monkhorst, J.D. Pack, Special Points for Brillouin-Zone Integrations, Phys. Rev. B, 13 (1976) 5188-5192.

[8] G. Henkelman, H. Jónsson, Improved tangent estimate in the nudged elastic band method for finding minimum energy paths and saddle points, The Journal of Chemical Physics, 113 (2000) 9978-9985.

[9]  G. Henkelman, B.P. Uberuaga, H. Jónsson, A climbing image nudged elastic band method for finding saddle points and minimum energy paths, The Journal of Chemical Physics, 113 (2000) 9901-9904.

[10] Kotri, E. El koraychy, M. Mazroui, Y. Boughaleb, Static investigation of adsorption and hetero-diffusion of copper, silver, and gold adatoms on the (111) surface, Surface and Interface Analysis, 49 (2017) 705-711.

Preferred adsorption site for atomic H on Pt (111) and coverage effects

Leave a reply

Author: Sharad Maheshwari

Introduction

In the following work, we aim to use plane wave Density Functional Theory (DFT) calculations to identify the preferred adsorption site for atomic H on (111) crystal surface of Platinum. We further seek to examine how the adsorption preferences change as the coverage of atomic H changes on the surface.

Adsorption sites on Pt (111)

Adsorption on (111) crystal surface of FCC metals can take place on 4 different sites with metal coordination as indicated in the parentheses, namely: atop (1), bridge (2), fcc (3) and hcp (3). The fcc and hcp site differ in the sub-surface stacking of the metal atoms. Figure 1 below illustrates all these binding sites with an H atom on these sites.

Figure 1. Different adsorption sites on Pt (111) crystal surface illustrated with H atom1. a) Atop b) Bridge c) FCC d) HCP. The white ball is H atom. Yellow balls indicate metal atoms in the first layer. Dark green balls are the metal atoms in the second layer. On FCC site, there is no metal atom in the second layer right below the site, whereas for HCP site, one can see a metal atom in the second layer exactly below the H atom.

 

Adsorption Energy

To estimate the strength of binding between adsorbate and the surface, we can evaluate the adsorption energy of the species on the surface.  For a reaction

A + * → A*

where, A is the adsorbate in the reference state, * indicates the surface and A* indicates the adsorbate bound on the surface, adsorption energy is can be given as

Eads = EA* – EA – E*                                       (1)

where Eads is the adsorption energy, EA* is the DFT energy of the adsorbate A bound to the surface, EA is the DFT energy of free A and E*  is the bare surface DFT energy It makes practical sense, however, to evaluate the adsorption energy with reference to a physically stable state of A. For diatomic molecule ( as in our case, H2), we can evaluate the reaction as

½ A2 + * → A*

and thus the corresponding adsorption energy is given as:

Eads = EA* – 0.5EA2 – E*                                       (2)

Coverage Effects

To evaluate the effect of surface coverage on the adsorption energy and preference to the adsorption site, we evaluate the adsorption energy per atom for different coverage. We define the coverage of adsorbate in terms of number of adsorbate atoms with respect to the number of metal atoms on the surface.  For e.g., 1/9 ML represents the 1 adsorbate atom per 9 surface metal atoms. The adsorption energy per metal atom is now evaluated as shown in eq. 3:

Eads = (EnA* – n/2EA2 – E*)/n                                     (3)

where n is the number of adsorbate atoms on the metal surface. We will use eq. 3 to report all the adsorption energies.

Calculation Details

Electronic structure calculations were performed using the Vienna Ab initio Simulation Package (VASP)[1,2] a plane wave basis set pseudo-potential code. We used the projector augmented wave (PAW)[3] method for core-valence treatment. The exchange and correlation energies were calculated using the Perdew, Burke, and Ernzerhof (PBE)[4] functional described within the generalized gradient approximation (GGA)[5]. A plane-wave basis set cutoff energy of 450 eV was used. The calculations were considered optimized when the force on every relaxed atom ( ionic convergence limit) was less than 0.02 eV Å-1. The electron densities were self consistently solved for the energy with the convergence limit (electronic convergence) set to 10-5 eV.

A 3 x 3 surface slab model was used to construct periodic surfaces of Pt (111) using an experimental bulk lattice constant of 2.775 Å. The slab models were comprised of 4 layers of metal atoms.  The top two layers of the slab were allowed to relax until convergence while the bottom two layers were kept fixed to imitate their bulk arrangement.   A vacuum region of 10 Å was inserted in the models to exclude periodic interaction between the slabs. Dipole corrections were also added in the direction normal to the surface.

The sampling of the Brillouin zone for all 3 x 3 surface cells was conducted with a k-point mesh of 5 x5 x 1 generated automatically using the Monkhorst Pack method [6]. K-point mesh of 1x1x1 was used for the isolated H2 molecule

Results

Figure 2 illustrates the low coverage (1/9 ML) adsorption energies for all different sites. At very low coverage limit, the fcc site is the most preferred site, followed by the atop site. The preference for hcp and bridge site is almost the same at low coverages. These results match well with the previously reported results.. [7]

Figure 2. Adsorption energy of atomic H on different binding sites on Pt (111).

To evaluate the coverage effect, we evaluate the adsorption energy per H atom for 1/9 ML, 2/9 ML, 3/9 ML and 1 ML on fcc, hcp and atop sites (Figure 3).

Figure 3. 2/9 ML, 3/9 ML and 1ML of H on a)atop b) FCC c) HCP site of Pt(111).

For higher coverages (2/9 ML, 3/9ML, 1ML), when the calculations were initiated from the bridge site, the H atoms moved to the fcc site during optimizations. Even after several tries, the optimized geometry for the adsorption on the bridge site at higher coverages could not be found.  Figure 4 shows the calculated adsorption energies per H atom for different coverages on atop, fcc and hcp sites.

Figure 4. Adsorption energy per H atom at 4 different coverages on atop, FCC and HCP sites on Pt (111).

As the coverage increases on the metal surface, we see that the adsorption energy per atom of adsorbate decreases. This is expected because as the coverage increases, it results in the repulsive interaction between the adsorbed species. We also see that as the coverage increases the preference for fcc and atop site becomes very similar. This trend is found to be true in the literature as well. However, the adsorption energies differ slightly as the energies reported in the literature are ZPVE (zero-point vibration energy) corrected whereas the numbers reported in this work do not correct for them [8,9].

References

[1] G. Kresse, J. Furthmuller, Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set, Comput. Mater. Sci., 6 (1996) 15-50.

[2] G. Kresse, J. Furthmuller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, Phys. Rev. B, 54 (1996) 11169-11186.

[3] G. Kresse, D. Joubert, From ultrasoft pseudopotentials to the projector augmented-wave method, Phys. Rev. B, 59 (1999) 1758-1775.

[4] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 77 (1996) 3865-3868.

[5] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Atoms, Molecules, Solids, And Surfaces – Applications of the Generalized Gradient Approximation for Exchange and Correlation, Phys. Rev. B, 46 (1992) 6671-6687.

[6] H.J. Monkhorst, J.D. Pack, Special Points for Brillouin-Zone Integrations, Phys. Rev. B, 13 (1976) 5188-5192.

[6] G. Papoian, J.K. Nørskov, R. Hoffmann, A Comparative Theoretical Study of the Hydrogen, Methyl, and Ethyl Chemisorption on the Pt(111) Surface, Journal of the American Chemical Society, 122 (2000) 4129-4144.

[7] I. Hamada, Y. Morikawa, Density-Functional Analysis of Hydrogen on Pt(111): Electric Field, Solvent, and Coverage Effects, The Journal of Physical Chemistry C, 112 (2008) 10889-10898.

[8] T.T.T. Hanh, Y. Takimoto, O. Sugino, First-principles thermodynamic description of hydrogen electroadsorption on the Pt(111) surface, Surface Science, 625 (2014) 104-111.

 

DFT Predictions for Lattice Parameters of Hafnium (Hf)

Leave a reply

by Sharad Maheshwari

Introduction

The aim of the following work is to use Density Functional Theory to predict the lattice parameters of the bulk crystal lattice of Hafnium (Hf). Experimentally, Hf crystallizes in a hexagonal close-packed structure with the ratio of lattice parameters, c/a = 1.58 [1]. Herein, we would use the DFT calculations to evaluate the c/a ratio and compare it with the above experimental value.

Figure 1. Top, side and perspective view of the primitive unit cell of Hafnium (Hf) in hexagonal closed packing with two lattice parameter, c and a.

 

Methodology

We use the plane wave DFT calculations to evaluate the binding energy of a primitive lattice of Hf as shown in the equation below. By minimizing the binding energy as a function of lattice parameters, we will evaluate the ratio of lattice parameters (c/a).

Binding Energy  = DFT Energy Calculated – 2*Atomic Energy of Hf

All the calculations were performed using CASTEP Simulation Package [2] in Material Studio. The exchange and correlation energies were calculated using the Perdew, Burke, and Ernzerhof (PBE) [3] functional described within the generalized gradient approximation (GGA) [4]. Electronic convergence tolerance of 2E-06 was used for all the calculations. The core  electrons were treated using “on the fly” generated (OTFG) ultra-soft pseudo-potential with core radius 2.096 Bohrs (1.109 Angstrom) generated  with panel of 26 valence electrons (4f14 5s2 5p6 5d2 6s2).

Cut-off Energy and K-Point Optimization

We use a plane wave basis set to perform our DFT calculations for the periodic lattice of Hf. In order to use the plane wave basis sets, it is essential to ensure the convergence of system energy with respect to the cut-off energy and k-point mesh. Therefore, we first perform cut-off energy and k-point optimization.  For this, we have used the hcp lattice with a = 3.1946 Å and c = 5.0511 Å which are the default dimensions for the Hf lattice in Material Studio and is very close to the experimental values[5, 6].

Cut-Off energy

In Fig. 2, we have plotted the absolute change in the consecutive total energy (|ΔE|) as the cut-off energy is varied. As a convergence criterion, we assume that the total energy is converged if |ΔE| is less than 0.003 eV. From Fig. 2, we can say that the cut-off energy of 500 eV sufficiently converged the total energy of the system.

Figure 2. Convergence of absolute change in consecutive total DFT energy with respect to cut-off energy. The blue horizontal line indicates the convergence criterion used (0.003 eV) and the plot is further zoomed in to illustrate the convergence

K-Point

In Fig. 3, we have plotted the absolute change in the consecutive total energy (|ΔE|) as the number of irreducible k-points, corresponding to different k-point grid used as shown in Table 1., is changed. As a convergence criterion, we assume that the total energy is converged if |ΔE| is less than 0.003 eV. From Fig. 3, we can say that k-point mesh generated automatically using the Monkhorst Pack method [7] of 10X10X6, corresponding to 42 irreducible k-points samples the Brillouin zone sufficiently well and converges the total energy of the system.

Figure 3. Convergence of total DFT energy with respect to irreducible k-points corresponding to different k-point grid. The orange horizontal line indicates the convergence criterion used (0.003 eV).

Table 1. K-point grid (xXyXz) and corresponding irreducible k-points
K-point gridIrreducible k-points
8X8X530
9X9X636
10X10X642
11X11X764
12X12X876
13X13X884

 

Hf Lattice Optimization

To optimize the hcp lattice of Hf, we minimize the binding energy of the system with respect to the two lattice constant, a and c.  To carry out this two-parameter optimization, we keep the ratio of c/a constant as the value of a is varied and then the same procedure is repeated for different values of c/a.

In order to ascertain the minima for any specific ratio of c/a, the binding energy is evaluated at 3 distinct values of a. A parabola is then fit through these points and the value of a corresponding to the minima of the parabola is evaluated using the equation of the parabola. The binding energy is then calculated at this evaluated value of a to check if this indeed is the minima. The binding energy is also evaluated for the points in the vicinity of this minimum to ensure that we find the least energy for the given c/a ratio.

Results

The above procedure for lattice optimization was repeated for different values of c/a. Fig. 4 plots the calculated binding energy with respect to a for different c/a ratio. Table 1 lists the minimum energy obtained for different c/a ratio and corresponding value of a.

Figure 4. DFT Energy of the crystal lattice of Hf as a function of one of the lattice parameter, ‘a’ for different values of the ratio c/a. An “Overall” curve passes through minimum for different “c/a” ratios and thus minimizing w.r.t. different c/a

Table 2. Minimum DFT binding energy calculated and the corresponding value for the lattice parameter a for different values of the ratio c/a
c/aBinding Energy (eV)
 a (Å)
1.581-17.487
3.204
1.481-17.649
3.282
1.781-17.693
3.082
1.381-17.648
3.355
1.681-17.541
3.141

A parabola is then plotted to pass through the minima obtained for different c/a ratio. Minimum of this parabola gives us the local minimum of binding energy with respect to c/a and thus, the optimum value for c/a.

The above optimization for binding energy minimization occurs at c/a = 1.581 and a = 3.204 Å (c =5.066 Å ). This value evaluated using the DFT calculation agree well with the experimental observation [6, 7].

Conclusion

Density Function Theory based calculations were used to ascertain the lattice parameters of the metal Hafnium. The binding energy of the primitive hcp lattice of Hf was minimized by altering the c/a ratio and the value of the lattice parameter a. The optimized lattice parameters obtained are reasonably coherent with the experimental results. This project thus helps illustrate that DFT can be used to predict lattice parameters reasonably well.

References

[1] D. Sholl, J.A. Steckel, Density Functional Theory: A Practical Introduction, Wiley2009.

[2] J. Clark Stewart, D. Segall Matthew, J. Pickard Chris, J. Hasnip Phil, I.J. Probert Matt, K. Refson, C. Payne Mike, First principles methods using CASTEP,  Zeitschrift für Kristallographie – Crystalline Materials, 2005, pp. 567.

[3] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 77 (1996) 3865-3868.

[4] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Atoms, Molecules, Solids, And Surfaces – Applications of theTHE Generalized Gradient Approximation for Exchange and Correlation, Phys. Rev. B, 46 (1992) 6671-6687.

[5] http://periodictable.com/Elements/072/data.html.

[6] https://www.webelements.com/hafnium/crystal_structure.html.

[7] H.J. Monkhorst, J.D. Pack, Special Points for Brillouin-Zone Integrations, Phys. Rev. B, 13 (1976) 5188-5192.