Significant Quantum Effects in Hydrogen Activation

Dissociation of molecular hydrogen is an important step in a wide variety of chemical, biological, and physical processes. Due to the light mass of hydrogen, it is recognized that quantum effects are often important to its reactivity. However, understanding how quantum effects impact the reactivity of hydrogen is still in its infancy. Here, we examine this issue using a well-defined Pd/Cu(111) alloy that allows the activation of hydrogen and deuterium molecules to be examined at individual Pd atom surface sites over a wide range of temperatures. Experiments comparing the uptake of hydrogen and deuterium as a function of temperature reveal completely different behavior of the two species. The rate of hydrogen activation increases at lower sample temperature, whereas deuterium activation slows as the temperature is lowered. Density functional theory simulations in which quantum nuclear effects are accounted for reveal that tunneling through the dissociation barrier is prevalent for H2 up to ∼190 K and for D2 up to ∼140 K. Kinetic Monte Carlo simulations indicate that the effective barrier to H2 dissociation is so low that hydrogen uptake on the surface is limited merely by thermodynamics, whereas the D2 dissociation process is controlled by kinetics. These data illustrate the complexity and inherent quantum nature of this ubiquitous and seemingly simple chemical process. Examining these effects in other systems with a similar range of approaches may uncover temperature regimes where quantum effects can be harnessed, yielding greater control of bond-breaking processes at surfaces and uncovering useful chemistries such as selective bond activation or isotope separation.

T he making and breaking of bonds involving hydrogen atoms at the surfaces of materials plays a major role in nature. It is also relevant to heterogeneous catalysis, sensing, energy storage, and chemical reactions in the interstellar medium. 1À6 Therefore, understanding the adsorption and reaction dynamics of hydrogen on solid surfaces is of widespread importance and central to process improvement and materials engineering. The light mass of hydrogen makes it subject to quantum behavior, which can strongly affect its reactivity. 7 It has been demonstrated both experimentally and theoretically that quantum tunneling plays a role in a variety of reactions, including conformational inversions, eliminations, and enzymatic reactions involving hydrogen. 7À10 While it is known that quantum effects can have a significant effect on the rate of reactions involving hydrogen, our understanding of the regimes in which these quantum effects are important is currently not well developed.
Heterogeneously catalyzed hydrogenation reactions often employ metals and alloys based on platinum, ruthenium, palladium, or nickel. These metals allow for the facile dissociation of molecular hydrogen due to a relatively low reaction barrier. In contrast, coinage metals such as Cu and Au exhibit a higher barrier for the dissociation of molecular hydrogen and are less reactive toward hydrogen. 11,12 Many surface science studies have examined the quantum behavior of hydrogen. 13À30 Experimental work has focused on the quantum tunneling-enabled diffusion of single hydrogen atoms on a variety of surfaces using field emission microscopy, 15 helium spinÀecho spectroscopy, 16 optical diffraction, 17À19 and scanning tunneling microscopy (STM). 20À22 Recently, quantum tunneling was found to allow 2D self-assembly of hydrogen atoms at low temperature (5 K) on Cu(111). 23 The key observable difference between classical and tunneling reaction pathways is the temperature dependence of the reaction rate, and divergence from classical behavior is often detected through careful examination of kinetic isotope effects in chemical reactions. 7,17,27,28,31 We recently reported that addition of just 1% of a monolayer of Pd to the surface of Cu can promote the dissociation of hydrogen. 32À34 These single-atom Pd sites lower the energy barrier to hydrogen uptake, allowing it to dissociate and then spill over onto the bare Cu terraces. 32À34 Despite careful use of density functional theory (DFT) calculations, a number of recent studies could not explain these results. 35À38 Here we use a well-defined model system and an array of sophisticated theoretical approaches that allow us to understand and quantify the precise role of quantum effects in hydrogen dissociation. We studied the uptake of hydrogen and deuterium using temperatureprogrammed desorption (TPD), DFT, path-integralbased DFT approaches to account for quantum nuclear effects, and kinetic Monte Carlo (KMC) simulations. The dissociative adsorption of molecular hydrogen and deuterium reveals vastly different rates and temperature dependence for the two species. The rate of hydrogen activation is higher at lower sample temperature, whereas deuterium activation slows down as the temperature is lowered. The dramatically different behavior observed in these experiments originates from the different extent of quantum tunneling for H 2 and D 2 . In both cases our calculations indicate that tunneling is relevant. KMC simulations indicate that in the case of H 2 the contribution from tunneling is so significant that the effective barrier is reduced to the point that H 2 uptake becomes dominated by thermodynamic effects. On the contrary, tunneling through the D 2 barrier is not as significant, and the D 2 uptake rate is dictated kinetically by the effective height of the barrier.

RESULTS AND DISCUSSION
The uptakes of hydrogen and deuterium on singleatom alloy Pd/Cu surfaces as a function of temperature were obtained using TPD experiments, as shown in Figure 1a and b. The measurements were performed on the 0.01 monolayer (ML) Pd/Cu(111) alloy surface, which consists of 1% Pd in the form of individual, isolated Pd atoms substituted into the Cu(111) lattice. 32À34 The surface was exposed to either hydrogen or deuterium at a selected sample temperature in the range of 75À193 K. The sample was subsequently warmed or cooled to 83 K, followed by the thermal desorption measurement. The 75À193 K temperature range reflects (i) the minimum temperature to which our sample could be cooled and (ii) the maximum temperature we could heat the sample to without desorption of gaseous hydrogen or deuterium. It is immediately clear that the uptakes of the two isotopes on the 0.01 ML Pd/Cu(111) surface are very different. Critically, the uptake of hydrogen decreases when the temperature is raised, while deuterium's uptake increases at higher sample temperature. Isotope effects that lead to different rates of H versus D reactions are common; however, the completely opposite behavior observed here as a function of temperature is exceptional. Figure 1c shows a natural log plot of the two uptake curves shown in Figure 1a and b as a function of inverse temperature. These plots were used to extract the apparent activation energy for the dissociative adsorption process. In the case of deuterium, the apparent activation energy was found to be 0.045 ( 0.005 eV, and for hydrogen it was À0.023 ( 0.005 eV. The small negative apparent activation barrier of hydrogen on 0.01 ML Pd/Cu(111) suggests that the dissociation process involves more than a single step. Typically, negative apparent activation barriers in Arrhenius plots are indicative of the presence of a precursor state; 26,39,40 however, we will show later in our analysis that it also relates to the thermodynamics of H uptake.
To elucidate the underlying mechanism and differences of the H 2 /D 2 dissociation process on Pd/Cu(111), we performed DFT calculations for the process. The specific exchangeÀcorrelation functional used was the optB88-vdW functional, 41 which accounts for van der Waals (vdW) dispersion forces within the nonlocal vdW-DF scheme. 42 Full details of the computational setup are given in the Materials and Methods section. The classical DFT potential energy barrier for H 2 dissociation was calculated on a Cu(111)-(3Â3) surface with a single Pd atom substitution ( Figure 2a). Along the minimum energy pathway of H 2 dissociation, a weak physisorption well was identified, with H 2 found 1.9 Å above Pd. Starting from this physisorbed state (with a binding energy of about À0.1 eV relative to a H 2 molecule in the gas phase), the HÀH bond dissociation occurs above the Pd atom with an energy barrier of 0.46 eV. The barrier obtained here with the optB88-vdW functional is slightly higher than that obtained in previous DFT calculations using a standard generalized gradient approximation functional, and both barriers are much too high to account for the experimental observations. 35 Bearing in mind that the computed DFT barrier is too high to account for the facile activation of H 2 below 100 K, as observed in experiments, we now consider how quantum nuclear effects may alter our physical picture of the dissociation process. The simplest and traditional first step is to consider how zero-point ARTICLE Figure 2. (a) Potential energy surface for H 2 dissociation from the physisorbed precursor state on the Pd/Cu(111) substrate. The total energy of the clean surface and the H 2 in the gas phase is used as the energy zero. Insets are top and side views of initial physisorbed precursor state (H 2 *), transition state (TS), and final state (2H*). White, brown, and cyan spheres indicate H, Cu, and Pd atoms, respectively. (b) Temperature dependence of the effective quantum energy barrier (relative to the gas phase) obtained from harmonic quantum transition state theory calculations that take into account tunneling through the chemisorption barrier and zero-point energies. It can be seen that below the quantum crossover temperature (260 K for H 2 , 190 K for D 2 ) the effective energy barrier decreases as the temperature is lowered. The orange/pink regions of the graph correspond to two different mechanisms for chemisorption below 80 K for H 2 or 50 K for D 2 . These two routes to chemisorption are the precursor-mediated process (c) or below 80 K for H 2 (50 K for D 2 ) a direct chemisorption for H 2 incident on a Pd site (d). ARTICLE energy (ZPE) will alter the energetics of the process. When we consider the effect of ZPE corrections, we do indeed find that they lower the energy barrier to dissociation. However, the reduction is small and the barriers remain too high at 0.37 eV for H 2 and 0.39 eV for D 2 , as compared to the experimentally observed apparent barriers (À0.023 and 0.045 eV for H 2 and D 2 , respectively). This indicates that the experimental results cannot be explained through a simple consideration of the ZPE effect and that instead quantum mechanical tunneling is likely to play an important role in this process.
To account for quantum tunneling effects in this system, we have used two path-integral-based approaches. Specifically, most calculations reported were obtained with harmonic quantum transition state theory (HQTST) (otherwise known as instanton theory), 43 although the more sophisticated ab initio path integral molecular dynamics (AIPIMD) 44,45 approach was also used. The HQTST method is an efficient path-integralbased approach that includes quantum tunneling by considering the spread of a "necklace" of beads over the top of the dissociation barrier. The beads are connected by mass-and temperature-dependent springs, such that the classical limit (with all beads contracted to a single point) is recovered at high temperature and/or high mass. As the temperature is lowered, these springs weaken and the beads sample lower energy states on either side of the classical saddle point, resulting in a lowering of the quantum energy barrier due to tunneling through the classical potential. The HQTST calculations were performed on an analytical one-dimensional potential constructed to resemble our DFT energy profile for dissociative adsorption ( Figure 2a). Figure 2b shows some of the key results of the HQTST calculations wherein the effective dissociation barriers for both H 2 and D 2 are plotted as a function of temperature. It can be seen that below the quantum-classical crossover temperature (which is the temperature below which classical and quantum mechanics diverge; in this system 260 K for H 2 , 190 K for D 2 ) the effective quantum energy barrier decreases as the temperature is lowered. A lower barrier will result in more facile chemisorption of H 2 from the physisorbed state following the precursor mechanism shown in Figure 2c. At 80 K (50 K) the effective quantum barrier for H 2 (D 2 ) to dissociatively chemisorb reaches zero (relative to H 2 or D 2 in the gas phase). This indicates that below these temperatures an incident molecule at a Pd site could undergo barrierless dissociation if, prior to doing so, it does not get trapped in the physisorbed state, as illustrated schematically in Figure 2d.
The HQTST calculations reveal that it is mainly tunneling that leads to the substantial reductions in the effective quantum dissociation barriers. To test this conclusion, we performed a separate set of AIPIMD simulations in which the quantum delocalization of the breaking H 2 bond was examined at various points along the classical dissociation pathway. In AIPIMD the forces on the atoms are computed "on the fly", and it does not rely on a predetermined potential energy surface. A snapshot from an 85 K AIPIMD simulation with the center of mass of the H 2 fixed at the classical saddle point is shown in Figure 3. Clearly there is a large spread of the path integral beads associated with the two hydrogen atoms. Delocalization of the beads along the dissociation reaction coordinate is a signature for quantum tunneling through the potential barrier. The extent of the delocalization observed in the AIPIMD simulations is similar to the spreading observed in the HQTST calculations at the same temperature ( Figure S1).
The path-integral-based DFT calculations reveal significant tunneling through the dissociation barrier at experimentally relevant temperatures, with a greater extent of tunneling for H 2 versus D 2 dissociation. To gain further insight into the experimentally observed differences in the H 2 /D 2 uptake, we performed a series of KMC simulations, making use of the data and insight gained from the various DFT calculations. With no adjustments to the DFT-derived parameters (energetics and reaction prefactors), these simulations ARTICLE predict no quantifiable uptake of H (D) adatoms for the two cases that employ (a) a DFT-derived (classical) activation energy barrier and (b) a temperaturedependent quantum tunneling barrier obtained via the HQTST simulations. Therefore, to capture the experimentally observed trends, systematic adjustments to the model parameters were performed until a reasonable agreement between the experimentally measured and model-predicted surface coverages was reached. These adjustments were small (0.20 eV or less) and were within the typical error bars associated with DFT calculations of processes at surfaces. 46 A more detailed account of these adjustments can be found in the Supporting Information ( Figure S4, Table S1). In particular we found that the coverages predicted by our KMC simulations are very sensitive to the binding energy of H 2 (D 2 ) on the Pd sites. Perturbations as small as 0.005 eV (stabilization) of H 2 (D 2 ) can cause the surface coverages to increase by as much as 20% from their original values, thereby indicating that the binding strength of the precursor state plays a critical role in the overall reaction. Accordingly, the correct treatment of the weak interaction of H 2 /D 2 with the metal surface 47,48 is highly desirable for the detailed understanding of the phenomena studied here. Figure 4 shows the results from the KMC simulations for H 2 /D 2 uptake, plotted against the experimental observations. The surface coverage of H (D) adatoms is found to decrease (increase) monotonically with an increase in temperature, consistent with the experimental trends. Our KMC results suggest that these opposite trends for the surface coverages of H/D adatoms with increasing temperature are due to the fact that H 2 dissociation is thermodynamically controlled, while D 2 dissociation is kinetically controlled over the entire temperature regime.
In the case of H 2 , our KMC simulations indicate that the decreasing trend of H adatom surface coverage with temperature is dictated by the thermodynamics, and not by the kinetics, of H 2 dissociation at surface Pd sites. In other words, the temperature-dependent activation energy barrier for this step is too small under all experimental conditions for it to be of great kinetic relevance, and the H surface coverage is governed by the thermodynamic driving force, i.e., the reaction free energy of this step. Since this is an exothermic step, an increase in temperature is expected to disfavor H 2 dissociation. In particular, the reaction free energy of this elementary step increases with an increase in temperature, which results in a monotonically decreasing equilibrium constant. The decreasing H surface coverage, and hence the negative apparent activation energy barrier, is a direct consequence of this phenomenon.
Contrary to the H 2 case just described, the monotonically increasing trend in D adatom surface coverage is due to the fact that D 2 dissociation is kinetically controlled. The temperature-dependent activation energy barrier for D 2 dissociation at surface Pd sites (Table S1, step 3) is large enough to ensure that it is the kinetically relevant rate-controlling step over the entire  ARTICLE temperature regime. Consequently, the D surface coverage is governed by the kinetics of this elementary step, thereby resulting in a monotonically increasing surface coverage with increasing temperature. In order to quantify this, we obtained the overall rate of D 2 dissociation from our model (total number of D 2 dissociation events, as obtained from our KMC model: "E net "). Plotting the natural log of E net against 1/T (where T = temperature (K)) yields the apparent activation energy barrier of the overall reaction, as shown in Figure 5. This apparent activation energy barrier, as obtained from this analysis, is 0.032 eV, consistent with its experimental value (0.045 ( 0.005 eV).
The results of our KMC simulations suggest that both H 2 and D 2 dissociation proceed through the same reaction pathway. The dissociation of H 2 (D 2 ) to give surface H(D) adatoms occurs exclusively over Pd sites; for the temperature range employed in our experiments, no dissociation is observed over Cu sites because of the prohibitively high activation energy barrier. The KMC simulations suggest that H 2 (D 2 ) must land directly on or adjacent to surface Pd atoms, as the precursor states in the vicinity of Cu atoms are extremely short-lived and prefer to desorb at all temperatures rather than diffuse across the surface to a Pd site. This was verified by reducing the rate constants for H 2 (D 2 ) diffusion steps by several orders of magnitude from the calculated values and observing that the final results were essentially independent of the diffusion rates.
Finally, we note that our calculations indicate that dissociation can occur through two distinct processes, each due to the combination of van der Waals interactions and quantum tunneling. The first process involves the trapping of H 2 in the physisorbed state, which acts as a precursor to the dissociative chemisorption process (Figure 2a and sketched in Figure 2c). The second mechanism may become possible at very low temperatures when quantum effects are so pronounced that the effective chemisorption barrier drops below the gas phase H 2 energy zero (sketched in Figure 2d). Our HQTST calculations predict that this can happen at less than ∼80 K for H 2 and less than ∼50 K for D 2 . This direct dissociation process occurs at lower temperatures than we have been able to probe in the current TPD experiments, but it is an intriguing implication of the calculations that we hope to examine in future studies. More generally, the results of this study suggest a novel approach for inducing selective reactions or isotope separations that utilize the vastly different rates of bond cleavage involving hydrogen and deuterium. For example, many important surfacecatalyzed reactions involve CÀH or OÀH bond cleavage. 49À52 If the C or O atom remains at the same surface site before and after the bond breaks and most of the motion is of the light H atom, similar temperature-dependent activation barriers may occur and novel ways to control bond breaking can be developed.

CONCLUSIONS
We have shown experimentally that the dissociation of H 2 and D 2 molecules at single Pd atom surface sites in a Pd/Cu alloy occurs at vastly different rates and with opposite trends as a function of sample temperature. Traditional DFT-based methods for modeling dissociation barriers are inadequate for explaining these results. Our experimental observations motivated us to take a theoretical approach using path-integral-based DFT simulations that account for the role of quantum nuclear effects in H 2 /D 2 dissociation, coupled with KMC to model the dissociation rates over the relevant range of temperatures. Our theoretical data reveal that the dissociation barriers are strongly temperature dependent and that H 2 /D 2 dissociation is dictated by totally different thermodynamic/kinetic parameters. Within our theoretical framework, we see a positive slope in Arrhenius plots for H 2 uptake, which originates from temperaturedependent quantum tunneling aided by a precursor state, and a negative slope for D 2 , which can be explained by a scheme in which thermally assisted barrier crossing dominates tunneling and the effect of the precursor state. The breaking of chemical bonds involving H atoms is ubiquitous in nature, and our results provide a unique insight into the highly temperature-dependent quantum effects that dictate these rates. We also hope to apply this refined theoretical approach to other light-atom systems in order to search for regimes where important chemistry may become more controllable once quantum effects are fully understood and quantified.

MATERIALS AND METHODS
Experimental Procedures. Measurements were carried out in two UHV chambers, both of which have been described in detail elsewhere. 23,33,34 The first UHV chamber was operated at a base pressure of <1 Â 10 À10 mbar and was equipped with a quadrupole mass spectrometer, used for temperatureprogrammed desorption measurements as well as LEED and AES capabilities. The Cu(111) crystal could be resistively heated to 775 K and cooled to 73 K. Normally the crystal could only be cooled to 83 K with liquid nitrogen, but lower temperatures were reached by bubbling helium gas through the cryogenic dewar. The sample was cleaned by cycles of Ar sputtering (1.5 kV, 11 μA) followed by annealing at 725 K. The second chamber, an Omicron Nanotechnology LT-UHV STM, was operated at a base pressure of <5 Â 10 À11 mbar. The instrument incorporated a preparation chamber in which sample cleaning, annealing, and Pd deposition were performed. Pd deposition was performed in both UHV systems using flux-monitored Omicron Nanotechnology Focus EFM 3 electron beam evaporators with the Cu(111) samples held at 380 K during deposition. Pd coverages were calibrated using both AES and TPD of CO. Hydrogen (99.999%, AirGas) and deuterium (99.999%, AirGas) ARTICLE were used in our study. Hydrogen coverage calculations were based on a saturation coverage of unity when hydrogen was adsorbed on a 5 ML Pd film, assuming that the film terminates as Pd(111).
Theoretical Calculations. Calculations of the classical barrier and pathway for H 2 /D 2 dissociation were performed using the VASP code 53 with projector augmented wave (PAW) potentials. 54 van der Waals dispersion forces were accounted for through the optB88-vdW functional, 41 which accounts for dispersion forces within the nonlocal vdW-DF scheme 42 and has been implemented in VASP by Klimes et al. 55 The electronic states were expanded using plane waves with an energy cutoff of 600 eV. The Pd/Cu(111) surface was modeled by a five-layer slab with a (3Â3) surface unit cell with a Cu atom at the top layer substituted by a Pd atom. A vacuum layer equivalent to six atomic layers was used. The Brillouin zone was sampled using a (8Â8Â1) k-point mesh based on the MonkhorstÀPack scheme. 56 The two bottom-most Cu(111) layers were fixed during relaxation. All structures were fully relaxed until the HellmannÀFeynman forces acting on the atoms were smaller than 0.01 eV/Å. The climbing image nudged elastic band (CI-NEB) method 57 was used to calculate the classical H 2 /D 2 activation energy barrier. ZPE corrections to the classical barrier were obtained simply by taking the difference between the sum of real-valued harmonic vibrational frequencies at the transition state and at the physisorbed (initial) state.
Quantum tunneling effects were accounted for using two path-integral-based approaches. Most results reported have been obtained from the harmonic quantum transition state theory method (also known as instanton) on an analytic potential specifically fitted to a DFT potential energy surface in which the atoms in the substrate were fixed at the positions they adopt prior to H 2 dissociation. Fixing the surface in this manner is an approximation, but it also gives a sensible reaction coordinate, as quantum effects frequently result in cornercutting of the classical minimum energy path during reactions.
Here the fixed surface results in a chemisorption barrier of 468 meV, only 6 meV higher than the barrier obtained when the surface is allowed to fully relax during the dissociation process. We used 200 beads in our instanton calculations, which were converged by comparing to 100 and 500 bead test calculations. AIPIMD simulations were also performed using 16 imaginary time-slices (known as beads) propagated by Langevin dynamics at 85 K and forces computed on-the-fly using DFT. For efficiency reasons we used a slightly less expensive computational setup than the barrier height calculations above. In particular, a threelayer slab was used with the bottom layer fixed, and the Brillouin zone was sampled with a (4Â4Â1) k-point mesh. These reduced settings result in a classical chemisorption barrier of 465 meV, which is a very good representation of the fully converged system.
The kmos 58 framework was used for the implementation of kinetic Monte Carlo simulations. From a given state of the system, the KMC algorithm determines the probability of all possible transitions to other states according to the user-defined rates of all available elementary reaction steps. A reaction step is then chosen in accordance with this probability distribution, and the system assumes the new state determined by that reaction step. The simulation has no memory; that is, the probability of undergoing a particular transition is determined only by the current state of the system and not by any prior configuration. The simulation time step after each transition is selected randomly from a Poisson distribution with a mean that equals the sum of the rates of all available transitions from the original state. Simulations are allowed to proceed until the desired simulation time is reached (50 s for H 2 , 200 s for D 2 , same as the dosing times in the experiments).
The code input comprises the lattice structure, binding energies and minimum energy configurations of all the species involved in the reaction mechanism, activation energy barriers and pre-exponential factors for all elementary steps, simulation time, and the reaction temperature and pressure. All simulations were performed on a Pd/Cu catalytic surface represented as a lattice in which all types of atop (Pd, Cu) and 3-fold hollow (Pd/Cu, Cu) sites were explicitly taken into account (see Figure S2). In particular, the lattice structure was represented by a 20Â20 Cu(111) grid with periodic boundary conditions. One percent of the surface Cu atoms were randomly substituted by Pd atoms with the restriction that they could not be nearest neighbors. This is in agreement with studies that have shown that Pd atoms disperse in the surface layer of Cu(111) 59 and that the most stable configurations for alloys with low Pd content (and coverages below X Pd = 0.4) are those where the surface Pd atoms are surrounded by Cu atoms. 32 This is also consistent with STM observations reported in past experimental studies. 60,61 H 2 (D 2 ) molecular precursors were only allowed to adsorb atop surface atoms, while H (D) adatoms occupied only hollow sites, consistent with our DFT calculations indicating that adsorption of H adatoms is ∼0.4 eV more stable on the hollow sites than the atop sites. Our DFT calculations show that adsorption of H (D) adatoms is thermodynamically favorable in up to three hollow sites around a single metal surface atom. Adsorption of a fourth, fifth, or sixth H (D) adatom around a metal surface atom is energetically unfavorable, with substantially weaker binding compared to that on other surface sites. The KMC simulations therefore allowed adsorption in up to three hollow sites around each metal surface atom. Simulations were performed at similar conditions as experiments: the surface coverage of H (D) was recorded after dosing with 10 À10 atm H 2 for 50 s (10 À9 atm D 2 for 200 s) at fixed temperatures in the range 75À195 K. Each KMC simulation was performed a total of 25 times to obtain an average value of final surface coverage for each experimental condition being simulated.
All the model parameters, with the exception of sticking coefficients for H 2 (D 2 ), were rigorously derived from our DFT calculations. The sticking coefficients for the H 2 (D 2 ) adsorption steps were assumed to be 1. The enthalpy estimates for the surface species and transition states were obtained by correcting the electronic energies, as obtained from our DFT calculations, for ZPE contributions and temperature variations. The temperature-dependent activation energy barriers for the H 2 /D 2 dissociation steps on Pd (Table S1, step 3) were obtained from the HQTST calculations. For our modeling purposes, these quantum tunneling barriers were approximated to be linear functions of temperature until they reach their classical values, as shown in Figure S3. The classical activation energy barrier for the H 2 (D 2 ) dissociation step on Cu (Table S1, step 4) was obtained using the CI-NEB method; 57 the barriers for all the surface diffusion steps (Table S1, steps 5À9) were found to be insensitive parameters for the KMC simulations and were fixed at 0.12 eV, in agreement with our earlier studies. 62 The entropy values were directly calculated from the vibrational frequencies of the respective states, and the pre-exponential factors were calculated from entropy differences between the initial and transition states of the respective elementary steps. The rate constants for reactions involving bond formation or bond scission were calculated using transition state theory; collision theory was employed to calculate the rate constants for the adsorptionÀdesorption steps. 63,64 Conflict of Interest: The authors declare no competing financial interest. Supporting Information Available: Supporting Information contains additional details about path integral DFT and KMC simulations. This material is available free of charge via the Internet at http://pubs.acs.org.