Benchmark First-Principles Calculations of ... - ACS Publications

Jan 24, 2018 - Benchmark First-Principles Calculations of Adsorbate Free Energies. Anshumaan Bajpai,. †,§. Prateek Mehta,. †,§. Kurt Frey,. †...
0 downloads 0 Views 5MB Size
Research Article Cite This: ACS Catal. 2018, 8, 1945−1954

pubs.acs.org/acscatalysis

Benchmark First-Principles Calculations of Adsorbate Free Energies Anshumaan Bajpai,†,§ Prateek Mehta,†,§ Kurt Frey,† Andrew M. Lehmer,† and William F. Schneider*,†,‡ †

Department of Chemical and Biomolecular Engineering, University of Notre Dame, Notre Dame, Indiana 46556, United States Department of Chemistry and Biochemistry, University of Notre Dame, Notre Dame, Indiana 46556, United States



S Supporting Information *

ABSTRACT: Adsorbate free energies are fundamental quantities in the microkinetic modeling of catalytic reactions. In first-principles modeling, finite-temperature free energies are generally obtained by combining density functional theory energies with standard approximate models, such as the harmonic oscillator, the hindered translator, or the twodimensional ideal gas. In this work, we calculate accurate free energies directly from first-principles potential energy surfaces combined with exact quantum mechanical solutions for the translational energy states to benchmark the reliability of common approximations. Through a series of case studies of monatomic adsorbates on metal surfaces, we show that no one free energy model performs satisfactorily in all cases. Moreover, even combinations of different approximations sometimes deviate significantly from the free energies calculated by our first-principles approach. Using observations from these case studies, we discuss how a full quantum mechanical approach can be extended to calculate accurate free energies for arbitrary adsorbate potential energy surfaces at computational cost similar to standard models. KEYWORDS: free energy, entropy, density functional theory, adsorption, harmonic oscillator, hindered translator, potential energy sampling

1. INTRODUCTION The past few decades have seen an explosion in the use of density-functional-theory-based models for computational catalyst screening, reaction mechanism discovery, and microkinetic modeling. This explosion has been met with an increasing focus on benchmarking the predictions made by standard DFT functionals against experiment and to higher levels of theory.1−3 For example, the development of functionals that incorporate Bayesian error estimation4,5 have allowed for uncertainty quantification in DFT calculated reaction rates.6−9 Compared to these efforts directed at 0 K energies, the performance of the models used to treat finite temperature contributions to the free energies that ultimately govern adsorption and reaction has been analyzed to a lesser extent. The free energy of an adsorbate-covered surface has a coverage-dependent, configurational contribution that arises from different arrangements of adsorbates among surface sites. A mean-field approach can be used to capture this configurational contribution at low-coverage or in the absence of strong adsorbate interactions, and accurate representations can be extracted from cluster expansion (or similar) techniques combined with Monte Carlo sampling.10−12 Adsorbate free energies also have a local contribution arising from internal rotations and vibrations of individual adsorbates, as well as © 2018 American Chemical Society

translational center of mass (c.o.m.) motion. Popular approaches to estimate the local contributions begin with an assumed functional form for the c.o.m. potential energy surface (PES). For example, in the harmonic oscillator model the PES is assumed to be separable into parabolic potentials along the three dimensions. The translational partition function then becomes a product of harmonic oscillator partition functions with characteristic frequencies taken from the second derivatives of the harmonic potentials. Recent publications highlight that the harmonic approximation underestimates c.o.m. free energies of certain adsorbates, especially at the temperatures at which surface adsorbates are catalytically active.13,14 In particular, evidence indicates that the free energies of physisorbed species may be better represented by a free translator (2D ideal gas) model15,16 or as a scaled fraction of the gas-phase entropy of the molecule.13 The hindered translator model was developed to bridge the gap between the harmonic oscillator and free translator models.17,18 In this model, the corrugated translational PES is approximated with a cosine function and parametrized based on the distance between adsorption sites and the hopping barrier between Received: October 9, 2017 Revised: December 25, 2017 Published: January 24, 2018 1945

DOI: 10.1021/acscatal.7b03438 ACS Catal. 2018, 8, 1945−1954

Research Article

ACS Catalysis sites. The hindered translator model is sometimes invoked when harmonic vibrational frequencies are computed to be smaller than some threshold, e.g., 150 cm−1.19 In this work, we benchmark these common translational free energy approximations against accurate free energies computed directly from first-principles. Using DFT, we compute the full in-plane adsorption PES for H, C, N, O, and S on Pt and Au (100) surfaces. We then solve the single-particle Schrödinger equation numerically to extract the full quantum-mechanical translational densities of states. We use the sum of states to compute the local contributions to the Helmholtz free energies and, by differentiation, adsorbate entropies. We find that the PESs of even these simple monatomic adsorbates are best described by different physical models, highlighting the need to exercise care in choosing a single model for all adsorbates in a system. Moreover, even combinations of different analytical models do not satisfactorily capture free energies for certain systems. Our first-principles approach outperforms conventional models without any parametrization at similar computational cost, and is extensible to arbitrary PESs, including those that describe the c.o.m. motion of polyatomic adsorbates.

Figure 1. Representative DFT-calculated PES for monatomic adsorption on a (100) metal surface facet. Grid points for vertical DFT relaxations indicated with black dots. The finite element meshes used for solving the Schrödinger equation for hollow and bridge site adsorption are shown in red and blue, respectively.

2. COMPUTATIONAL METHODS 2.1. DFT Calculations. Plane-wave, supercell DFT calculations were performed using the Vienna Ab-initio Simulation Package (VASP 5.4.1).20−22 Plane waves were included up to a cutoff of 420 eV. Electron exchange and correlation were described with the Perdew−Wang 91 (PW91) implementation of the generalized gradient approximation (GGA).23 Core electrons were modeled with the projector augmented wave20,22,24 (PAW) method. An order 1 MethfesselPaxton electron smearing of width 0.15 eV was used to set partial band occupancies. Electronic energies were considered converged when the difference between subsequent selfconsistent iterations was less than 10−5 eV, and atomic positions were relaxed until the force on each atom was less than 0.03 eV Å−1. The bulk lattice constants for Pt and Au within these approximations are 3.988 and 4.182 Å, respectively. The surfaces were modeled using five metal layers separated by ∼16 Å of vacuum. A 10 × 10 × 1, Γ-centered kpoint mesh was used for the 2 × 2 × 1 metal surface. 2.2. Construction of Potential Energy Surfaces. Adsorbate-free metal surfaces were constructed in 2 × 2 supercells of the (100) facet by relaxing the top three surface layers of a five layer slab, with the remaining layers fixed in their bulk positions. Interactions between adsorbates and their periodic images are expected to be of the order 0.01 eV,12,25 small enough not to influence the results here. In-plane adsorbate PESs were computed by rastering an adsorbate on a 11 × 11 grid within the cell, as illustrated in Figure 1. At each grid point, the adsorbate was relaxed in the z direction at fixed xy coordinates. The periodicity and symmetry of the supercell is used to obtain the complete in-plane PES. We assume that the potential normal to the surface is separable, similar to previous works.17,18,26 Further, in a series of calculations we found that the potential normal to the surface was well described by a parabola for all adsorbates (Figure S1). Therefore, we treated the surface normal degree of freedom within the harmonic limit (see section 2.3). The normal potential was obtained by placing the adsorbate at the minimum in the in-plane PES and displacing along the z axis in steps of 0.1 Å within a range of ±0.3 Å about the minimum energy location above the surface.

Generally the coupling between adsorbate and metal surface degrees of freedom are small, and accordingly we fixed the metal surface in constructing all adsorbate PESs. We validated this assumption by comparing surface vibrational frequencies with and without the adsorbate. The comparison reveals that the total energy of the vibrational modes of the surface atoms are effectively unperturbed (∼1 meV/atom) in the presence of an adsorbate. 2.3. Adsorbate Free Energy Calculations. We seek here to compute the free energy Aads relative to the ground-state of noninteracting adsorbates at a concentration of one per site. This quantity appears in expressions commonly used to describe the free energies as a function of surface concentration, and we discuss this connection in section S1 of the Supporting Information (SI). Temperature-dependent adsorbate Helmholtz free energies Aads (eq 1) and entropies (eq 2) are calculated using five different approaches: first-principles potential energy sampling (PESq), semiclassical potential energy sampling (PESsc), the harmonic oscillator (HO) approximation, the free translator (FT), and the hindered translator (HT). For each model, we construct the canonical partition function as the product of a model-specific 2D partition function, qxy and a model-independent harmonic oscillator partition function qHO for motion along the surface z normal. A ads = −kBT ln(qxyqzHO) Sads = −

∂A ads ∂T

(1)

(2)

In the PESq approach, the in-plane partition function is calculated as qxyPESq =



Ei ⎞ ⎟ ⎝ kBT ⎠

∑ exp⎜− i=1

(3)

where the accessible energy states, Ei, are the eigenvalues of the time-independent Schrödinger equation. 1946

DOI: 10.1021/acscatal.7b03438 ACS Catal. 2018, 8, 1945−1954

Research Article

ACS Catalysis −h2 2 ∇ Ψ = (E − V )Ψ 8π 2m

(4)

q HO =

The eigenvalues were computed numerically using the MATLAB finite element module, pdeeig27 (an example MATLAB file is included in the SI). The in-plane potential, V, is calculated at each point on the finite element mesh by cubic interpolation of the explicit DFT calculated PES at the grid points shown in Figure 1. We center the computational box around the PES minimum and set its area to that of a unique adsorption site, i.e., a 1 × 1 cell for hollow site 1 1 adsorption and a 2 × 2 cell for bridge site adsorption (see Figure 1). Wave functions are assumed to vanish at the site boundaries. Similar boundary conditions are used for the FT and HT models (described in more detail below). For computational efficiency, we only include eigenvalues smaller than 0.5 eV. The contribution of states beyond this energy cutoff to the free energy were found to be negligible in the temperature range considered in this work. Free energies can also be extracted from the PES by treating the in-plane adsorbate motion as a classical continuous system. This approach has been used in a recent study of CO oxidation by Jørgensen et al.,26 where in-plane partition function was calculated as qxyPESsc =

2πmkBT h2



⎛ − V (x , y ) ⎞ exp⎜ ⎟ dx dy ⎝ kBT ⎠

∏ i=x ,y,z

−hνiHO 2kBT

( ) 1 − exp( ) exp

−hνiHO kBT

(6)

We calculated vibrational frequencies νHO from the second i derivatives of the parabolic energy profiles at the energy minima, obtained by displacing the adsorbate from the minimum energy position along x and y directions with a step size of ±0.01 Å. The hindered translator model relaxes the assumption that the in-plane potential is entirely parabolic by using a corrugated, cosine-like potential (see dashed line in Figure 2), V (x , y) = V0 +

y V bx ⎛ 2πy ⎞ 2πx ⎞⎟ Vb ⎛⎜ ⎜1 − cos ⎟ + 1 − cos 2 ⎝ a ⎠ 2 ⎝ a ⎠

(7)

Here a is the distance between neighboring sites, V0 is the potential energy at the minima, and Vxb and Vyb are the diffusion barriers in the x and y directions. Following Sprowl and co-workers,18,28 we use a ZPEcorrected 2D partition function for the hindered translator:

qxyHT

( )exp⎡⎣⎢−

M =

(5)

πrx Tx

rx ⎤ ⎡ 1 ⎤ 2 rx ⎥exp⎣⎢ − Tx ⎦⎥I0 2Tx Tx ⎦ 2

(1 − exp⎡⎢⎣− ⎤⎥⎦) 1 Tx

( ) exp⎡⎢ ⎣

⎤ 2 ⎥ (2 + 16rx)Tx ⎦

The quantum harmonic oscillator approximation was used to compute qz. Because qxy goes to zero at T = 0, this semiclassical model does not recover the adsorbate zero-point energies. It is expected to converge to the quantum mechanical solution h when the thermal wavelength, λ = , is much smaller

(8)

This partition function scales linearly with the number of sites M and is a function of a dimensionless barrier height, rx = Vb/ HT HT hνHT is a x , and translational energy, Tx = kT/hνx . Here ν characteristic frequency which, from the second derivative of eq

2πmkBT

than the dimensions of the box. As such, the accuracy of this approach is dependent on the temperature, the mass of the adsorbate, and the size of the box. We will now briefly summarize free energy models that assume an analytical expression for the PES (illustrated in Figure 2). In the harmonic oscillator approximation, the adsorbate experiences a parabolic potential in the three dimensions, whose quantum mechanical partition function is given by

7 at the energy minima, is given by ν HT =

Vb 2ma 2

. For 4-fold

hollow site adsorption, diffusion barriers are located along the symmetric normal vibrational coordinates, such that Vxb = Vyb = Vb. For bridge site adsorption, the potential is asymmetric around the adsorption site, and the smallest diffusion barriers do not occur along the normal vibrational coordinates (see discussion in section 3.4). In eqn 8, the last exponential factor is a zero-point energy correction to Hill’s hindered translator

( ) is the zero order modified Bessel

formulation,17 and I0

rx 2Tx

function of the first kind. The HT model has been shown to approach the HO model for kBT ≪ Vb and FT model for kBT ≫ Vb for small organic molecules on Pt(111).18 The free translator (or 2D gas) partition function (eq 9b) can be computed as a Boltzmann sum of the quantum mechanical energy states Enx,ny (eq 9a) available to the adsorbate in a zero potential 2D box of edge length a. In the literature,17,26,29 it is common to use the 2D gas model in the classical limit, where the energy levels are treated as continuous. The summation in eq 9b can then be replaced by an integration to obtain an area (A)-dependent closed form analytical expression.17,29

Figure 2. One dimensional potential energy surfaces for in-plane adsorbate motion as approximated by the harmonic oscillator (HO), the hindered translator (HT), and the free translator (FT) models. a is the nearest neighbor distance between surface atoms for fcc (100) surfaces and Vb is the diffusion barrier.

Enx , ny = 1947

h2 (nx2 + ny2) 8ma 2

nx , ny = 1, 2, 3, ...

(9a)

DOI: 10.1021/acscatal.7b03438 ACS Catal. 2018, 8, 1945−1954

Research Article

ACS Catalysis

Figure 3. In-plane PESs for adsorption of H, C, N, O, and S on Au(100) and Pt(100). The grid points indicate the supercell for solving the Schrödinger equation. Magenta and yellow markers on the supercell grids depict the normal modes for adsorbate vibration. Energies are referenced to the energy minima for the respective adsorbate/metal system.

qxyFT =

⎛ −Enx , ny ⎞ 2πmkBT ⎟⎟ ≈ A h2 ⎝ kBT ⎠

maximum at the bridge site. Harmonic potentials (solid lines in Figure 4b) fit the energy profile quite well at displacements up to 0.5 Å from the energy minima. At larger displacements, the true energy profile diverges from the parabolic fit. In-plane and surface normal vibrational frequencies are computed to be 403, 403, and 213 cm−1, respectively. Helmholtz free energies as calculated by the various models are plotted as a function of temperature in Figure 4c. In this figure, positive free energies in the limit T → 0 K reflect the contribution of zero-point energy. HO free energies are found to almost exactly match the first-principles PESq free energies over the entire temperature range considered here. The deviation between the HO and PESq models is approximately 5 meV at 1200 K. This agreement between the two models is unsurprising, due to the close match between the harmonic fit and the DFT based PES for energies as high as ∼0.45 eV (Figure 4b). The semiclassical PESsc free energy curve converges to the quantum mechanical one at temperatures greater than 400 K. As noted in section 2.3, the semiclassical approach fails to capture the zero point energy and free energies at low temperatures, as the integral approximation is no longer valid because the energy states are not continuous. Free energies predicted by the HT model fall below the PES curves over the entire temperature range. Zero-point energies are 25 meV lower in the HT model compared to the PESq or HO models. At 1200 K, AHT deviates from APES by 140 meV. The origin of the underprediction of free energies by the HT model is related to the failure of the cosine potential to capture the curvature of the DFT PES. The hindered translator potential is less steep around the energy minimum: its vibrational frequencies of 208 cm−1 are 195 cm−1 less than the corresponding HO frequencies. The zero potential free translator model provides a lower bound on the free energies and is clearly not an appropriate choice for the steep potentials of group 1. The temperature-dependent adsorbate entropies predicted by the various models are presented in Figure 4d. Similar to the free energy trends, the harmonic oscillator approximates the PESq entropy very well; deviations are 0.45 eV 4-fold preferred, Vb < 0.15 eV bridge site preferred

magnitude of the diffusion barriers. The members of group 1, including C/Au, N/Au, S/Au, C/Pt, and S/Pt, prefer 4-fold adsorption and have diffusion barriers >0.45 eV. Group 2 includes O/Au and N/Pt, which also prefer 4-fold adsorption but have diffusion barriers 0.03 eV), the error exceeds 100 meV at temperatures above 800 K. The error is smaller at lower T and flatter potentials. As we found in section 3.3, using the width of the deeper well, ab (see Figure 8a) may result in a better agreement between the HT and PESq models. Accordingly, in Figure 8c, b we plot ΔAaHT against temperature for the scaled potentials. The HTab model shows the opposite behavior to the standard b is close to zero. HT model. For larger diffusion barriers, ΔAaHT The error becomes increasingly positive for flatter potentials, where the adsorbate can move easily across the barrier, especially at higher temperatures. The above discussion suggests that there is an optimal width opt = 0 for any value of Vb . We between a and ab at which ΔAaHT

Schrödinger equation. The in-plane partition function can then be approximated as qxy = qxqy (10) Free energies computed using the separability assumption for a group 1 system (C/Au) are compared against the full 2D PESq free energies in Figure 9a. For understanding the trends in

Figure 9. Parity plot comparison of free energy predicted by models with the separability assumption and that predicted by the twodimensional model PESq model for (a) scaled C/Au and (b) scaled N/ Pt.

model performance, one low (400 K) and one high temperature (1200 K) case is presented, and the potential is scaled similar to the approach described above. The HT free energies are also presented as triangles for comparison. We find a near perfect agreement between free energies computed by the 1D (APES1D2) and 2D (APESq) PES sampling approaches at both temperatures and across all vibrational frequencies. In contrast, the HT model shows large deviations from the parity line at higher frequencies and temperature. b A similar comparison of APES1D2, AHT, AaHT is made with APESq for the double welled potentials of N/Pt (group 2) in Figure 9b. At low temperatures, deviations from the parity line are not pronounced, but all three 1D models deviate from the parity line at high temperatures. APES1D2 and AHT deviate from APESq at b larger vibrational frequencies. AaHT shows the opposite trend, deviating from APESq at low vibrational frequencies (as also shown in Figure 8). Figure 9 shows that, if potentials are assumed to be separable (as is commonly assumed within the HO and HT models), numerical solution of the one-dimensional Schrö dinger equations will always be more reliable than any simplified model. The computational cost for the HT and the PES1D2

kBT

find the relationship between

aopt a

( ) to be sigmoidal

and ln

Vb kBT

(Figure 8d). This plot shows that the application of the HT model for double well potentials is not straightforward. In general, at smaller Vb the nearest-neighbor distance should be kBT

used, and smaller values are more appropriate with increasing Vb . kBT

Next, we address the question of separability of the in-plane potential. In the HO and HT models, it is common to assume that the potential is separable in the planar dimensions. We can make a similar approximation in the first-principles sampling approach, by using the one-dimensional potentials along the normal mode coordinates to compute eigenvalue solutions for 1952

DOI: 10.1021/acscatal.7b03438 ACS Catal. 2018, 8, 1945−1954

ACS Catalysis



methods is similar, nudged elastic band31,32 type calculations that are used to compute Vb in the HT model also provide the full energy profile along the diffusion coordinate. Rather than simply using the diffusion barrier and expressing the potential as an analytical expression (eq 7), using the complete 1D energy profile and solving for the eigenstates numerically provides more accurate free energies. For potentials that are not separable, both the HT and the PES1D2 approaches produce similar results, and a full 2-D solution is required to get accurate free energies, especially at higher temperatures. The same potential energy sampling and state summing approach is applicable to lower-symmetry facets or sites, but such potential energy surfaces may not be as readily separable as are the (100) facet potentials. Lastly, we note that the free energy models developed here can be applied to energy landscapes computed with any DFT functional or other means. We expect the adsorbate energy landscapes and thus free energy corrections to be less sensitive to this choice than are absolute binding energies. The models developed here provide a mechanism to test these and similar questions.

Research Article

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acscatal.7b03438. Detailed discussion of configurational contribution to free energy, potential energy profiles along the surface normal, free energy plots for the remaining systems in group 1, group 2, and group 3, density of states for all the systems. (PDF) Folder with matlab scripts to solve Schrödinger equation in 1 and 2 dimensions for a given PES data set. Includes sample data and README files. (ZIP)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Anshumaan Bajpai: 0000-0003-3567-7217 Prateek Mehta: 0000-0001-6233-8072 William F. Schneider: 0000-0003-0664-2138 Author Contributions §

A.B. and P.M. contributed equally to this work.

5. CONCLUSIONS

Notes

We have described an approach for accurate calculation of translational adsorbate free energies, by using quantum mechanical energy states obtained by numerically solving the single-particle Schrödinger equation on complete first-principles potential energy surfaces. We present a series of case studies of monatomic adsorbates where this approach is used to benchmark free energies calculated using standard analytical approximations. Our overall conclusion is that the performance of standard models depends on the shape and separability of the adsorbate PES. The commonly used harmonic approximation is reliable when the adsorbate is confined within a deep potential well, especially at lower temperatures. The HO model completely breaks down for flatter PESs, as energy states available to the adsorbate are no longer localized within the adsorption site. On the other hand, the hindered translator model works well for relatively flat PESs, but exhibits increasing error for steeper potentials. Our analysis identifies an intermediate regime where both the HO models and HT models have errors that are similar in magnitude but of opposite sign, suggesting that the practice of switching between the two models at some threshold vibrational frequency could produce up to 1 order of magnitude error in calculated rate or equilibrium constants. Further, we show that for complex and/or nonseparable PESs (e.g., double-well PESs or bridge site PESs) application of the HT model is not straightforward and may require reparameterization depending on the temperature and the magnitude of the diffusion barrier. In contrast, potential energy sampling based approaches (both quantum mechanical and classical) can be used to evaluate the free energies accurately for arbitrary PESs. The advantage of the quantum mechanical approach is that it is accurate across all temperatures, while the classical approach fails for low temperatures. Finally, we show that for separable potentials, calculation of the complete PES can be avoided, and one-dimensional eigenvalue problems can be solved numerically to obtain accurate free energies at the same computational cost as that of the hindered translator method.

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the U.S. Department of Energy, Office of Science, Basic Energy Sciences, Catalysis Science Program, under Award DE-FG02-06ER15839. Computational resources were provided by the Notre Dame Center for Research Computing. P.M. acknowledges support through the Eilers Graduate Student Fellowship for Energy Related Research.



REFERENCES

(1) Wellendorff, J.; Silbaugh, T. L.; Garcia-Pintos, D.; Nørskov, J. K.; Bligaard, T.; Studt, F.; Campbell, C. T. Surf. Sci. 2015, 640, 36−44. (2) Bligaard, T.; Bullock, R. M.; Campbell, C. T.; Chen, J. G.; Gates, B. C.; Gorte, R. J.; Jones, C. W.; Jones, W. D.; Kitchin, J. R.; Scott, S. L. ACS Catal. 2016, 6, 2590−2602. (3) Sharada, S. M.; Bligaard, T.; Luntz, A. C.; Kroes, G.-J.; Nørskov, J. K. J. Phys. Chem. C 2017, 121, 19807−19815. (4) Wellendorff, J.; Lundgaard, K. T.; Møgelhøj, A.; Petzold, V.; Landis, D. D.; Nørskov, J. K.; Bligaard, T.; Jacobsen, K. W. Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 85, 235149. (5) Lundgaard, K. T.; Wellendorff, J.; Voss, J.; Jacobsen, K. W.; Bligaard, T. Phys. Rev. B: Condens. Matter Mater. Phys. 2016, 93, 235162. (6) Medford, A. J.; Wellendorff, J.; Vojvodic, A.; Studt, F.; AbildPedersen, F.; Jacobsen, K. W.; Bligaard, T.; Norskov, J. K. Science 2014, 345, 197−200. (7) Christensen, R.; Hansen, H. A.; Vegge, T. Catal. Sci. Technol. 2015, 5, 4946−4949. (8) Christensen, R.; Hummelshøj, J. S.; Hansen, H. A.; Vegge, T. J. Phys. Chem. C 2015, 119, 17596−17601. (9) Deshpande, S.; Kitchin, J. R.; Viswanathan, V. ACS Catal. 2016, 6, 5251−5259. (10) Piccinin, S.; Stamatakis, M. ACS Catal. 2014, 4, 2143−2152. (11) Bray, J. M.; Schneider, W. F. ACS Catal. 2015, 5, 1087−1099. (12) Bajpai, A.; Frey, K.; Schneider, W. J. Phys. Chem. C 2017, 121, 7344−7354. (13) Campbell, C. T.; Sellers, J. R. V. J. Am. Chem. Soc. 2012, 134, 18109−18115. 1953

DOI: 10.1021/acscatal.7b03438 ACS Catal. 2018, 8, 1945−1954

Research Article

ACS Catalysis (14) Campbell, C. T.; Á rnadóttir, L.; Sellers, J. R. V. Z. Phys. Chem. 2013, 227, 1435−1454. (15) Savara, A.; Schmidt, C. M.; Geiger, F. M.; Weitz, E. J. Phys. Chem. C 2009, 113, 2806−2815. (16) Réocreux, R.; Ould Hamou, C. A.; Michel, C.; Giorgi, J. B.; Sautet, P. ACS Catal. 2016, 6, 8166−8178. (17) Hill, T. L. An Introduction to Statistical Thermodynamics, 1st ed.; Addison-Wesley: Boston, MA, 1960. (18) Sprowl, L. H.; Campbell, C. T.; Á rnadóttir, L. J. Phys. Chem. C 2016, 120, 9719−9731. (19) Choksi, T.; Greeley, J. ACS Catal. 2016, 6, 7260−7277. (20) Kresse, G.; Furthmüller, J. Comput. Mater. Sci. 1996, 6, 15−50. (21) Kresse, G.; Furthmüller, J. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 11169−11186. (22) Kresse, G.; Joubert, D. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758−1775. (23) Perdew, J. P.; Chevary, J. A.; Vosko, S. H.; Jackson, K. A.; Pederson, M. R.; Singh, D. J.; Fiolhais, C. Phys. Rev. B: Condens. Matter Mater. Phys. 1992, 46, 6671−6687. (24) Blöchl, P. E. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953−17979. (25) Schmidt, D. J.; Chen, W.; Wolverton, C.; Schneider, W. F. J. Chem. Theory Comput. 2012, 8, 264−273. (26) Jørgensen, M.; Grönbeck, H. J. Phys. Chem. C 2017, 121, 7199− 7207. (27) MATLAB and Partial Differential Equation Toolbox. 2015b. (28) Campbell, C. T.; Sprowl, L. H.; Á rnadóttir, L. J. Phys. Chem. C 2016, 120, 10283−10297. (29) Mcquarrie, D. A.; Simon, J. D. Physical Chemistry, A Molecular Approach, 1st ed.; University Science Books: Sausalito, CA, 1997. (30) Ayala, P. Y.; Schlegel, H. B. J. Chem. Phys. 1998, 108, 2314. (31) Henkelman, G.; Uberuaga, B. P.; Jónsson, H. J. Chem. Phys. 2000, 113, 9901−9904. (32) Henkelman, G.; Jónsson, H. J. Chem. Phys. 2000, 113, 9978− 9985.

1954

DOI: 10.1021/acscatal.7b03438 ACS Catal. 2018, 8, 1945−1954