A flat surface study of the Eley-Rideal dynamics of recombinative desorption of hydrogen on a metal surface

Mats Persson

Department of Applied Physics, Chalmers University of Technology and Göteborg University, S-412 96, Göteborg, Sweden

Bret Jackson

Department of Chemistry, University of Massachusetts, Amherst, Massachusetts 01003, USA (January 19, 1995)

Abstract

The dynamics of a direct reaction between a gas phase H atom and an adsorbed H atom, often referred to as an Eley-Rideal mechanism, is explored using a fully three dimensional flat surface model for Cu(111). The model is based on a flat-surface approximation for a single electronically adiabatic potential energy surface (PES). This reduces the inherently six dimensional reactive scattering problem (for a rigid surface) to a three dimensional problem by introducing three constants of motion. The resulting scattering problem is treated quantum mechanically by solving the time-dependent Schrödinger equation, and also by quasi classical trajectory calculations. We have considered four different model PESs which are all more or less compatible with available knowledge about the interactions between hydrogen atoms and molecules and Cu(111). We have studied the dependence of the reactive cross section and product translational energy and rovibrational state distributions on the kinetic energy and angle of the incident H atom, and the vibrational state of the adsorbed H atom.

I. INTRODUCTION

The detailed mechanisms behind gas-surface reactions are in many cases most convincingly revealed by dynamical studies. A particularly interesting class of such reactions is the direct reaction between a gas phase species and an adsorbed species, with the product molecule immediately desorbing. This so called Eley-Rideal (ER) reaction mechanism has been speculated to occur on surfaces for a long time, but only recently has strong evidence been provided for this mechanism. This evidence comes from two characteristic dynamical features of the reaction mechanism; i) the large exothermicity that leads to unusually large translational and internal energies in the desorbing products, and ii) their sensitivity to the initial state of the gas-phase species. This behaviour should be contrasted with the other class of idealized reactions which proceed via the so called Langmuir-Hinschelwool (LH) mechanism, where the reactants are both initially adsorbed on the surface and in thermal equilibrium with the substrate. These reactions are often thermoneutral because the strong molecular bond which is formed generally has an energy close to the sum of the weaker adsorbate-substrate bonds. Moreover, since the reactants are assumed to be thermally equilibrated with the surface there is no memory of the state of the incident gas-phase species in the state of the product molecule.

The first dynamical evidence for an ER mechanism was provided by measurements of the internal state distributions of product molecules. In two separate studies of H atom sources by Hall et. al. [1] and Eenshustra et. al. [2], in which H_{2} is dissociated on a hot Tungsten filament in a cell, highly vibrationally excited H_{2} was observed. The amount of vibrational excitation was so large that it could only be consistent with the energetics of an ER recombination, where one of the H atoms was adsorbed on the cold metal wall. However, in one of the experiments the observed vibrational energy was more than available for this mechanism, making the interpretation less clear. The molecular beam study by Rettner [3], examining the reaction of gas phase H(D) with D(H) adsorbed on Cu(111), was performed under more well-defined conditions. An analysis of the HD product translational energy and angular distributions demonstrated their sensitivity to the incident energy as well as to which isotope was incident. The spread of the translational energy distribution was found to range from zero to about the maximum available energy that is released in the ER mechanism for this system.

For the reaction of a gas phase H with Cl adsorbed on Au(111), Lykke and Kay [4] were the first to suggest, from a quasi-effusive beam study of the product HCl rovibrational distributions, that this reaction proceeded by both an ER and a LH mechanism. The rotational distributions in the ν = 0 channel obeyed Boltzmann statistics at the various surface temperatures, and were consistent with a LH mechanism for the recombination of reactants that are thermally equilibrated with the surface. The rotational distribution in the ν = 1 channel, on the other hand, were fairly insensitive to the surface temperature, and the scattered flux was more sharply peaked towards the surface normal than for ν = 0, suggesting an ER mechanism. This conclusion has been confirmed by Rettner and Auerbach [5] in a detailed molecular beam study of the correlations between the initial state of the reactants and the final state of the product. Specifically, the two mechanisms are reflected in the bimodality of the measured velocity and angular distributions of the HCl product, where the ER components of these distributions depend on the kinetic energy of the incident H atom. A similiar correlation between the translational energy of the product molecule and the incident energy of the gas phase reactant was first established in the protonation of hyperthermal N(C_{2}H_{4})_{3}N molecules at a Pt(111) surface, which demonstrated that atom abstraction can occur via an ER mechanism [6]. Most recently, Rettner [7] has complemented his molecular beam study with detailed and extensive measurements of the HCl internal state distributions that not only confirm the findings of Lykke and Kay but also bring up several new and interesting issues about ER dynamics, such as the disposal of the reaction energy, parallel momentum conservation, and the variation of the mean translational energy with final angle. Gates and co-workers have examined the reaction between gas phase H(D) and D(H) adsorbed on Si(100) and Si(111), and have suggested an ER mechanism based on kinetics studies and measured activation energies [8,9].

These experimental studies of ER dynamics have stimulated a few groups to perform theoretical studies based on fully quantum mechanical models [10–13]. All of the models are based on a restricted collinear geometry involving two degrees of freedom, and use a simple London-Eyring-Polanyi-Sato (LEPS) description of the potential energy surface. Both Kratzer and Brenig [10] and Jackson and Persson [11] (henceforth referred to as JP) found unusually large vibrational excitation of the H_{2} product from the reaction of a gas phase H with an adsorbed H, in qualitative agreement with the results from the gas-cell experiments [1,2] mentioned above. Our early work in JP and subsequent papers [12] has also stressed the important effects of isotopic substitution in hydrogen recombination, and their rationalization in terms of a simple kinematical model. Later experimental observations [3] of the isotope effects on the translational energy of the product molecule were in agreement with our results. A specific result of our isotope studies was that the reaction probability drops when the mass of the adsorbed species is increased relative to that of the incident atom. The ER mechanism for the reaction between gas phase H and Cl adsorbed on Au(111) is in this respect most interesting and our recent study [13] reveals that despite the large mass difference the reaction is still possible, and involves some interesting dynamics. For example, a long-lived reaction intermediate is formed. The relatively low vibrational excitation of the product HCl could be understood in terms of a PES where a major part of the potential energy release is in the exit channel. Also, the reaction probability was found to increase strongly with the initial Cl-metal vibrational quantum number. This leads to a substrate temperature dependent reaction probability that can account for the low apparent thermal activation energy of about 60 meV for the ER component of the reaction of gas phase H with adsorbed Cl, as observed in the beam experiments [4,7] mentioned above. This effect could also be the reason for the low activation energies in the range 70-90 meV observed by Yates and co-workers [14] in their kinetic study of the reaction between H and different adsorbed halogens on a Si(100) surface.

The conclusions that have been drawn from these low dimensional dynamical studies of ER reaction dynamics have been valuable, as shown above, but are limited in range and should only be viewed as qualitative and as a guide for further studies. For example, the collinear model cannot address the role of product molecule rotational motion or nonzero impact parameter, and reaction cross sections cannot be computed.

In this paper we develop and explore a fully three dimensional model for the ER recombination dynamics of a gas-phase atom with an atom adsorbed on a flat surface. This model enables us to calculate rovibrational distributions of the product molecules and reaction cross sections even at off-normal incidence. lt is based on a flat-surface approximation for a single electronically adiabatic potential energy surface (PES), which reduces the six dimensional reactive scattering problem (for a rigid substrate) to a three dimensional problem by introducing three constants of motion; the normal component of the angular momentum and the lateral components of the linear momentum. The resulting scattering problem is treated quantum mechanically by solving the time-dependent Schrödinger equation using a new extension of the pseudospectral method to curvilinear coordinates for the propagation of the wave function. We have applied this method to the ER reaction of a gas phase H atom with an H atom adsorbed on Cu(111). This system has several attractive features; i) the interaction between hydrogen and copper has turned out to be a model system for the study of the dynamics of activated dissociative adsorption [15–17] as well as physisorption [18], which in turn has generated a wealth of detailed knowledge about the H- and H_{2}-Cu(111) interaction, which can be exploited in the modelling of the PES, ii) the light mass of an H atom might introduce some interesting quantum effects, and iii) the possibility of studying ER dynamics by isotopic substitution. We have considered four different model PESs which are all more or less compatible with the available knowledge regarding H and H_{2} interactions with Cu(111). We have studied the dependence of the reactive cross section and translational energy and rovibrational state distributions on the kinetic energy and angle of the incident H atom and the vibrational state of the adsorbed H atom. In order to elucidate the role of quantum effects on the dynamics we compare our results with quasi-classical trajectory calculations.

The paper is organized as follows. Sec. II presents a detailed account of our model and the different techniques used in the calculations. We begin with a description of our “flat-surface” approximation and our LEPS model for the different PESs in Sec. II A. The wave packet propagation scheme and the final state analysis, together with some computational aspects, are detailed in Sec. II B while the quasi-classical trajectory calculations are described in Sec. II C. The results from our quantum mechanical wave packet and quasi-classical trajectory calculations are presented and discussed in Sec. III. These results include reactive and non-reactive cross sections (Sec. III A), translational energy (Sec. III B) and rovibrational state (Sec. III C) distributions as a function of different initial conditions such as the energy and angle of the incident gas-phase H atom and the vibrational state of the adsorbed H atom. Sec. IV summarizes our work and contains some concluding remarks.

II. THEORY

A. Model

In this section, we define the flat-surface approximation for the potential energy surface (PES) and show how it effectively reduces the six-dimensional problem for the ER recombination of two atoms on a rigid surface to a three-dimensional problem in cylindrical coordinates. We then describe our LEPS construction of such a PES.

1. Flat-surface approximation

We first introduce some necessary notation. The location of the two atoms are specified by two position vectors r_{1} and r_{2} where the subscripts “1” and “2” henceforth refer to the atom initially adsorbed on the surface and the incident atom, respectively. Since we are primarily interested in the state of the product molecule, it is convenient to make a transformation to the center-of-mass position vector R = (m_{1}r_{1} + m_{2}r_{2})/M and the relative position vector r = r_{2} − r_{1} where M = m_{1} + m_{2} is the total mass and the reduced mass is denoted by μ = (1/m_{1} + 1/m_{2})^{−1}.

The flat-surface approximation for the PES reduces the dimensionality of the wave function by introducing three constants of motion. In this approximation the interaction between the atoms and the surface is assumed to be translationally invariant along the surface. This introduces an axial symmetry, since the interaction is also rotationally invariant around an axis normal to the surface. These symmetries imply that the potential energy, V (r_{1}, r_{2}), in terms of the center-of-mass and relative position vectors has the form,

V (r_{1}, r_{2}) = V (ρ, z, Z) ,  (2.1)

where ρ and z together with an azimuthal angle ϕ are cylindrical coordinates for the relative position vector, r = ρ(cos ϕx̂+sin ϕŷ)+zẑ. The translational and axial symmetries introduce two sets of constants of motion; the lateral momentum, P_{∥}, of the center-of-mass motion and the z-component, J_{z}, of the angular momentum of the relative motion. The corresponding “good” quantum numbers are then given by the wavevector K_{∥} = P_{∥}/h̄ and the azimuthal quantum number m_{j} = J_{z}/h̄. For definite values of these quantum numbers, the full sixdimensional wavefunction Ψ(r_{1}, r_{2}) of the two atoms is given by,

Ψ(r_{1}, r_{2}) = Ψ(ρ, z, Z) exp(iK_{∥} · R_{∥})\sqrt{\frac{1}{2π}} exp(im_{j}ϕ),  (2.2)

where \sqrt{1/2π} is a normalization constant associated with the integration over ϕ. Hence, the original six-dimensional problem reduces to a set of decoupled three-dimensional problems parametrized by K_{∥} and m_{j}. In fact, for normal incidence, we need only consider a single three-dimensional wave function, since K_{∥} = 0 and m_{j} = 0 for the initial state. For offnormal incidence the initial state is a superposition of states with different m_{j} and we have to propagate a set of independent three-dimensional wave functions in time.

An effective Hamiltonian for the three-dimensional wave function in (2.2) can now be derived as follows. The kinetic energy operators for the center-of-mass motion and relative motion are given in terms of Z, the constants of motion, and the cylindrical coordinates (ρ, z) as,

K̂_{CM} = −\frac{h̄^{2}}{2M} \frac{∂^{2}}{∂Z^{2}} + \frac{h̄^{2}K_{∥}^{2}}{2M}  (2.3)

K̂_{rel} = −\frac{h̄^{2}}{2μ} (\frac{1}{ρ} \frac{∂}{∂ρ} (ρ\frac{∂}{∂ρ}) + \frac{∂^{2}}{∂z^{2}} ) + \frac{h̄^{2}m_{j}^{2}}{2μρ^{2}} .  (2.4)

The mixing of ρ and ∂/∂ρ in K̂_{rel} is an effect of the curvature of the cylindrical coordinates and is avoided by a standard rescaling of the wavefunction by the square root of the cylindrical volume element:

Φ(ρ, z, Z) = ρ^{1/2}Ψ(ρ, z, Z) .  (2.5)

This results in a quadratic form for the effective kinetic energy operator of the relative motion:

K̂_{rel}^{eff} = −\frac{h̄^{2}}{2μ} (\frac{∂^{2}}{∂ρ^{2}} + \frac{∂^{2}}{∂z^{2}} ) .  (2.6)

The attractive “curvature” term, −h̄^{2}/8μρ^{2}, which results from this rescaling is combined with the kinetic energy term, h̄^{2}m_{j}^{2}/(2μρ^{2}), in (2.4) to give a centrifugal term in the effective potential energy:

V^{eff}(ρ, z, Z) = V (ρ, z, Z) + \frac{h̄^{2}(m_{j}^{2} − 1/4)}{2μρ^{2}} .  (2.7)

The total effective Hamiltonian, Ĥ^{eff},for Φ(ρ, z, Z), excluding the conserved kinetic energy of the lateral center-of-mass motion, is then given by,

Ĥ^{eff} = K̂_{zCM} + K̂_{rel}^{eff} + V^{eff}(ρ, z, Z) ,  (2.8)

where K̂_{zCM} = K̂_{CM}−h̄^{2}K_{∥}^{2}/(2M). An important feature of Ĥ^{eff} is that the simple quadratic form of the kinetic energy operator allows one to use an efficient scheme for the propagation of the wavefunction, as shown in Section II B 1, but with the price that the effective potential is singular at ρ = 0.

2. Model potential energy surface

We assume that the interaction between the surface and the atoms is governed by a single electronically adiabatic PES. Recent progress in total energy calculations based on the local density approximation have made it possible to calculate with reasonable accuracy such a PES for the interaction of hydrogen with metal surfaces. Efforts have so far been devoted to the calculation of the parts of the PES that are relevant to dissociative adsorption, but it should also be possible to calculate the parts relevant to an ER reaction. In this study, we use a LEPS model for V (ρ, z, Z), as in our earlier collinear studies [11]. This PES is defined as,

V (ρ, z, Z) = U_{a}(z_{1}) + U_{a}(z_{2}) + U_{m}(r) −
\sqrt{Q_{m}(r)^{2} + (Q_{a}(z_{1}) + Q_{a}(z_{2}))^{2} − (Q_{a}(z_{1}) + Q_{a}(z_{2}))Q_{m}(r)} ,  (2.9)

where r = (ρ^{2}+z^{2})^{1/2} is the inter-atomic distance, and z_{1} = Z−m_{2}z/M and z_{2} = Z+m_{1}z/M are the z coordinates for atom “1” and “2”, respectively. The potential terms labelled by the subscripts a and m are associated with the single atom and the molecule, respectively, and are given by,

U_{i}(x) = \frac{D_{i}}{4(1 + Δ_{i})}[(3 + Δ_{i}) exp (−2α_{i}(x − x_{i})) − (2 + 6Δ_{i}) exp (−α_{i}(x − x_{i}))]  (2.10)
Q_{i}(x) = \frac{D_{i}}{4(1 + Δ_{i})}[(1 + 3Δ_{i}) exp (−2α_{i}(x − x_{i})) − (6 + 2Δ_{i}) exp (−α_{i}(x − x_{i}))], (2.11)

where i stands for a or m. For ρ = 0 this PES reduces to the collinear model used in JP. This LEPS-PES reduces to a particularly simple form either when the atoms are far away from each other, r → ∞, or when they are both far away from the surface, z_{1}, z_{2} → ∞. In the former situation Q_{a}(z_{1}) + Q_{a}(z_{2}) ≫ Q_{m}(r) while in the latter situation Q_{m}(r) ≫ Q_{a}(z_{1}) + Q_{a}(z_{2}) . In these two limits, V (ρ, z, Z) has the asymptotic forms,

V (ρ, z, Z) = \{ \begin{matrix}V_{a}(z_{1}) + V_{a}(z_{2}) + δV_{m}(r), & if r → ∞ \\ V_{m}(r) + δV_{a}(z_{1}) + δV_{a}(z_{2}), & if z_{1}, z_{2} → ∞ \end{matrix} ,  (2.12)

where V_{i}(x) is a Morse potential with well-depth D_{i}, equilibrium position x_{i} and repulsive exponent 2α_{i}. The asymptotic interaction term is given by,

δV_{i}(x) ≃ \frac{−5D_{i}(Δ_{i} − 0.2)}{8(1 + Δ_{i})} exp(−α_{i}(x − x_{i})) .  (2.13)

Hence, when r → ∞ the interaction is dominated by a Morse potential, V_{a}, for each atom with the surface, while for z_{1}, z_{2} → ∞ the interaction is dominated by a molecular Morse potential V_{m} between the atoms. The Sato parameters Δ_{a} and Δ_{m} describe the form of V (ρ, z, Z) in the interaction region between the atoms and the surface, as is evident from the asymptotic interaction term, δV_{i}(x), in (2.13).

Our motivation for choosing the different potential parameters in Eq.(2.9) is the same as in JP. A good description of the hydrogen-hydrogen interaction is given by V_{m}(r) with the parameter values given in Table I. Our parameters for V_{a}(z_{1,2}) are the same as in JP, which give a hydrogen chemisorption energy, D_{a} ≈ 2.45eV, and vibrational frequency, h̄ω_{a} ≈ 130 meV, in agreement with experimental [19] and theoretical studies [20] of H on Cu(111). In our choice of values for Δ_{a} and Δ_{m}, we are guided by what is known from experimental [15,16] and theoretical studies [17] of the topology of the dissociative part of the PES. These studies suggest that the barrier height V^{∗} ≈ 0.7 eV and is late in the entrance channel. From this gross knowledge, we are not able to determine uniquely Δ_{a} and Δ_{m}. We have therefore considered four different sets that are more or less compatible with this knowledge, and the resulting model PES’s are labelled by the suffixes -I to -IV. There is a caveat in this procedure introduced by our neglect of surface corrugation since the site for the dissociation is not necessarily the same as the equilibrium site for a hydrogen atom on Cu(111). The flat-surface approximation is in part justified by the relatively low diffusion barriers of about a few tenths of an eV for H atoms on metal surfaces. The model PES-I has the same Δ_{a} and Δ_{m} as in the collinear model used in JP, and its dissocative section, V (ρ, z = 0, Z), is shown in Fig. 1. This figure shows that there is a late “activation” barrier located at (ρ^{∗}, Z^{∗}) ≈ (1.16, 0.12) Å with height, V^{∗} ≈ 0.6 eV. Our other model potentials have a similiar topology of the dissociative section, as shown by the corresponding V^{∗} and (ρ^{∗}, Z^{∗}) in Table I.

A characteristic feature of the ER mechanism for recombinative desorption is the large potential energy release, D_{m} −D_{a}, in this reaction. Where this potential energy drop occurs in the PES is revealed by the two-dimensional section, V (ρ, z, Z = z_{0} + m_{2}z/M), where atom “1” is fixed in its initial equilibrium position, z_{0} = x_{a}. As shown in Fig. 1, there is no potential barrier for atom “2” to approach atom “1”, a sizeable part of the potential energy release is in this entrance channel, and the “activation” barrier for dissociation shows up as a saddle point displaced laterally away from atom “1”. A measure of the potential energy release in the entrance channel is given by V_{exo} = −V (ρ = 0, z = r_{0}, Z = z_{0} +m_{2}r_{0}/M)−D_{a}, and is tabulated in Table I for PES-I-IV.

B. Solution of the time-dependent Schrödinger equation

1. Wave-packet propagation scheme

The representation of a reactive scattering problem on a finite grid requires some care since the pseudospectral technique imposes periodic boundary conditions on the wave function. In the asymptotic coordinate directions where we have open channels, this problem is handled by introducing absorbing boundaries, defined by optical potentials. The boundary at ρ = 0, which requires special treatment, is avoided by making an odd extension of the wavefunction to negative ρ:

Φ(ρ, z, Z; t) = −Φ(−ρ, z, Z; t) .  (2.14)

This extension ensures that Φ(ρ = 0, z, Z; t) = 0 and that the grid boundaries of ρ are in the asymptotic directions. The specific form of the grid and the optical potentials will be discussed in Section II B 4.

The initial wave function, Ψ(r_{1}, r_{2}; t = 0), should represent an adsorbed atom, “1”, and an incident atom, “2”, far away from the surface. We assume that atom “1” is in a bound state n of V_{a}(z_{1}) with zero lateral momentum, p_{1∥} = 0. Atom “2” has an angle of incidence θ_{i} and kinetic energy E_{i}. It is represented by a plane wave parallel to the surface with a lateral wave vector k_{2∥} = \sqrt{2m_{2}E_{i}} sin θ_{i}/h̄ and a Gaussian wave function, G(z_{2}), normal to the surface. Ψ(r_{1}, r_{2}; t = 0) has then a definite lateral momentum P_{∥} = h̄k_{2∥} for the centerof-mass motion, and using the well-known cylindrical wave expansion of a two-dimensional plane wave it can now be expanded as,

Ψ(r_{1}, r_{2}; t = 0) = ∑_{m_{j}=−∞}^{∞} ρ^{−1/2}Φ_{m_{j}}(ρ, z, Z; t = 0)i^{|m_{j}|}\sqrt{\frac{1}{2π}} exp(im_{j}ϕ) exp(iK_{∥} · R_{∥}), ρ > 0,

(2.15)

where

Φ_{m_{j}}(ρ, z, Z; t = 0) = ρ^{1/2}u_{n}(z_{1})G(z_{2})J_{|m_{j}|}(k_{∥}ρ)F (ρ), ρ > 0.  (2.16)

J_{|m_{j}|}(k_{∥}ρ) is a cylindrical Bessel function and k_{∥} = μk_{2∥}/m_{2} is the lateral wave vector for the relative motion of the atoms. Note that for off-normal incidence it is sufficient to consider only positive m_{j} while for normal incidence k_{∥} = 0 and we need only to consider a single wavefunction Φ_{m_{j}=0}. u_{n}(z_{1}) is a bound state wavefunction of V_{a}(z_{1}) and we use the analytical form for a Morse potential.

The Gaussian wave function G(z_{2}) is given by,

G(z_{2}) = \frac{1}{\sqrt{πΔz}} exp[−ik_{z}(z_{2} − z_{i}) − \frac{(z_{2} − z_{i})^{2}}{2Δz^{2}} ]  (2.17)

where h̄k_{z} is the average momentum of atom “2” normal to the surface, z_{i} is its average initial position, and Δz is its position spread with a corresponding (root-mean-square) normal kinetic energy spread ΔE_{z} ≃ h̄\sqrt{E_{z}/(2m)}/Δz. The value of k_{z} has been corrected for the small potential energy in the entrance channel via k_{z} = \sqrt{2m_{2}(E_{iz}− < V^{eff} >)}/h̄, where E_{iz} = E_{i} cos^{2} θ_{i} is the average incident energy at infinity and < V^{eff} > is the average of V^{eff} over the initial wave packet. The spatial variation of V^{eff} over the wave packet also introduces an additional energy uncertainity, but this is typically much smaller than ΔE_{z} and can be neglected.

A cut-off function F (ρ) is introduced in order to localize the lateral extension of the wavefunction within the grid and is choosen to have the form,

F (ρ) = \sqrt{\frac{π}{A^{eff}c(m_{j}; k_{∥})}} \{ \begin{matrix}1 & if |ρ| < ρ_{c} \\ exp(−(\frac{|ρ|−ρ_{c}}{λ})^{2}) & otherwise \end{matrix} ,  (2.18)

where A^{eff} is the effective area of the lateral spatial extension of the initial wave function [21]. c(m_{j}; k_{∥}) is the relative contribution from component m_{j} and is defined in such a way that \sqrt{π/A^{eff}c(m_{j}; k_{∥})} is the normalization constant for ρ^{1/2}J_{|m_{j}|}(k_{∥}ρ)F (ρ) over ρ. This definition implies that

c(0; k_{∥}) + 2 ∑_{m_{j}=1}^{∞} c(m_{j}; k_{∥}) = 1 .  (2.19)

Note that c(m_{j} ≠ 0; k_{∥}) ≠ 0 only when k_{∥} ≠ 0. F (ρ) will introduce some kinetic energy, E_{ρ}^{cut}, in the lateral relative motion of the atoms, but for the choice of parameters used in the calculation it is negligible [21].

The time-evolution of each m_{j} component of the wave function, Φ_{m_{j}}, by Ĥ^{eff}, as defined in (2.8), is done numerically using the standard split-time propagator [22]. In this spectral technique, Φ_{m_{j}} is computed at a later time using the Trotter formula with an approximate form for the short time propagator as,

Φ_{m_{j}}(ρ, z, Z; t) = [exp(−i\frac{K̂^{eff}Δt}{2h̄}) exp(−i\frac{V̂^{eff}Δt}{h̄}) exp(−i\frac{K̂^{eff}Δt}{2h̄})]^{N}Φ_{m_{j}}(ρ, z, Z; t = 0) , (2.20)

where Δt = t/N is the time step. The term exp(−iK̂^{eff}Δt/2h̄) has been computed in a standard manner using fast Fourier transform techniques. The potential terms V (ρ, z, Z) and h̄^{2}m_{j}^{2}/(2μρ^{2}) of V̂^{eff} are represented directly on the grid while the singular attractive curvature term −h̄^{2}/(8μρ^{2}) requires special treatment. We find that a good grid represention of this latter term is obtained by keeping the original form,

−\frac{h̄^{2}}{8μρ^{2}} = \frac{h̄^{2}}{2m}[ρ^{−1/2} \frac{d^{2}ρ^{1/2}}{dρ^{2}} ] .  (2.21)

and representing the deratives on the grid via discrete Fourier transforms, where we have defined ρ^{1/2} = sign(ρ)\sqrt{|ρ|}. The quality of this represention will be discussed later in Section II B 3. This kind of method should also be useful for handling the propagation of wave packets in other curvilinear coordinate systems.

2. Final state analysis

The final rovibrational and translational energy distributions of the product molecules have been computed using a method developed independently by Gates and Darling [23], which is a special case of the projection operator approach by Neuhauser and Baer [24]. In this section, we briefly describe how this method has been applied to our problem, and how the cross sections have been calculated.

In this method the reacted part of the wave packet, which is absorbed by the optical potential near the grid boundary of the coordinate Z, is analyzed on a separate asymptotic grid where the time-evolution is governed by the free molecular Hamiltonian. More specifically, at every time step the absorbed part of each m_{j} component of the wave packet is put on an asymptotic (2+1)-dimensional grid in the spatial relative coordinates (ρ, z) and the wave vector, k_{Z}, for the center-of-mass motion in the Z-direction. The Hamiltonian on this grid is the free space molecular Hamiltonian, Ĥ_{m}^{eff} + K̂_{zCM} where now the effective potential energy is V_{m}^{eff} = V_{m}(r)+h̄^{2}(m_{j}^{2} −1/4)/(2μρ^{2}). Only positive values for k_{Z} are allowed on the asymptotic grid in order to represent molecules that propagate away from the surface. At each time step, the wave function on the asymptotic grid, Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; t), is propagated by the free molecular Hamiltonian and combined with the absorbed part of Φ_{m_{j}}(ρ, z, Z; t) in the following way:

Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; t + Δt) = exp(−i\frac{Ĥ_{m}^{eff}Δt}{h̄}) exp(−i\frac{E_{Z}(k_{Z})Δt}{h̄})Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; t)
+ ∫ (1 − d(Z)) exp(−ik_{Z}Z)Φ_{m_{j}}(ρ, z, Z; t)dZ .  (2.22)

Here, d(Z) is the damping function used in the absorbing boundary for the Z-direction. In the propagation we have used the fact that the kinetic energy operator K̂_{zCM} is diagonal on the asymptotic grid with energy eigenvalues E_{Z}(k_{Z}) = h̄^{2}k_{Z}^{2}/(2M).

After propagating the wave function for sufficiently long time T so that most of the reacted part of the wave function resides on the asymptotic grid, the different state transition probabilities can be extracted and resolved in the following way. First, for a homonuclear molecule like H_{2} the para (p) and ortho (o) components of the wave function are obtained by combining the exchange symmetry adapted components of Φ_{asymp} as,

\{ \begin{matrix}Φ_{m_{j}}^{asymp,(p)} \\ Φ_{m_{j}}^{asymp,(o)} \end{matrix} \} (ρ, z, k_{Z}) = \frac{1}{\sqrt{2}} [Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}) ± (−1)^{m_{j}}Φ_{m_{j}}^{asymp}(ρ, −z, k_{Z})]  (2.23)

The reaction probability P_{r}(m_{j}) for p- or o-H_{2} is then given by the norm of the corresponding components of Φ_{m_{j}}^{asymp} times the weight c(m_{j}; k_{∥}). The total reactive cross section, σ_{r}, is now obtained as,

σ_{r} = A^{eff} [c(0; k_{∥})P_{r}(0) + 2 ∑_{m_{j}=1}^{∞} c(m_{j}; k_{∥})P_{r}(m_{j})] ,  (2.24)

where the effective area, A^{eff}, and the weights, c(m_{j}; k_{∥}), are defined as in Eq.(2.18). The factor 2 accounts for the degeneracy with respect to ±m_{j}. The non-reactive cross section for the exchange of atom “1” with atom “2” is obtained in a similiar manner from the norms of the different m_{j} components of the wave function absorbed at the grid boundary in the direction z → −∞.

Now we turn to the computation of the internal state distributions of the product molecules. For each m_{j} channel the probability P (E_{Z}; m_{j})dE_{Z} of finding a molecule with kinetic energy E_{Z}, for center-of-mass motion normal to the surface, is computed from the norm of the spatial part of Φ_{m_{j}}^{asymp} on the asymptotic grid for each k_{Z} as,

P (E_{Z}; m_{j})dE_{Z} = (∫ ∫ dρdz|Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; t = T)|^{2}) \frac{h̄^{2}k_{Z}}{M} dk_{Z} .  (2.25)

The rovibrational state distribution for each m_{j} channel is computed from an autocorrelation function for Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; t) in the following way. The first step is to compute the timedependent wave functions Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; t + T) using the Trotter formula in (2.20), but with a short-time propagator based on the free molecular Hamiltonian Ĥ_{m}^{eff} with zero kinetic energy in the center-of-mass motion. The auto-correlation function C(t) is calculated from the overlap integral of this time-dependent wave function with its initial value at t = T as,

C(t) = ∫ ∫ dρdz ∫_{0}^{∞} \frac{dk_{Z}}{2π} Φ_{m_{j}}^{asymp∗}(ρ, z, k_{Z}; t + T)Φ_{m_{j}}^{asymp}(ρ, z, k_{Z}; T) .  (2.26)

A resolution of the asymptotic wave function in terms of the rovibrational states, labelled by the vibrational level ν and angular momentum j, shows that the time-dependence of C(t) is determined by the rovibrational energies E_{(ν,j)}, and gives directly the probabilities P (ν, j; m_{j}) of finding the molecule in the corresponding state. These energies and probabilities are extracted from the Fourier transform of C(t), which is given by,

C̃(E) = ∫_{0}^{∞} w(t)C(t) exp(−i\frac{Et}{h̄})dt
= ∑_{(ν,J)} w̃(E − E_{(ν,j)})P (ν, j; m_{j}) ,  (2.27)

where w̃(E) = ∫_{0}^{∞} dtw(t) exp(−iEt/h̄) is the half-sided Fourier transform of an appropriate window function w(t). The function C̃(E) consists thus of a series of peaks located at the different E_{(ν,j)}, whose line shapes are determined by w̃(E) and with amplitudes proportional to P (ν, j; m_{j}). How these different peaks are assigned to different rovibrational states, (ν, j), is described in the next section. The total translational and rovibrational distributions are obtained by averaging the corresponding distributions for each m_{j} component of the wavefunction with the same weights, c(m_{j}; k_{∥}), as used in Eq. (2.24).

3. Rovibrational energies and wavefunctions

A crucial point in the final state analysis is the assignment of the peaks in the calculated spectra, C̃(E), to different rovibrational states of the molecule. An unambious assignment has been accomplished by comparing the peak positions with those obtained from a calculation of the autocorrelation functions of wave functions with definite values of j and m_{j}, and in a few instances from their nodal structure.

The most direct method for determining the rovibrational states and energies of a diatomic molecule is to view it as a central potential problem. The angular part of a boundstate wave function is then a spherical harmonic, Y_{l}^{m_{j}}, and its radial part is determined from a one-dimensional Schrödinger equation with an effective potential V_{m}(r)+h̄^{2}j(j+1)/(2μr^{2}). In our wave-packet calculation we have determined the rovibrational energies by performing a separate calculation for the autocorrelation function, C(t), for a two-dimensional wave function propagated in time by H_{m}^{eff} with the initial form,

Φ(ρ, z; t = 0) = N \frac{ρ^{1/2}}{r} g(r)P_{j}^{|m_{j}|}(\frac{z}{r}), ρ > 0 ,  (2.28)

where g(r) is a Gaussian wave function, P_{j}^{|m_{j}|} is an associated Legendre function and N a normalization constant. With appropriate parameters for g(r), the peak positions in the corresponding spectra, C̃(E), give the energies of the different low-lying vibrational states with definite values of (j, m_{j}). The “exact” rovibrational energies, on the other hand, have been obtained by solving numerically the time-independent radial Schrödinger equation. The calculated E_{(ν,j)} (for m_{j} = 0) are typically within a few meV of the “exact” E_{(ν,j)}, except for the highly rotationally excited states (e.g. n = 0, j = 18) where the difference is about 9 meV. The lower accuracy for these states is not unexpected since they have the most rapid spatial variation, particularly in the angular direction.

The proposed wave propagation scheme is also able to handle states with azimuthal quantum numbers m_{j} ≠ 0. The degeneracy of E_{(ν,j)} with respect to m_{j} is only approximately true in the wave-packet propagation scheme in cylindrical coordinates, and the magnitude of the deviations depends on the numerical accuracy of the calculation. The calculated energies are found to be nearly degenerate with deviations being typically less than their deviations from the “exact” values. In fact they tend even in some cases to be in better agreement at larger m_{j}. The behaviour of E_{(ν,j)} with m_{j} is not easy to anticipate since the number of nodes in the polar angle decreases with increasing m_{j}, and the range of the polar angle decreases due to localization of the wavefunction near z = 0 by the increasingly repulsive centrifugal potential.

4. Some computational details

The construction of the grid requires some consideration of the energy-scales of interest and the nature of the interaction potential. The size of the grid-spacing, Δx, along a coordinate x, is determined by the maximum kinetic energy, K_{max} = h̄^{2}π^{2}/(2mΔx^{2}), that must be represented along this coordinate as well as the spatial scale for the variation of the potential. We have used Δρ = Δz =0.1Å and ΔZ =0.05Å , corresponding to K_{max} ≈ 4eV for each coordinate, which is sufficent to accomodate the energies of interest for the reaction of the hydrogen atoms. The difference of a factor of 4 between the masses M and μ for a homonuclear system makes it necessary to decrease ΔZ by a factor of 2 as compared to Δρ and Δz. The number of gridpoints along x, N_{x}, is determined by the size of the interaction region. This region is relatively small for our model potentials, and we find that N_{x} = 128 for each coordinate is enough to accomodate a sufficiently large part of the asymptotic region in the different scattering channels. N_{ρ}^{asymp} = N_{z}^{asymp} = 32 is sufficient for the spatial part of the asymptotic grid since it has only to cover the outgoing molecular channel, while for the part of the grid that represents k_{Z}, N_{Z}^{asymp} is set by the main grid to be N_{Z}/2. When calculating the detailed translational and rovibrational distributions, we use twice as many gridpoints, N_{Z}, in the Z-direction in order to have a sufficiently dense grid in k_{Z}. The optical potential at the absorbing boundaries has been chosen to be linear in x, which corresponds to a simple exponential form for the resulting damping function, d(x). We have used a decay parameter of 1Å^{−1} and a range of about 1.5Å for d(x).

There are some numerical considerations associated with the time-propagation of the wavefunction using the split-time operator, as defined in Eq.(2.20). The necessary time step, Δt, is set by the kinetic energy range of interest in the reaction and not by the maximum total kinetic energy allowed on the grid, which in this case is about 12eV. We have used Δt = 0.25fs which corresponds to a faithful treatment of kinetic energies up to h̄/(Δt) ≈ 2.5eV. We find that it takes less than 4000 timesteps for the reacted part of the wavefunction to leave the main grid and to build up on the asymptotic grid. The number of timesteps in the final state analysis, N_{asymp}, is determined by the energy resolution required to resolve the rovibrational energies of the molecule. We use a standard Hann window function with an energy-resolution ΔE_{rovib} = 2πh̄/(NΔt), and ΔE_{rovib} ≈4meV, corresponding to N = 4000, is found to be sufficient. The rovibrational energies and probabilites have been determined from the calculated C̃(E) by fitting each peak position with a simple quadratic polynomial, which is able to account for 98-99% of the total spectral strength.

In this section we also comment on the treatment of the curvature term. A direct grid representation of this term gives, for example, a (n, j) = (0, 0) wave function that is depleted close to ρ = 0, while the proposed representation, as defined in eq. (2.21), gives a much better wave function. Another illustration of the merits of the latter representation is obtained by considering the expectation value of the kinetic energy of the groundstate wavefunction of a three dimensional harmonic oscillator. For such a wavefunction with vibrational energies h̄ω_{z} = h̄ω_{x} = h̄ω_{y} = 0.2eV centered at the ρ = 0 axis the expectation value of the kinetic energy differs only by about 1% from its exact value, h̄(ω_{z} + ω_{x} + ω_{y})/4 while the direct representation gives an error of about 40%.

The calculations have been performed on a CM200 data parallel machine with 8K processors. The most computationally intensive parts are the time-evolution of the wave function on the main grid and the successive transformations of the reacted part onto the asymptotic grid. For example, on the large 128^{2} × 256 grid as used in the calculation of the detailed product distributions, each time step takes about 6.0s, while the time-propagation, which includes two FFT’s, takes about 2.7s.

C. Quasi-classical trajectory calculations

The calculated reaction cross sections and translational energy and rotational distributions have also been compared with the results from quasi-classical trajectory studies using the same model PES and scattering conditions in order to elucidate the role of quantumeffects in the reaction dynamics. In this section we give a brief description of the calculations behind these quasiclassical results. The prefix “quasi” refers to the fact that it is necessary to take into account the zero-point vibrational motion of the initially adsorbed atom.

The classical counterpart to the Hamiltonian in (2.8) governs the time-evolution of the coordinates ρ, z and Z, and is given by

H^{eff} = \frac{p_{ρ}^{2} + p_{z}^{2}}{2μ} + \frac{p_{Z}^{2}}{2M} + \frac{P_{∥}^{2}}{2M} + V (ρ, z, Z) + \frac{J_{z}^{2}}{2μρ^{2}} ,  (2.29)

where p_{ρ}, p_{z} and p_{Z} are the momenta conjugate to these coordinates. P_{∥} and J_{z} are the lateral center-of-mass momenta and the perpendicular component of the angular momentrum of the relative motion, respectively, and both are conserved in the flat-surface approximation. Note that in this classical Hamiltonian the quantum-mechanical “curvature” term is absent in the centrifugal potential. Given appropriate initial conditions, the equations of motion for the coordinates generated by the Hamiltonian in (2.29) have been integrated in a standard manner by the velocity Verlet algorithm [25]. We have prefered to integrate the individual atomic coordinates z_{1} and z_{2} instead of z and Z.

We have only considered normal incidence and the adsorbed atom is assumed initially to have no lateral momentum, which implies, as in the quantum mechanical case, that J_{z} = 0 and P_{∥} = 0. The strict classical limit of the initial conditions used in the quantummechanical case corresponds to the initial atom “1” being at rest and atom “2” being incident with velocity v_{i} = \sqrt{2E_{i}/m_{2}}, with a distribution of impact parameters, b = ρ(t = 0), with a relative weight b. The vibrational motion of the initially adsorbed atom in a quantummechanical vibrational state ν_{i} has been mimicked in the quasi-classical calculations by having a micro-canonical distribution P (z_{1}, v_{1}; E) for the initial position z_{1} and velocity v_{1} of atom “1” with an energy E equal to the corresponding bound state energy E_{ν}. For the Morse potential V_{a}(z_{1}), P (z_{1}, v_{1}; E) can be sampled effectively because of the analytical form of the corresponding cumulative distribution function. Specifically, the proper initial distribution of z_{1} and v_{1} is obtained by sampling an uniform distribution of ξ in the range [0, 1] as,

z_{1} = −\frac{1}{α_{a}} ln(1 + c sin(2 tan^{−1}(\frac{ς − 1}{ς + 1}))),  (2.30)
ς = \sqrt{\frac{1 − c}{1 + c}} tan(\frac{π}{2} ξ),  (2.31)
v_{1} = ±\sqrt{2(E − V_{a}(z_{1}))/m_{1}}  (2.32)

where c = \sqrt{1 + E_{ν}/D_{a}}.

The final state analysis of the product molecule is straightforward in the classical case. For each impact parameter, b, the reaction probability, P_{r}(b), is determined by the relative number of reactive trajectories. The corresponding quasi-classical reaction cross section is given by

σ_{r} = 2π ∫_{0}^{∞} bdbP_{r}(b) .  (2.33)

For each reactive trajectory far away from the surface, the translational energy, E_{Z} = MŻ^{2}/2, and the angular momentum, J = μ(zρ̇ − ρż) with corresponding quantum number j = J/h̄, are calculated readily for the product molecule. The vibrational distributions have been obtained from the corresponding classical action variable as

h̄(ν + \frac{1}{2}) = \frac{1}{π} ∫_{r_{1}}^{r_{2}} \sqrt{2μ(E_{int} − V_{m}(r) − \frac{J^{2}}{2μr^{2}} )}dr  (2.34)

where r_{1} and r_{2} are the classical turning points for the radial distance r, for a molecule with internal energy E_{int} = μ(ρ̇^{2} + ż^{2})/2 + V_{m}(r) + J_{z}^{2}/(2μρ^{2}). Here we have obtained the vibrational quantum numbers, ν, by equating the action variable with h̄(ν + \frac{1}{2}). Continous distributions for E_{Z}, j and ν have been obtained by introducing a finite resolution via a triangular resolution function. For appropriate values of the resolutions [26], we find typically that these state distributions are converged for about 10000 reactive trajectories.

III. RESULTS AND DISCUSSION

We have performed extensive wave packet and classical trajectory studies of the dependence of the reactivity and the final state of the product H_{2} molecule on the initial state of the reacting H atoms, for four different models of the PES (PES-I,-II,-III, and -IV). The final state of the H_{2} molecule has been characterized with respect to the translational energy and rovibrational distributions. The most detailed studies have been based on the model PES-I. In most of the studies the incoming H atom is at normal incidence, but we have also done a calculation for off-normal incidence. We have varied the kinetic energy of the incident atom and the vibrational state of the adsorbed atom. Our quantum mechanical results for reactive cross sections and final state distributions have also been compared with the corresponding (quasi)-classical results in order to identify and clarify the role of quantum effects in the reaction dynamics.

A. Cross sections

1. Reactive scattering

We find no preference for the formation of p- or o-H_{2} molecules in our wave packet calculations for any of the PES models. The corresponding cross sections, σ_{r}^{(p)} and σ_{r}^{(o)}, are about equal, as shown in Table II. We will henceforth make no distinction between p- or o-H_{2} molecules and just consider the total cross section, σ_{r} = σ_{r}^{(p)} + σ_{r}^{(o)}. The quantum mechanical results in Fig. 2 for σ_{r}, as calculated using the model PES-I, demonstrate that σ_{r} decreases slowly with incident kinetic energy, E_{i}, and is sensitive to the initial vibrational state, ν_{i}, of the adsorbed H atom. The corresponding quasi-classical result gives a good account of the magnitude of the quantum mechanical σ_{r} but misses the enhancement at low E_{i}. The rather weak dependence of σ_{r} on E_{i} in the range considered is not unexpected for this kind of a PES, where a large part of the exothermicity is in the entrance channel. The other PES models also have this topology, and we have therefore limited our investigations of these other models to a single energy, which has been choosen to be E_{i}=75meV. The results for these surfaces, as depicted in Table II, shows that quasi-classical results are in general not in good agreement with the quantum mechanical results. For instance, the quantum and quasi-classical results for PES-III differ by about a factor of 10. The calculated values for the reactive cross sections are markedly smaller than the approximate values 4±2.7Å^{2} and 5.4±1.3Å^{2} which are estimated [27] from the observed reaction probabilities between H(D) and saturated overlayers of D(H) on Cu(111) [3]. Our values are more in line with the measured cross section of about 0.3Å^{2} for the ER component of the reaction of H with Cl on Au(111) at a substrate temperature of 100K [7]

The dependence of σ_{r} on ν_{i}, for PES-I, is more complex than for E_{i} ; it increases first when ν_{i} is changed from 0 to 1 and then decreases for ν_{i} = 2. This trend with ν_{i} is nicely reproduced in the quasi-classical calculations, as can be seen in Fig. 2. The importance of ν_{i} is further illustrated by the fact that a strict classical calculation, in which the adsorbed atom is initially at rest, gives a σ_{r} = 0.24Å^{2}, which is roughly half of the corresponding quasi-classical cross section for ν_{i} = 0. This increase in reactivity with ν_{i} was seen, but in much more dramatic way, in an earlier theoretical study of H(g) + Cl/Au [13], and may be the source of the increased yields observed at higher surface temperatures for several Eley-Rideal reactions [4,7,14].

In a classical description of the scattering it is possible to define a reaction probability, P_{r}(b), as a function of the impact parameter b. We find that this quantity gives useful information about the origin of the calculated magnitudes of σ_{r}. The result for P_{r}(b) using the model PES-I, as shown in Fig. 3a, demonstrates that the reaction can take place in two different ranges of impact parameters; one where the atom is incident more or less on top of the preadsorbed atom, b ≈ 0, and a second where the atom is incident with b ≈ 0.7Å, which is close to the equlibrium distance, r_{0} ≡ x_{m}, of the H_{2} molecule. The dominant contribution to σ_{r} (see Eq. (2.33)) comes from the latter range of b since its area is much larger than for b ≈ 0. Note that the low weight for small b to σ_{r} shows that there is a qualitative distinction between the flat-surface model and the collinear model. Comparison of the quasi-classical result for P_{r}(b) with the corresponding classical result shows that the effect of the initial zero-point vibration is to extend the range of reactive b around r_{0} with a concomitant decrease of P_{r}(b). The net effect is an increase of σ_{r}.

The quasi-classical results for the PES models II,III, and IV, as shown in Fig. 3b, demonstrate that by modifing PES-I, the reactive ranges of b can be changed in different ways. For PES-II, we still have the same bimodal distribution of reactive b as for PES-I, but the range around r_{0} is much more narrow, resulting in a smaller σ_{r} than for PES-I. This reactive range of b around r_{0} disappears completely for PES-III, and the reaction takes place only for b close to 0, resulting in a very small cross section in relative terms, as shown in Table II. For PES-IV, no reaction takes place in this range, and instead the reactive distribution is located around b ≈ 0.5Å. In fact, these different behaviors of the reactive impact distributions of the model PES’s has been the underlying motivation to study them. We shall see how they influence the internal state distribution of the product molecules in Sec. III C.

2. Non-reactive scattering

Several outcomes other than reaction are possible in this model. These non-reactive events includes reflection and trapping of the incident atom on the surface, exchange with the preadsorbed atom, and collision-induced desorption. The collision-induced desorption process is energetically forbidden since its energy threshold is equal to the binding energy (≈ 2.4 eV) of the adsorbed atom and is well above the incident energies of interest here. We find a large cross section for trapping of the atom, while the cross section for the exchange process is found to be relatively small and is typically less than 3 × 10^{−3}Å^{2}. In our wave packet calculations we have not attempted to make a detailed final state analysis of the trapping but it is straightforward, however, to obtain these cross sections from a final state analysis of the trajectories in the quasi-classical calculation. We find that for the PES-I the cross section for trapping of the incident atom is about 18 Å^{2}, at E_{i}=75 meV and ν_{i} = 0. Thus the trapping of a gas phase H atom via energy and momentum transfer to an adsorbed H atom is very efficient. This process can possibly be a source of ”hot” H adatoms [28] that can react with other adsorbed H atoms before equilibrating with the surface. Such a reaction mechanism is intermediate between ER and LH.

B. Final state energy distributions

A central question in this work is how the initial state of the reacting H atoms affects the relative distribution of the energy released in this reaction, Q, between translational and internal energy, E_{Z} and E_{int}, respectively, in the product H_{2} molecule. In our flatsurface model these two energy distributions are directly related to each other via energy conservation; E_{int} + E_{Z} = Q. Q in turn depends on the potential energy release, D_{m} − D_{A}, E_{i}, and ν_{i};

Q = D_{m} − D_{a} + E(ν_{i}) + E_{i}.  (3.1)

Thus it is sufficient to consider the translational energy distribution, P (E_{Z}), only.

We begin by presenting the results for the dependence of the average translational energy, ⟨E_{Z}⟩, on E_{i}, ν_{i} and the shape of the model PES. The calculated results in Fig. 4 for ⟨E_{Z}⟩ using PES-I show that ⟨E_{Z}⟩ increases almost linearly with E_{i} and decreases monotonically with ν_{i}. The slope of ⟨E_{Z}⟩ with E_{i} shows that the efficiacy for transfer of E_{i} to ⟨E_{Z}⟩ is about 0.5. The decrease of ⟨E_{Z}⟩ with vibrational energy, E(ν_{i}), with about 0.13 eV between consequetive ν_{i}, shows that the corresponding efficiacy is negative and approximately -1. In turn this implies directly from energy conservation that the efficiacy for transfer of E_{i} to ⟨E_{int}⟩ is also about 0.5, while the corresponding efficiacy of ν_{i} is positive and is about 2. The quasi-classical results for ⟨E_{Z}⟩ agree rather well with the quantum mechanical results with regard to trends with E_{i} and ν_{i}. The most notable difference is its overestimate of the effect of changing ν_{i} from 0 to 1 which probably reflects the inherent defiencies associated with a quasi-classical description of the initial vibrational state of the adsorbed atom. When considering our different PES models, we find, as seen in Table II, that the quantum mechanical values for ⟨E_{Z}⟩ are rather insensitive to the model PES while there are large differences among the corresponding quasi-classical values. Thus, the quasi-classical results are in general not as accurate as suggested by Fig. 4. Finally, we note that the calculated values for ⟨E_{Z}⟩ are consistent with the measured values of about 0.85 eV and 1.1 eV for the HD molecules produced in the reactions H(g)+D(ads) and D(g)+H(ads), respectively [7].

A more critical comparison between the quantum mechanical and the quasi-classical results is obtained by studying the detailed translational energy distributions, P (E_{Z}). The quantum mechanical result for PES-I, as shown in Fig. 5, gives a broad distribution for E_{Z} of the product molecule. The peaks in this distribution reflect via energy conservation the discrete rovibrational energy levels of the product H_{2} molecule. The gross features of this distribution are well reproduced by the corresponding quasi-classical result, except that it overestimates P (E_{Z}) at large E_{Z}. That the distribution is sensitive to the vibrational state of the adsorbed atom is demonstrated in the same figure by the strict classical result, where P (E_{Z}) is sharply peaked around 1.0 eV, close to the quasi-classical value of ⟨E_{Z}⟩. A common behaviour of P (E_{Z}) among the other PES models, as shown in Fig. 6, is that the quantum mechanical results are broader in energy than the corresponding quasi-classical results and have, in particular, a pronounced high energy tail. This behavior is most notable for the PES-II model.

C. Internal state distributions

An interesting feature of our flat-surface model is the ability to calculate, as detailed in Sec. II B 2, proper rovibrational distributions, P (ν, j), for a 3-D rotor. One should keep in mind that in this model, at normal incidence only molecules which undergo cartwheel motion (m_{j} = 0) will be produced. Before presenting full results for P (ν, j; m_{j} = 0), we present the rotational and vibrational distributions separately; i.e., P (j) = ∑_{ν} P (ν, j; m_{j} = 0) and P (ν) = ∑_{j} P (ν, j; m_{j} = 0), respectively. We end with a presentation of the corresponding results for off-normal incidence, where molecules will also undergo helicopter motion (m_{j} ≠ 0).

1. Rotational distributions

In Figs 7 and 8 we make a direct comparison between quantum mechanical and (quasi)- classical rotational distributions, P (j), of the product H_{2} molecule both as a function of ν_{i} and the shape of our different PES’s. We find no dramatic variations of P (j) with ν_{i} for the model PES-I, as can be seen in Fig. 7. A characteristic feature is the bimodal form of P (j), which is pronounced for ν_{i} = 1. One change with increasing ν_{i} is that the strength of the low-j distribution increases at the expense of the high-j states, as seen in the P (j) distributions for ν_{i} = 0 and ν_{i} = 1, 2. The quasi-classical P (j) distribution reproduces well the width of the quantum mechanical distribution, but fails to reproduce the bimodal behaviour. The maximum quantum numbers for rotational excitation, j_{max}, are set by the energy, Q, released in the reaction and are listed together with Q in Table III. Note that j_{max} is slightly larger in the quasi-classical treatment since it allows for a vibrational state ν = −1/2 that is lower in energy than the quantum mechanical ground state ν = 0. The importance of the initial vibrational state of the adsorbed H atom on P (j) is demonstrated by the strict classical result in Fig. 7a. This distribution is very narrow compared to the quasi-classical result and is peaked around j ≈ 13. This value is close to the quantum number j ∼ 10, corresponding to the angular momentum of an incident H atom with an impact parameter in the range b ≈ 0.75 Å, which has been accelerated by the potential energy release in the entrance channel.

The bimodal behavior of the quantum mechanical P (j) is also present for most of our different PES models, as shown in Fig. 8. The only exception is PES-III, which has the unique feature that in the quasi-classical calculations the reaction takes place only for b close to 0. This behavior suggests that the bimodal form of P (j) is a characteristic feature of a PES where classically the reaction takes place with a sizeable impact parameter away from b = 0. The quasi-classical calculations do not reproduce the bimodal distributions, which suggests that they have a quantum mechanical origin. The deviations are most notable in the P (j) for PES-II, as can be clearly seen in Fig. 8a. In this case the quasi-classical distribution of J is much more narrow than the corresponding quantum distribution and is similiar in behavior to the strict classical P (j) for the model PES-I, as shown in Fig. 7. In the case of the model PES-III, on the other hand, the quasi-classical results reproduce very well the calculated quantum distribution P (j). This distribution P (j) is as broad as for the other PES models, which is rather unexpected given that the (quasi-classical) reactive range of b is close to 0, and that the incident atom should have a negligible J with respect to the adsorbed atom. However, large values for j can be produced by an initially small impact parameter as a results of the deflection of the incident atom from the adsorbed atom when the former is accelerated by the potential energy release in the entrance channel.

2. Vibrational distributions

The quantum mechanical and (quasi)-classical vibrational distributions, P (ν), of the product H_{2} molecules as a function of ν_{i} are presented in Figs 9 and 10 for our different PES models. A clear trend of P (ν) for the model PES-I with ν_{i}, as can be seen in Fig. 9, is the increase of the vibrational excitation with ν_{i}; the maximum P (ν) is located at ν=2,3 and 4 for ν_{i}=0,1 and 2, respectivily. The quasi-classical result reproduces well the observed behavior, particularly for ν_{i}=1 and 2, but misses the high ν tail for ν_{i}=0. The maximum possible vibrational excitation, ν_{max}, is set by Q and is listed for each ν_{i} in Table III. The importance of describing the zero-point motion is once again illustrated in Fig. 9 by the large deviation between the strictly classical result, where P (ν) is strongly peaked for ν around about 0.3, and the quasi-classical P (ν).

It is also interesting to compare the calculated P (ν) with the corresponding result obtained in the collinear model, as done in Fig. 9a. That the collinear model gives much larger vibrational excitation than the flat surface model is not unexpected, given that there is no competing rotational channel to transfer energy to, and that the force exerted by the incident atom on the adsorbed atom is always in the bond direction.

The quantum mechanical vibrational distributions for our different PES models are quite similiar and show no distinct features. For instance, they all decay from about 0.25 for ν above 2 up to ν_{max} = 4. In particular, P (ν) for the models PES-I and -II are almost identical. The quasi-classical results for P (ν) reproduce the gross features of the quantum P (ν) on our different PES models, except once again for the model PES-II. In this latter case the quasiclassical distribution P (ν) is strongly peaked for ν around 0.3, in a very similiar manner as P (j) for the strictly classical result for PES-I. The most notable deviation between the quasi-classical and quantum mechanical results for PES-III and -IV is the underestimation in the low ν region.

D. Rovibrational distributions

A more complete picture of the underlying dynamics behind the rotational and vibrational distributions is revealed in Figs 11 and 12 by the barcharts for the full quantum mechanical rovibrational distributions, P (ν, j). A common behavior is the anticorrelation between rotational and vibrational excitation; high j states are found among the low ν states and low j states among the high ν states. This behaviour is in part a consequence of the kinematical constraint imposed by energy conservation, which prevents, as can be seen in Table III, high j states from being produced with high ν. However, the suppression of states with both low j and ν = 0, as present in all our cases except ν_{i} = 2 for the model PES-I, is a dynamical effect. Another consequence of this kinematical constraint is the effect of ν_{i} on P (ν, j), as shown for the model PES-I in Fig. 11. As ν_{i} increases, the number of states with high ν and low j increases at the expense of states with low ν and high j. The net effect is the decrease in rotational excitation with ν_{i} seen in Fig. 7.

The quantum mechanical P (ν, j) for all our PES models and ν_{i} = 0, except for the model PES-III, are quite similiar, as demonstrated in Figs 11 and 12. This was anticipated by the similiarities of P (ν) and P (j) among these PES models. The characteristic bimodal behavior of these similar P (j), as discussed above, is an effect of the bimodal behavior of P (ν, j) for ν = 0 and 1 with j and the relatively low rotational excitation for ν > 1, which in part is the same kinematical effect discussed above. The bimodal behaviour is most pronunced for ν = 1 and there is only a weak low j contribution to P (ν = 0, j). For the model PES-III, which sticks out of the pack, there is no distinct bimodal behaviour of P (ν, j) for fixed ν. Almost all rovibrational excitations that are kinematically possible are present in P (ν, j) except for the suppresion of P (ν, j) for low j and ν.

1. Off-normal incidence

In our flat-surface approximation, molecules executing helicopter motion, m_{j} ≠ 0, are produced only at off-normal incidence. Furthermore, the conservation of J_{z} in this model implies that the amount of helicopter motion is determined from the reactivity of each m_{j} channel. At typical incident energies, E_{i}, and angles of incidence, θ_{i}, only low |m_{j}| states are produced. This fact is illustrated in Fig. 13 from the relative reaction probabilities, P_{r}(|m_{j}|) at E_{i}=100 meV, θ_{i} = 30^{∘}, and ν_{i} = 0 for the model PES-I, where only m_{j} = 0 and m_{j} = ±1 are appreciably excited. In fact, this behavior can be simply understood from the reactive range of the impact parameters presented in Sec. III A 1. The probability P_{r}(|m_{j}|; b) of finding the incident particle with an impact parameter b in a state with definite |m_{j}| is according to Eq. (2.16) given by g(|m_{j}|)|J_{|m_{j}|}(k_{∥}b)|^{2} where h̄k_{∥} is the lateral momentum of the relative motion of the atoms and g(|m_{j}|) is the degeneray of the state m_{j} (1 for m_{j} = 0 and 2 otherwise). In Sec. III A 1, the reactive range of b was found to be centered around b ≈ 0.7 Å and for the initial conditions E_{i}=100 meV and θ_{i} = 30^{∘} (k_{∥} ≈ 1.7 Å^{−1}), this gives P_{r}(|m_{j}|; b) (plotted as solid diamonds) that are close to the calculated P_{r}(|m_{j}|), as can be seen from Fig. 13.

The rovibrational distributions are found to be sensitive to m_{j}. The rotational excitations tend to increase with |m_{j}| while the vibrational excitation decreases as shown by their average values in Table IV. The partial rotational distributions, P (j; m_{j}) = ∑_{ν} P (ν, j; m_{j}), for both |m_{j}| = 0 and 1 show bimodal behavior, but with the nodes at different j, resulting in a much less pronounced bimodal behavior for P (j), as shown in Fig. 14. This is in contrast to P (j) at normal incidence in Fig. 7a, for the same incident energy normal to the surface. A comparison of the vibrational distributions, P (ν), for normal and off-normal incidence in Figs. 9a and 15 shows that they are quite similiar. The most notable difference is that the off-normal P (ν) peaks at a lower ν than the one at normal incidence.

IV. CONCLUDING REMARKS AND SUMMARY

A flat-surface model has been developed and used to explore the dynamics of the Eley-Rideal recombination of a gas phase H atom with an adsorbed H atom on a metal surface. This model is based on a flat-surface approximation for a single electronically adiabatic potential energy surface, which reduces the six dimensional reactive scattering problem on a rigid surface to a tractable three dimensional problem by introducing three constants of motion; the angular momentum normal to the surface and the lateral components of the linear momentum. The three dimensional reactive scattering problem is treated quantum mechanically by evolving the time dependent Schrödinger equation, using an extension of the pseudospectral method to cylindrical coordinates for the propagation of the wavefunction. In contrast to previous studies based on two dimensional collinear models, this model enables us to calculate reactive cross sections and full rovibrational distributions for a three dimensional rotor. In addition, this flat-surface model should also be applicable to other reactive and nonreactive scattering problems at surfaces.

We have considered four different models for the PES which all are more or less compatible with available knowledge about atomic and molecular hydrogen interactions with Cu(111). We have studied the dependence of the reactive cross section, translational energy distributions and rovibrational state distributions on the kinetic energy and angle of the incident H atom and the vibrational state of the adsorbed H atom. In order to elucidate the role of quantum effects on the dynamics we have compared our results with quasi-classical trajectory calculations. Our main conclusions are:

i) We find that the Eley-Rideal recombination takes place with a cross section ranging up to 0.5 Å^{2}, depending on the PES, and para- and ortho-hydrogen molecules are produced in about equal amounts. The quasi-classical calculations show that the dominant contribution to the cross-section comes from a narrow range of impact parameters centered around a value close to the H_{2} equilibrium separation except for one of the PESs where this range is close to collinear incidence, resulting in a small cross section. We find that the cross section for non-reactive exchange scattering is negligible, ∼ 3 × 10^{−3} Å^{2} while it is large, ∼ 20 Å^{2}, for trapping of the incident H atom via energy and momentum transfer to the adsorbed H atom.

ii) A major part of the large exothermicity of this reaction of ∼2.5eV is channelled into the rovibrational distribution of the molecule, as expected for an attractive PES where a large part of this energy is released in the entrance channel. The rotational and vibrational distributions are broad with quantum numbers up to the maximum energetically possible, j ∼ 20 and ν ∼ 4 − 5, with an anticorrelation between rotational and vibrational excitation; the high j-states are found among the low ν-states. The broad rovibrational energy distribution results in a translational energy distribution which ranges from 0 to ∼ 2eV, with an average of about 0.8-0.9 eV. A characteristic feature among the PESs where the reactive cross section is dominated quasi-classically by non-zero impact parameters is the bimodal rotational distribution. However, the amount of rotational excitation is not simply related to the range of reactive impact parameters.

iii) The most important initial condition effecting the Eley-Rideal reaction dynamics is the vibrational quantum number, ν_{i}, of the adsorbed H atom. The reactive cross section increases non-monotonically with increasing ν_{i}. The efficiacy of energy transfer from the adsorbate vibration to the internal energy of the molecule is about 2. The amount of rotational excitation tends to decrease with increasing ν_{i}, while the amount of vibrational excitation increases. The effect of the incident atom’s kinetic energy on the dynamics is small, as expected for an attractive PES. At off-normal angles of incidence the bimodality of the rotational distribution becomes less pronounced and the rotational distribution consists not only of molecules in cartwheel motion but also in helicopter motion. The amount of helicopter motion is determined simply by the range of reactive impact parameters and the lateral momentum of the relative motion of the two atoms.

iv) The quasi-classical trajectory calculations are able to account for the gross features of the quantum mechanical results, but with some notable exceptions for some of the PESs. Inclusion of the zero-point vibrational motion of the adsorbed atom is found to substantially broaden the widths of the translational, vibrational and rotational distributions, and a strictly classical calculation (no initial zero point energy) is not in general good enough,resulting in too narrow distributions for some of the PESs.

We are presently applying the flat-surface model to different isotopic combinations of the H atoms in order to investigate, as suggested by our previous collinear studies, any interesting mass effects on the ER dynamics. This work will also enable us to make a more direct comparison with available data from molecular beam experiments, and hopefully also stimulate measurements of the internal state distribution of the product molecules. The wealth of dynamical data from molecular beam experiments on the ER reaction between gas phase H atoms and Cl atoms adsorbed on Au(111) makes this system a most interesting case for the application of the flat-surface model, particularly given the observation that the angular distribution of the HCl product seems to be consistent with parallel momentum conservation. There are certainly several limitations of the low-dimensional flat-surface model that require more attention. The importance of the vibrational state of the preadsorbed atom suggests that its lateral vibrational motion might require a better description, by introducing surface corrugation. Also, large changes in the electronic configurations associated with the making and breaking of bonds at metal surface suggests that the assumption of a single electronically adiabatic potential energy surface might have to be relaxed. Ultimately, one would need to include energy dissipation into the substrate excitations in order to be able to treat the role of possible reaction mechanisms intermediate between the idealized ER and LH; i.e., mechanisms involving “hot” adatoms.

ACKNOWLEDGMENTS

This work has been supported by the Swedish Natural Science Research Council (NFR). Bret Jackson acknowledges support from the Division of Chemical Sciences, Office of Basic Energy Sciences, Office of Energy Research, US Department of Energy, under Grant # DE-FG02-87ER13744. Allocation of computer time by the Center for Parallel Computers (PDC) is also gratefully acknowledged.

REFERENCES

[1] R.I. Hall, I. Cadez, M. Landau, F. Pichou, and C. Schermann, Phys. Rev. Lett. 60, 337 (1988).

[2] P.J. Eenshustra, J.H.M. Bonnie, J. Los, and H.J. Hopmann, Phys. Rev. Lett. 60, 341 (1988).

[3] C.T. Rettner, Phys. Rev. Lett. 69, 383 (1992).

[4] K.R. Lykke and B.D. Kay, in Laser Photoionization and Desorption Surface Analysis Techniques, Ed. N.S. Nogar, (SPIE, Bellingham, WA, 1990), vol. 1208, p. 18.

[5] C.T. Rettner and D.J. Auerbach, Science 263, 365 (1994).

[6] E.W. Kuipers, A. Vardi, A. Danon, and A. Amirav, Phys. Rev. Lett. 66, 116 (1991).

[7] C.T. Rettner, J. Chem. Phys. 101, 1529 (1994).

[8] D. D. Koleske, S. M. Gates, and J. A. Schultz, J. Chem. Phys. 99, 5619 (1993).

[9] D. D. Koleske, S. M. Gates, and B. Jackson, J. Chem. Phys. (in press).

[10] P. Kratzer and W. Brenig, Surf. Sci. 254, 275 (1991).

[11] B. Jackson and M. Persson, J. Chem. Phys. 96, 2378 (1992).

[12] B. Jackson and M. Persson, Surf. Sci. 269/270, 195 (1992).

[13] B. Jackson, M. Persson, and B.D. Kay, J. Chem. Phys. 100, 7687 (1994).

[14] C.C. Cheng, S.R. Lucas, H. Gutleben, W.J. Choyke, and J.T. Yates, Jr., J. Am. Chem. Soc. 114, 1249 (1992).

[15] H.A. Michelsen, C.T. Rettner, and D.J. Auerbach, in Surface Reactions, Ed. R.J. Madix, (Springer, Berlin 1993)

[16] B.E. Hayden in Dynamics of Gas-Surface Reactions, Eds C.T. Rettner and M.N.R.

Ashfold (Royal Soc. Chemistry, Cambridge, 1991), p.137.

[17] S. Holloway in Dynamics of Gas-Surface Reactions, Eds C.T. Rettner and M.N.R. Ashfold (Royal Soc. Chemistry, Cambridge, 1991), p.88.

[18] M. Persson and S. Andersson, Surf. Rev. and Lett. 1, 187 (1994).

[19] E.M. McCash, S.F Parker, J. Pritchard, and M.A. Chesters, Surf. Sci. 215, 363 (1989).

[20] K. Gundersen, B. Hammer, K.W. Jacobsen, J.K. Nørskov, J.S. Lin, and V. Milman, Surf. Sci. 285, 27 (1993).

[21] A^{eff} = π(ρ_{c}^{2} + λ^{2}/2 + \sqrt{π/2}ρ_{c}λ) and E_{ρ}^{cut} ≃ h̄^{2}(\sqrt{2π}/(4ρ_{c}λ))/(2μ) < 1meV for ρ_{c} = 3.2Å and λ = 1Å.

[22] M.D. Feit, J.A. Fleck, Jr., and A. Steiger, J. Comp. Phys. 47, 412 (1982).

[23] G.A. Gates, G.R. Darling, and S. Holloway, Chem. Phys. 179, 199 (1994).

[24] D. Neuhauser and M. Baer, J. Chem. Phys. 91, 4651 (1989).

[25] See, e.g. p. 81 in “Computer Simulations of Liquids”, M.P. Allen and D.J. Tildesley, (Clarendon Press, Oxford 1987)

[26] The base of the triangular resolution function have been choosen both for the rotational and vibrational distributions to be equal to 1. Quantities like the reactive cross-section, average center-of-mass energy and angular momentum change less than 1% while the change in the most detailed internal state distribution (angular momentum) is less than 10% when the number of reactive trajectories is decreased by a factor of 2 either by changing the sampling over the impact parameter or the vibrational phase of the adsorbed atom.

[27] The cross sections have been estimated by taking the measured reaction probabilities in Ref. [3] and dividing these values with the estimated saturation coverage of about 0.15Å^{−2}.

[28] J. Harris and B. Kasemo, Surf. Sci. 105, L281, (1981). 37

FIGURES

FIG. 1. Dissociative and pick-up sections of the model PES-I. The upper panel depicts the dissociative section, V (ρ, z = 0, Z), of the model PES-I as a function of ρ and Z. The lower panel depicts the pick-up section, V (ρ, z, Z = z_{0} + m_{2}z/M), of the model PES-I, but now as a function of ρ and z. The contour spacing is 0.2 eV.

FIG. 2. Quantum mechanical (open symbols) and quasi classical (solid symbols) reaction cross sections, σ_{r}, as a function of incident energy, E_{i}, and the initial vibrational state of the adsorbed atom, ν_{i}, at normal incidence and ΔE_{z}= 25 meV. The results are based on the model PES-I. ν_{i} = 0 (circles), 1 (squares), and 2 (diamonds).

FIG. 3. Classical reaction probabilities, P_{r}(b), as a function of impact parameter, b. The upper panel presents results for the model PES-I, for both the adsorbed H-atom at rest (dashed line) and the quasi-classical result (solid line) for ν_{i} = 0. The lower panel depicts the quasi-classical results for ν_{i} = 0 using the potential energy surfaces PES-II (solid line), PES-III (short-dashed line), and PES-IV (long-dashed line). Both panels correspond to E_{i}=75meV and normal incidence.

FIG. 4. Quantum mechanical (open symbols) and quasi-classical (solid symbols) average translational energies, ⟨E_{Z}⟩, as a function of E_{i} and ν_{i} as calculated from the model PES-I at normal incidence and ΔE_{z}= 25 meV. ν_{i}=0 (circles), 1 (squares) and 2 (diamonds).

FIG. 5. Quantum mechanical and (quasi)-classical translational energy distributions, P (E_{Z}), as calculated from PES-I, at E_{i}=75meV, ΔE_{z}= 50 meV and normal incidence. ν_{i}=0 for the quantum mechanical (solid line) and quasi-classical (long-dashed line) P (E_{Z}), while the adsorbed atom is initially at rest for the classical result (short-dashed line).

FIG. 6. Quantum mechanical (solid lines) and quasi-classical (dashed lines) translational energy distributions, P (E_{Z}), as calculated from (a) PES-II, (b) PES-III, and (c) PES-IV, for ν_{i}=0, E_{i}=75meV, ΔE_{z}= 50 meV and normal incidence.

FIG. 7. Quantum mechanical (solid lines) and quasi-classical (long-dashed lines) rotational distributions, P (j), as calculated from PES-I, for E_{i}=75meV, ΔE_{z}= 50 meV and normal incidence.

(a) ν_{i}=0 and strictly classical result (short-dashed line), (b) ν_{i}=1, and (c) ν_{i}=2.

FIG. 8. Quantum mechanical (solid lines) and quasii-classical (dashed lines) rotational distributions, P (j), as calculated from (a) PES-II, (b) PES-III, and (c) PES-IV, for ν_{i}=0, E_{i}=75meV, and normal incidence.

FIG. 9. Quantum mechanical (open symbols) and quasi-classical (solid symbols) vibrational distributions, P (ν), as calculated from PES-I for E_{i}=75meV, ΔE_{z}= 50 meV and normal incidence. (a) ν_{i}=0, strictly classical result (short-dashed line) and two-dimensional collinear result (diamonds), (b) ν_{i}=1, and (c) ν_{i}=2.

FIG. 10. Quantum mechanical (open symbols) and quasi-classical (solid symbols) vibrational distributions, P (ν), as calculated from (a) PES-II, (b) PES-III, and (c) PES-IV, for ν_{i}=0, E_{i}=75meV, ΔE_{z}= 50 meV , and normal incidence.

FIG. 11. Quantum mechanical rovibrational distributions, P (ν, j), as calculated from PES-I for (a) ν_{i}=0, (b) ν_{i}=1, and (c) ν_{i}=2, for E_{i}=75meV, ΔE_{z}= 50 meV and normal incidence.

FIG. 12. Quantum mechanical rovibrational distributions, P (ν, j), as calculated from (a) PES-II, (b) PES-III, and (c) PES-IV, for ν_{i}=0, E_{i}=75meV, ΔE_{z}= 50 meV and normal incidence.

FIG. 13. Quantum mechanical reaction probability as a function of azimuthal quantum number, P (|m_{j}|), at off-normal incidence for the model PES-I. E_{i}=0.1 eV, θ_{i} = 30^{∘} and ν_{i} = 0. The results (solid diamonds) from a simple model discussed in Sec. III D 1 are also included.

FIG. 14. Quantum mechanical rotational distribution, P (j), at off-normal incidence for the model PES-I. E_{i}=0.1 eV, ΔE_{z}= 50 meV, θ_{i} = 30^{∘} and ν_{i} = 0.

FIG. 15. Quantum mechanical vibrational distribution, P (ν), at off-normal incidence for the model PES-I. E_{i}=0.1 eV, ΔE_{z}= 50 meV, θ_{i} = 30^{∘} and ν_{i} = 0.

TABLES

TABLE I. Parameters for the four different model potential energy surfaces. D_{i}, α_{i} and x_{i} are V_{i} parameters for the H_{2} molecule (i = m) and H-Metal (i = a) interactions. Δ_{m} and Δ_{a} are Sato parameters for the different model PESs. V^{∗} and (ρ^{∗}, Z^{∗}) are the corresponding heights and locations of the activation barriers of the dissociative section. V_{exo} gives a measure of the potential energy release in the entrance channel of the different model PESs.

H_{2}  H-Metal
D_{i} (eV)  4.745  2.45
α_{i} (Å^{−1})  1.943  1.0
x_{i} (Å)  0.741  0.0
PES-I  PES-II  PES-III  PES-IV
Δ_{m}  -0.2  -0.2  -0.1  -0.3
Δ_{a}  0.2  0.1  0.0  0.1
V^{∗} (eV)  0.59  0.75  0.68  1.05
(ρ^{∗}, Z^{∗}) (Å)  (1.16,0.12)  (1.08,0.2)  (0.96,0.34)  (1.16,0.2)
V_{exo} (eV)  1.7  1.3  1.1  1.1

TABLE II. Quantum mechanical and quasi-classical reaction cross sections σ_{r}, average translational, and internal energies, ⟨E_{int}⟩ and ⟨E_{Z}⟩, respectively. They are calculated at normal incidence with the incident energy E_{i} = 0.075 meV and initial vibrational state ν_{i} = 0 for the four different potential energy surfaces. The quasi-classical results are within the parenthesis.

PES-I  PES-II  PES-III  PES-IV
σ_{r}^{(p)}/σ_{r}^{(o)}  0.95  0.94  0.94  0.90
σ_{r}(Å^{2})  0.47(0.42)  0.11(0.054)  0.024(0.003)  0.063(0.049)
⟨E_{Z}⟩(eV)  0.88(1.00)  0.85(0.59)  0.82(0.59)  0.83(0.78)
⟨E_{int}⟩(eV)  1.56(1.43)  1.59(1.85)  1.62(1.85)  1.61(1.66)

TABLE III. Exothermicity and maximum possible rovibrational excitations. Q is the energy released in the reaction (see Eq. 3.1) as a function of the initial vibrational state, ν_{i}, of the adsorbed H atom for E_{i} = 75meV, and j_{max} and ν_{max} are the corresponding maximum rotational and vibrational quantum numbers allowed energetically by this Q. The values within the parenthesis are the corresponding quasi-classical results. The lower panel gives the maximum rotational excitation, j_{max}(ν), allowed energetically for each individual vibrational state ν.

ν_{i}  Q(eV)  j_{max}  ν_{max}
0  2.44  18  (20.8)  4  (4.8)
1  2.57  19  (21.5)  5  (5.1)
2  2.70  20  (22.2)  5  (5.5)
j_{max}(ν)
ν_{i}  ν=0  1  2  3  4  5
0  18  16  13  10  6  –
1  19  17  14  11  8  2
2  20  18  15  12  9  5

TABLE IV. Average internal state distributions at off-normal incidence. Average vibrational and rotational quantum numbers, < ν > and < j >, respectively, and internal energy, < E_{int} >, of the product H_{2} molecule as a function of azimuthal quantum number |m_{j}| at an angle of incidence, θ_{i} = 30^{∘}, and incident energy, E_{i}= 0.1 eV.

|m_{j}|=  0  1  2  3
< ν >  1.4  1.1  1.0  0.5
< j >  8.5  9.6  10.1  12.8
< E_{int} >(eV)  1.60  1.57  1.59  1.75