Author Archives: Hoang Tran

Study of O-H bond stretching frequency at different adsorption environments of OH adsorbed on Pt(111)

Author: Hoang (Bolton) Tran

Problem Statement

Previously (Post 2), the energetic and tilt angle of hydroxyl (OH) adsorption on Pt (111) was studied. In this study (Post 3), the OH bond stretching vibration \((\nu_{OH})\) is of primary interest. Firstly, the zero-point energy (ZPE) will be considered for adsorption energy. Secondly, the adsorption energy \((E_{ad})\) and \(\nu_{OH}\) will be computed and analyzed at different OH adsorption environments:

  • Adsorption sites: Atop vs. bridge (see Fig. 1 in Post 2)
  • Surface coverage: 1.0 ML vs. 0.25 ML (Monolayer – number of adatom per surface atom)
  • Water solvation: introduce water molecules above the adsorbed molecule.

Methodology

Software:

DFT was implemented in both CASTEP and DMol3, a plane wave basis set software embedded inside BIOVIA Material Studio (MS). [1]

Calculation:

Geometric optimization (CASTEP): Please see Post 2

Vibrational analysis (DMol3):

  • Exchange-Correlation Functional: Generalized Gradient approximations – Perdew-Burke-Ernzerhof (GGA-PBE) [2]
  • Self-Consistent Field tolerance: 1E-05 eV
  • Monkhost-Pack [3] grid size: 10x10x1 (50 irreducible kpoints) for 1×1 supercell & 5x5x1 (13 irreducible kpoints) for 2×2 supercell.
  • Self-consistent dipole correction in z-direction (normal to surface)

The procedure of optimizing, cleaving and adsorbing OH onto Pt(111) surface are the same as outlined in Post 2. This procedure was carried out again for each adsorption environments outlined below. The optimized structure of each adsorption environment was then used as input for the vibrational frequency calculation.

ZPE correction:

The adsorption energy is defined as:

\begin{equation} E_{ad}=E_{OH/slab}-E_{OH}-E_{slab} \end{equation}

To calculate the ZPE correction for \(E_{ad}\), ZPE was calculated for both the free OH radical and the adsorbed OH. Since the free OH only has O-H stretching vibration, the ZPE of adsorbed OH also considered only the O-H stretching \(\nu_{OH}\).

Adsorption environments:

1.  Adsorption sites:

The two adsorption sites of interest are the high symmetry atop and bridge sites. Values of \(E_{ad}\) at 1.0 ML coverage were already calculated from previous study. See Post 2 for more details.

2. Surface coverage:

Adsorption energy at the atop site for two difference OH coverage values are considered: 1.0 ML and 0.25 ML.  In Post 2, \(E_{ad}\) was already calculated for 1.0 ML. To achieve 0.25 ML, 2×2 supercell with one absorbed OH was built (Fig. 1).

1.0 ML
0.25 ML

Figure 1. Visualization of OH surface coverage of 1.0 ML and 0.25 ML (click on tab to view)

3. Water solvation:

Solvation effect, or the interaction of water molecules near surface with the adsorbed OH was explored. Atop sites at 0.25 ML was used with 1 or 2 water molecules floating just above the adsorbed OH (Fig. 2).

1 water
2 water

Figure 2. Optimized structure of water floating above adsorbed OH (click on tab to view)

The adsorption energy for the system with extra water molecules were calculated differently:

\begin{equation} E_{ad}=E_{OH/slab/water}-E_{OH}-E_{slab/water} \end{equation}

In this case, \(E_{slab/water}\) was calculated by optimizing a bare slab with floating water(s).

Results and discussion

1. ZPE corrections:

The adsorption energies, calculated across different adsorption environments, with and without ZPE are shown in table 1 below.

Table 1. Comparing \(E_{ad}\) with and without ZPE:

Adsorption environments\(E_{ad}\) without ZPE (eV)\(E_{ad}\) with ZPE (eV)
1.0 MLBridge-2.14-2.15
Atop-3.02-3.03
0.25 MLAtop 0 water-2.90-2.90
Atop 1 water-2.99-2.99
Atop 2 waters-3.21-3.22

The difference in ZPE are all 0.01 eV or less, which is a lot smaller than the difference between \(E_{ad}\) at different environment. Therefore, it is concluded that ZPE does not affect the adsorption energy of OH on Pt (111) significantly.

This result also has significant implementation later.

2. Adsorption environments:

Calculated values of \(\nu_{OH}\) and \(E_{ad}\) for each environments are shown in Figure 3 and 4 below:

Figure 3. Frequency of OH stretching at different adsorption environments and of free OH.

Figure 4. Adsorption energy of OH at different environments.

The hypothesis of this study is: The adsorption energy is stronger (more negative) when the adsorption environment (site, coverage, solvation) allows the adsorbed OH to form stronger hydrogen bonds.

The strength of hydrogen bond is characterized by \(\nu_{OH}\), weaker OH stretching (lower \(\nu_{OH}\)) means stronger hydrogen bonds. The presented results put this hypothesis to the test.

Adsorption sites:

Looking first at the 1.0 ML (blue section) in Fig. 3, it is observed that the OH vibrations at atop and bridge sites are both smaller than that of free OH radical, indicating some degree of hydrogen bonds formed at the sites. The hydrogen bond strength is greater at the atop site. Looking now at 1.0 ML section in Fig. 4, atop adsorption is stronger than that of the bridge site (results of Post 2).

Surface coverage:

Results at “1.0 ML Atop” and “0.25 ML Atop 0w” represent the direct comparison for surface coverage. \(\nu_{OH}\) at 0.25 ML is about the same as free OH, indicating the absence hydrogen bonds. Fig. 4 then shows that OH adsorbs weaker at 0.25 ML.

Water solvation:

Focusing now at the 0.25 ML sections (orange sections). Surprisingly, adding one water (Atop 1w) does not lower \(\nu_{OH}\), indicating that there are still no hydrogen bonding. This, however, makes sense by observing the position of this water molecule in Fig. 2 above. The placement of this water molecule did not seem to produce hydrogen bonding, which would be water’s O atom pointing toward OH’s H atom.

Nevertheless, OH adsorption is still slightly stronger with one water. This can be attributed to the solvation effect on the slab itself, because this effect is the only difference left between the two equations for \(E_{ad}\) (\(E_{slab}\) and \(E_{slab/water}\)).

Considering now the case with two water, \(\nu_{OH}\) is now lower, indicating formation of hydrogen bonding, which is also confirmed by looking at position of the second added water in Fig. 2. Finally, the adsorption with 2 waters is much stronger than with both 0 water and 1 water by looking at Fig. 4.

Conclusion

Despite the above analysis seems to indicate that stronger hydrogen bonds lead to stronger adsorption, the aforementioned hypothesis is false. The most obvious evidence is from the ZPE correction. It was already shown that ZPE does not contribute significantly to the adsorption energy. The ZPE correction only comes from the difference in frequency of OH stretching \(\nu_{OH}\). Therefore, \(\nu_{OH}\), despite varying between different environments, does not contribute significantly to \(E_{ad}\).

In other words, the most that one can conclude from these results is that there is a correlation between hydrogen bond strength \(\nu_{OH}\) and adsorption energy \(E_{ad}\). The reason behind this correlation can be thought of as hydrogen bonds formation lowering the energy of adsorbed OH \(E_{OH/slab}\) or \(E_{OH/slab/water}\), which subsequently lowers adsorption energy \(E_{ad}\). However, it would be wrong to conclude that \(\nu_{OH}\) is solely responsible for the difference in \(E_{ad}\), because the magnitude does not check out.

Finally, another interesting takeaway is that water solvents increase the adsorption tendency, not merely by interacting directly with the adsorbed species (through hydrogen bonding), but potentially through other electronic interaction with the surface as well. This is why understanding of the solvation effects on solid-liquid interfacial processes (e.g. electrocatalysis) is of great interest.

References

[1] “First principles methods using CASTEP”, Zeitschrift fuer Kristallographie 220(5-6) pp. 567-570 (2005) S. J. Clark, M. D. Segall, C. J. Pickard, P. J. Hasnip, M. J. Probert, K. Refson, M. C. Payne

[2] J. P. Perdew, K. Burke, M. Enzerhof (1996). Generalized Gradient Approximation Made Simple, Phys, Rev. Lett., 77, 3865.

[3] H. J. Monkhorst and J. D. Pack (1976), Special points for Brillouin-zone integrations, Phys. Rev. B 13, 5188.

Exploring process energetic and tilt angles of adsorption of OH on Pt(111).

Author: Hoang (Bolton) Tran

Problem Statement [1]

In this study, the preferred hydroxyl (OH) adsorption configuration on Pt (111) surface was determined by using Density functional theory (DFT) [2].

Firstly, adsorption site was considered. The top mono-layer of Pt (111) had 4 distinct high symmetry adsorption sites, namely atop, bridge, hcp hollow and fcc hollow (Fig. 1). Secondly, the angle of O-H bond with respect to surface normal (z-axis) was optimized and analyzed for each adsorption site.

Atop Site
Bridge Site
Hollow HCP
Hollow FCC

Figure 1. Visualization of OH adsorbed at each high symmetry sites. Blue is Platinum, red is Oxygen and white is Hydrogen. Cell size is extended for better visualization.

Methodology

Software:

DFT was implemented in CASTEP, a plane wave basis set software embedded inside BIOVIA Material Studio (MS). [3]

Calculation:

1. Calculation setting

  • Exchange-Correlation Functional: Generalized Gradient approximations – Perdew-Burke-Ernzerhof (GGA-PBE) [4]
  • Self-Consistent Field tolerance: 2E-06 eV
  • Pseudo potential: On-the-fly Generated (OTFG) ultrasoft (Vanderbilt-type) [5].
    • Pt: Core radius: 2.4 a.u.. Core orbitals: 1s2 2s2 2p6 3s2 3p6 3d10 4s2 4p6 4d10. Atomic Pseudo Energy: -13042.3015 eV
    • O: Core radius: 1.1 a.u.. Core orbitals: 1s2. Atomic Pseudo Energy: -431.8855 eV
    • H: Core radius: 0.6 a.u.. Core orbitals: None. Atomic Pseudo Energy: -12.4592 eV
  • Monkhost-Pack [6] grid size: 10x10x1 (50 irreducible kpoints)
  • Plane-wave energy cutoff: 450 eV
  • Self-consistent dipole correction in z-direction (normal to surface).
  • Geometric Optimization parameter:
    • Convergence criteria: 1E-05 eV/atom energy | 0.03 eV/Å force |0.001 Å displacement.
    • No cell optimization.
    • BFGS algorithm [7]

2. Pt (111) slab adsorption:

A primitive cell of FCC Platinum was first optimized for lattice constants using DFT. The (111) slab was cleaved from the optimized cell with 4 layers thickness. A vacuum of 12Å was set between each periodic surface in the z direction (surface normal) . No supercell was created for this slab. This means that with one OH attached to any site inside the unit cell, the surface coverage is effectively 1 ML.

Firstly, the bare slab was geometrically optimized to determine its relaxed energy \(E_{slab}\). The bottom two layers of the slab was fixed in position, in order to simulate bulk behavior. Both Michaelides [8] and Seitsonen [9] performed almost identical study and reported that 3-layer-slab with only one free top layer is sufficient for convergence of adsorption energies.

Secondly, a OH radical was optimized alone to determine \(E_{OH}\). Although it might not be physical in term of having a OH radical in the gas/vacuum space, it was merely used as a reference energy to calculate the adsorption energy  \(E_{ad}\) (one could simulate dissociation of water with H co-adsorption or with dissociation energy instead).

Finally, the optimized OH was set on each of the 4 adsorption sites. Aside from the bottom 2 layers of the slab, the O atom was also constrained in space to avoid transitioning to other adsorption sites (although this was proven to not work as expected later on). The H atom was allowed to be relaxed in terms of position, which meant the O-H angle w.r.t surface normal \((\theta)\) would also be relaxed. Geometric optimization was carried out to calculate the energy of the OH/slab system \(E_{OH/slab}\).

The adsorption energy was calculated as

\begin{equation} E_{ad}=E_{OH/slab}-E_{OH}-E_{slab} \end{equation}

This geometric optimization would yield the local minimum of \(\theta\) and also the Pt-O bond length  \((L_{Pt/O})\), since surface Pt atom could still move. Therefore, initial values of \(\theta\)  and  \(L_{Pt/O}\) are important (Fig. 2).

Figure 2. Visualization of tilt angle \((\theta)\) and Pt-O bond length \((L_{Pt/O})\) of OH at atop site

  • \(L_{Pt/O}\) obviously must not be too close to overestimate the energy due to repulsion (since O atom is fixed could not bounce off). However, it should be close enough to form a bond. It has been experimentally reported [8] that optimal Pt-O bond is around 2.11Å. Therefore, \(L_{Pt/O}\) was initially set to be around that value.
  • \(\theta\) ranges from 0º to 90º. However, 0º was avoided since it is quite symmetric (not perfectly because the adjacent Pt sites are different) in the xy-plane, which means it maybe a local minimum and \(\theta\) may not change during optimization. To sample the angles and ensure that there is only one optimal \(\theta\) between 0º to 90º (though in theory, there should be only one), two initial \(\theta\) which are 30º and 60º, were used at each adsorption site.

All energy was calculated as cohesive energy, which is the total 0 K energy given in CASTEP subtracted by the total pseudo-potential energy.

Results and discussion

The cohesive energy of Pt(111) bare slab was -33.35 eV and that of bare OH radical was -7.32 eV. The results are presented in Table 1 below.

Table 1. Adsorption energy, tilt angle and Pt/O bond length at different adsorption sites

Adsorption SitesInitial \(\theta\)Adsorption Energy (eV)Optimized \(\theta\)Optimized \(L_{Pt/O}\) (Angstrom)Stable site?
Atop30\(^o\)-3.0273.4\(^o\)1.975Yes
60\(^o\)-3.0273.5\(^o\)1.975Yes
Bridge30\(^o\)-2.1468.3\(^o\)2.210Yes
60\(^o\)-2.1468.7\(^o\)2.210Yes
FCC Hollow*30\(^o\)-1.0945.6\(^o\)-No
60\(^o\)-1.0945.2\(^o\)-No
HCP Hollow**30\(^o\)-1.2056.3\(^o\)-No
60\(^o\)-1.2056.9\(^o\)-No

*FCC Hollow site has 3 \(L_{Pt/O}\) between O and 3 surrounding Pt fixed at 2.140, 2.134 and 2.133 Å

**HCP Hollow site has 3  \(L_{Pt/O}\) between O and 3 surrounding Pt fixed at 2.244, 2.240 and 2.272 Å

For the two hollow adsorption sites, attempts to geometrically optimize actually resulted in Pt surface atoms shifting and create a bridge adsorption site. Therefore the hollow sites are not stable. To really simulate the hollow configurations then, both O atom and the three adjacent surface Pt were fixed, effectively creating fixed bonds. Therefore optimization for the hollow configuration was solely in terms of angle \(\theta\).

Above point brings up a potential pitfall from fixing Oxygen atom in all 3 directions. Surface Pt atoms were still allowed to relaxed and might protrude upward spuriously. However, optimized \(L_{Pt/O}\) showed that the distance that Pt moved was within 0.1 Å which is considered to be normal for surface relaxation [1]. Retrospectively, it is clear that Oxygen atom should be fixed in only the x and y direction and allowed to move in the z-direction to avoid this pitfall.

The atop site, with lowest adsorption energy, is the predicted preferred adsorption site for OH. This preferred atop site for 1 ML coverage is supported by Michaelides [8] through DFT calculation and Seitsonen [9] through both DFT calculations and experimental LEED measurements.

For adsorption energy, Michaelides reported value of -2.5 eV. For \(L_{Pt/O}\), Michaelides reported 1.99 Å while Seitsonen reported 2.03 Å through DFT calculations (2.11 Å from experiment). For \(\theta\),  Michaelides reported \(73^o\) while Seitsonen reported \(75.8^o\). The literature values are in decent agreement with this study (-3.02 eV, 1.98 Å , \(73.5^o\)). The discrepancy can be attributed to different irreducible kpoints, energy cutoff, slab layers and relax-able layers used in the two papers. Fixing oxygen atom in all 3 directions in this study may also contribute to this discrepancy.

There are two interesting points of discussion from this result:

1. Pt-O bond and O hybridization:

In atop configuration, the O atom is sharing electrons with H atom and only 1 Pt atom. This fact, together with \(\theta\) being very similar to that of water (Pt-O-H of \(104^o\)), suggest that oxygen atoms are in sp3 hybridization state while being adsorbed.

At bridge site, \(\theta\) is different, but the Pt-O-H angle is still roughly \(104^o\) (Fig. 3), suggesting that oxygen also forms sp3 hybridization here.

At both the hollow sites, however, one would expect \(\theta\) to be 0 in order to form 4 bonds (with H and with 3 surrounding Pt) as is characteristic of sp3. Results clearly shows that \(\theta\) is not 0 (though it may stay at 0 if initially built that way), because H wants to form hydrogen bonds with neighboring O. That may suggests the reason why OH does not adsorb favorably on hollow sites, because it cannot both form sp3-type bonds and have non-zero (\theta\) at the same time.

Another interesting comparison can be made with the adsorption of a bare O atom. Without having to share an electron with H atom, the bare O atom actually prefers to form bonds with 3 other Pt atoms (hollow sites) [see Angela’s post].

2. Hydrogen bond network:

In any configuration, H atom wants to form hydrogen bonds with other neighboring O atoms. It wants to move down in z-direction (let’s call it “nodding”) and effectively, lower \(\theta\) because it gives the shortest distance with respect to neighboring O atoms (illustrated in fig. 3). Ideally, it would want to be perfectly parallel to the surface \((\theta=90^o)\) to form the strongest hydrogen bonds it could.

However, back to table 1, it can be seen that \(\theta\) is biggest in the atop site, but not quite 90º. This again goes back to sp3 hybridization of oxygen, due to steric hindrance of Oxygen’s electron pairs, it wants to maintain about \(104^o\) between each pair (i.e. the Pt-O-H angle), so it cannot “nod” all the way down.

Figure 3. Visualization of the hydrogen bonds (dotted red lines) between H atom and neighboring O atoms at optimized atop site (top) and bridge site (bottom). In both cases, Pt-O-H angle is roughly \(104^o\) indicating sp3 hybridization.

In conclusion, it is suggested that OH prefers to be on the atop site because it can form the strongest hydrogen bond network (highest \(\theta\)) without breaking its sp3 hybridization.

Limitations

There are many limitations to the scope of this study in terms of how “real” the system is as compared to the many interesting electrochemical or heterocatalytic systems.

  1. This system modeled the adsorption energy based on the adsorption of a OH radical from gas phase, which was clearly not very realistic. As mentioned, a more realistic approach would be to consider the dissociation of water which leads to co-adsorption of both H and OH. The adsorption of both OH and water molecule is also possible [8,9].
  2.  The environment of the slab was greatly oversimplified. There was no acidic/alkaline solvent (i.e. solvation effects) or applied surface electrostatic potential being considered, which would potentially change the preferred adsorption site of OH [10].
  3. Being a DFT calculation, thermal effect and entropy are often neglected. Adsorption normally leads to negative changes in entropy, which are more significant at high temperature (\(T\Delta S\)). That leads to weaker adsorption at higher temperature.
  4. 1 ML coverage may not be very probable in real system. Lower coverage would lead to different hydrogen bond effect and consequently different adsorption preferences [8,9].

Future study (3rd post)

Possible ideas that could be explored using DFT in later post:

  • Measure frequency of O-H bond stretching at different \(\theta\) and use it to quantify the strength of hydrogen bonds.
  • Lower surface coverage (0.25 ML) to observe its effect on hydrogen bond strength, \(\theta\) and energetic.
  • Float a few water molecules/cluster above the adsorbed OH and observe its effects on hydrogen bond strength, \(\theta\) and energetic.

References

[1] Sholl, D., Steckel, J., & Sholl. (2011). Density Functional Theory (p. 46). Somerset: Wiley.

[2] R.O Jones (2015), Density functional theory: Its origins, rise to prominence, and future, Rev. Mod. Phys. 87, 897.

[3] “First principles methods using CASTEP”, Zeitschrift fuer Kristallographie 220(5-6) pp. 567-570 (2005) S. J. Clark, M. D. Segall, C. J. Pickard, P. J. Hasnip, M. J. Probert, K. Refson, M. C. Payne

[4] J. P. Perdew, K. Burke, M. Enzerhof (1996). Generalized Gradient Approximation Made Simple, Phys, Rev. Lett., 77, 3865.

[5] D. Vanderbilt (1990), Soft self-consistent pseudopotentials in a generalized eigenvalue formalism. Phys. Rev. B 41 (11), 7892-7895.

[6] H. J. Monkhorst and J. D. Pack (1976), Special points for Brillouin-zone integrations, Phys. Rev. B 13, 5188.

[7] J. D. Head and M. C. Zerner (1985), A Broyden-Fletcher-Goldfarb-Shanno optimization procedure for molecular geometries, Chem. Phys. Letters, 122 (3), 264-270.

[8] A. Michaelides and P. Hu (2001), A density functional theory study of hydroxyl and the intermediate in the water formation reaction on Pt, J. Chem. Phys. 144, 513.

[9] A. P. Seitsonen, Y. Zhu, K.Bedurftig and H. Over (2001), Bond mechanism and atomic geometry of an ordered hydroxyl overlayer on Pt(111), J. Am. Chem. Soc. 123, 7347 – 7351.

[10] D. Strmcnik, M. Uchimura, C. Wang, R. Subbaraman, N. Danilovic, D. Vlilet, A.P. Paulikas, V.R. Stamenkovic and N.M. Markovic (2013), Improving the hydrogen oxidation reaction rate by promotion of hydroxyl adsorption, Nature Chem. 5, 300 – 306.

Predicting optimal crystal structure of Pt and its lattice constants using DFT- CASTEP

Author: Hoang (Bolton) Tran

Problem statement [1]

Platinum (Pt) crystal structure are hypothesized to take the form of either a simple cubic (SC), a face-centered cubic (FCC) or a hexagonal closed pack (HCP) structure.

Density functional theory [2] (DFT) calculation is used to first, find the optimal lattice constant for each of the proposed structure. The optimal lattice constant is one that yields the lowest cohesive energy (E) per atom in the unit cell.  Then, the lowest cohesive energy among the three structures would indicate which structure is optimal for Pt.

Methodology

Software:

DFT calculation is carried out using CAmbridge Serial Total Energy Package (CASTEP), a plane wave basis set tool that embedded inside BIOVIA Material Studio (MS) software. [3]

Calculation:

1. DFT energy calculation setting:

  • Exchange-Correlation Functional: Generalized Gradient approximations – Perdew-Burke-Ernzerhof (GGA-PBE) [4]
  • Self-Consistent Field tolerance: 2E-06 eV
  • Pseudo potential: On-the-fly Generated (OTFG) ultrasoft (Vanderbilt-type) [5]. Core radius: 2.4 a.u.. Core orbitals: 1s2 2s2 2p6 3s2 3p6 3d10 4s2 4p6 4d10

Cohesive energy is calculated by subtracting the atomic pseudo energy (-13042.3015 eV) from the estimated 0 K energy. All energy is in per atom basis.

2. Lattice constant optimization

Crystal structures of Pt are generated in MS with editable lattice constant.

The lattice constant for SC and FCC is a which is side length of the cubic unit cell. a is varied and the cohesive energy (E) is calculated. Data of E vs. a is then fitted to the Birch-Murnaghan (BM) equation of state [6]. Then the minimum from the fitted function is considered as the optimal lattice constant.

The lattice constants for HCP are a and c, the basal side length and height respectively. All lattice constants are defined and optimized in its primitive cell.

Simple cubic
Face centered cubic
Hexagonal closed pack

3. Convergence

The energy cutoff for plane wave set is varied for a particular set of lattice constant optimization. As energy cutoff is incremented, the minimum cohesive energy is calculated until the difference with respect to result at highest energy cutoff is less than 0.01eV . This optimal energy cutoff is then used throughout all calculations.

The kpoints, which represent Monkhorst-Pack sampling space [7] for the plane wave, is optimized for each crystal structure. The grid parameters for k space is varied for each set of lattice constant optimization and then optimized similarly to energy cutoff, i.e. until minimum cohesive energy does not change by more than 0.01 eV with respect to highest kpoints result.

Results

1. Energy cutoff optimization:

Simple cubic structure is chosen to be optimized in term of energy cutoff at constant kpoints of  120 (15x15x15 grid). The lattice constant a is varied from 2 to 3.5 Å and cohesive energy is extracted and fitted to BM equation. Similar procedure is then repeated for different values of energy cutoff. Figure 1 illustrate such fitting for two values of energy cutoff. Python’s scipy.optimize module [8] is used to find the minimum energy for each fitted function.

Fig 1. Finding minimum cohesive energy for different energy cutoffs, SC structure.

Table 1 below shows that the minimum cohesive energy goes down as energy cutoffs increases. Difference with respect to highest energy cut off (i.e. 450 eV in this case) is used as the criteria of convergence. Therefore, 400 eV point would be considered as converged because it differs from 450 eV point by less than 0.01 eV. But as to ensure good convergence throughout, 450 eV will be used as the optimal energy cutoffs for all following calculations.

Table 1. Minimum total energy at different energy cutoffs – SC

Energy cutoff (eV)Minimum total energy (eV/atom)Difference with respect to highest cutoff (eV)
270-7.6260.539
300-8.0030.161
350-8.1510.013
400-8.1640
450-8.1640

2. Simple cubic

To determine to optimal lattice constant, similar process of varying a and calculating cohesive energy and fitting it to BM equation to find the minimum point is repeated, now with different kpoints. The results are shown in table 2.

Table 2. Minimum cohesive energy at different kpoints – simple cubic

kpointsMinimum cohesive energy (eV/atom)Difference with respect to highest kpoints (eV)
35-8.194-0.040
56-8.176-0.022
84-8.209-0.055
120-8.157-0.003
286-8.1540.003

It is seen that kpoints of 120 (15x15x15 grid) is the optimal value which gives the optimal lattice constant a as 2.637 Å and minimum cohesive energy as -8.157 eV/atom

3. Face centered cubic

Kpoints convergence is again performed for FCC structure and results are shown in table 3.

Table 3. Minimum cohesive energy at different kpoints – FCC

kpointsMinimum cohesive energy (eV/atom)Difference with respect to highest kpoints (eV)
56-8.659-0.004
84-8.6430.012
110-8.6550.000
120-8.6550.000

The optimal kpoints as seen from the table, is 110 (10x10x10 grid) which correspond to lattice constant as 2.811 Å and minimum cohesive energy as -8.655 eV/atom.

4. Hexagonal closed pack

With two lattice constants a and c, the optimal values are found by fixing a ratio of c/a then perform similar curve fitting for minimum cohesive energy per atom.

First, ratio c/a of 1.50 is chosen to optimize in terms of kpoints.

Table 4. Minimum cohesive energy per atom at different kpoints – HCP

kpointsMinimum cohesive energy (eV/atom)Difference with respect to highest kpoints (eV)
3-8.691-0.026
14-8.681-0.016
36-8.6650
76-8.6640.001
136-8.6650

With similar interpretation as above sections, kpoints of 36 (9x9x6 grid) is optimal in performing calculation for HCP structure. This kpoints will be used for different c/a ratio.

Figure 2 illustrate the BM equation fitting for different c/a ratio. The minimum for each fitted curve is then plotted against the c/a ratio in Figure 3.

Fig. 2: Fitting of BM equation – HCP, varying c/a ratio

Fig. 3: Minimum cohesive energies at different c/a ratio for HCP structure

It is observed from Figure 3 that c/a ratio of 1.65 has the lowest cohesive energy/atom which is -8.674 eV/atom and lattice constants of = 2.817 Å and = 4.648 Å.

Conclusions

It is interesting to see that the minimum cohesive energy per atom of HCP structure is lower than that of FCC structure, which suggest that Pt prefers the HCP structure and that does not agree with what is well known [9].

It is noted, however, that the minimum cohesive energy for each configuration is found by curve fitting to the BM equation. Hence to reaffirm the result, re-calculation was done in CASTEP for each of the crystal structure, using the determined optimal kpoints and lattice constants to reevaluate this discrepancy (table 5).

Table 5. Re-calculation for energy of each crystal structure

SCFCCHCP
a (A)2.6372.8122.817
c (A)--4.468
Fitted minimum E (eV/atom)-8.157-8.655-8.674
Re-optimized E (eV/atom)-8.207-8.664-8.605

From re-calculation results, cohesive energy per atom of FCC structure is the lowest, indicating the expected result of FCC structure being the stable structure for Pt.

Retrospectively, the convergence criteria for both kpoints and energy cutoff is ±0.01 eV, therefore a meaningful comparison can only be made if difference in energy between two structures is more than 0.02 eV. The difference between fitted minimum E of HCP and FCC is 0.019 eV (third row table 5), while that of the re-calculated energy is 0.059 eV (fourth row table 5). Therefore more confidence can be given to the re-calculated results.

A method that does not require re-calculation is to fit the BM equation multiple times to converge on a good cohesive energy value. This study, instead, calculated the minimum cohesive energy after fitting the BM equation only once to roughly determine the optimal lattice constant, then re-calculate the energy at that optimal lattice constant to use for comparison between crystal structures.

Lastly, the optimal lattice constant reported in literature [9] is 3.92 Å. For this study, converting primitive to conventional FCC yields a lattice constant of 3.98 Å. The relative error is 1.5% which is a very good agreement.

Reference.

[1] Sholl, D., Steckel, J., & Sholl. (2011). Density Functional Theory (p. 46). Somerset: Wiley.

[2] R.O Jones (2015), Density functional theory: Its origins, rise to prominence, and future, Rev. Mod. Phys. 87, 897.

[3] “First principles methods using CASTEP”, Zeitschrift fuer Kristallographie 220(5-6) pp. 567-570 (2005) S. J. Clark, M. D. Segall, C. J. Pickard, P. J. Hasnip, M. J. Probert, K. Refson, M. C. Payne

[4] J. P. Perdew, K. Burke, M. Enzerhof (1996). Generalized Gradient Approximation Made Simple, Phys, Rev. Lett., 77, 3865.

[5] D. Vanderbilt (1990), Soft self-consistent pseudopotentials in a generalized eigenvalue formalism. Phys. Rev. B 41 (11), 7892-7895.

[6] Franis Birch (1947), Finite Elastic Strain of Cubic Crystals, Phys. Rev. 71, 809.

[7] H. J. Monkhorst and J. D. Pack (1976), Special points for Brillouin-zone integrations, Phys. Rev. B 13, 5188.

[8] Jones E, Oliphant E, Peterson P, et al. SciPy: Open Source Scientific Tools for Python, 2001.

[9] Hermann, K. (2011). Crystallography and surface structure (Appendix E: Parameter Tables of Crystals). Wiley‐VCH Verlag GmbH & Co. KGaA.