Kinetics of NH3 Desorption and Diffusion on Pt: Implications for the Ostwald Process

We report accurate time-resolved measurements of NH3 desorption from Pt(111) and Pt(332) and use these results to determine elementary rate constants for desorption from steps, from (111) terrace sites and for diffusion on (111) terraces. Modeling the extracted rate constants with transition state theory, we find that conventional models for partition functions, which rely on uncoupled degrees of freedom (DOFs), are not able to reproduce the experimental observations. The results can be reproduced using a more sophisticated partition function, which couples DOFs that are most sensitive to NH3 translation parallel to the surface; this approach yields accurate values for the NH3 binding energy to Pt(111) (1.13 ± 0.02 eV) and the diffusion barrier (0.71 ± 0.04 eV). In addition, we determine NH3’s binding energy preference for steps over terraces on Pt (0.23 ± 0.03 eV). The ratio of the diffusion barrier to desorption energy is ∼0.65, in violation of the so-called 12% rule. Using our derived diffusion/desorption rates, we explain why established rate models of the Ostwald process incorrectly predict low selectivity and yields of NO under typical reactor operating conditions. Our results suggest that mean-field kinetics models have limited applicability for modeling the Ostwald process.


S2. Determination of the thermal sticking coefficient of NH3 at Pt(111)
The thermal sticking coefficient on Pt(111) is derived from the experimentally obtained translational energy distributions of desorbing NH3 molecules using the principle of detailed balance. From detailed balance, the relation between the desorbing flux distribution ( , S ) , the Maxwell-Boltzmann thermal (3D) distribution MB ( , S ) , and the initial sticking probability function 0 ( ) is: ( , S ) ∝ 0 ( ) MB ( , S ), (S1) We observe that the desorbing NH3 is subthermal, which indicates a decreasing sticking coefficient with increasing translational energy. The experimentally derived translational energy distributions of NH3 desorbing from Pt(111) at 473-873 K are fitted globally to an empirical temperature independent 0 ( ). After fitting, the resulting sticking probability is obtained.  This approach gives the shape but not the absolute value of the sticking probability curve. We establish the absolute value using 0 ( = 0) = 1, an excellent assumption for a non-activated adsorption system like this one 1 . Figure S3 Figure S4: The kinetic traces of NH3 desorption are shown for two Pt(111) crystals with different step densities: 0.15 ± 0.05 % (+) and 0.4 ± 0.2 % (×). At reduced step density the trace resembles more closely a single exponential. This shows that the slow component of the bi-exponential kinetic trace arises due to NH3 interactions with steps.

S4. The Diffusion-Desorption
We recognize that, within our approximations, expression S7 is nothing else but the TST expression for direct desorption from steps. We notice that this approximation is only valid for our setup of the kinetic model as we assume that the TS for terrace hopping is the same as for hopping between terrace and step sites.

S4.2. Kinetic rate equations in a non-uniform spatial grid
Our approach follows closely the work reported in Ref. 3 , where it was shown how to account for surface diffusion by sequential hopping events between adjacent binding sites. The method employs periodic boundary conditions. For a step density of 0.4% this requires solving 250 coupled differential equations numerically, which is time consuming especially when employed within an optimization routine designed to obtain fits to the experimental data. However, in order to describe diffusion accurately, we do not need to treat each terrace site individually. The only unique binding sites that must be treated individually are the steps-the individual binding sites at the terraces can be coarse-grained into spatial regions that include multiple atomic rows, where the diffusion equation is solved. This dramatically reduces the number of coupled differential equations that need to be solved. We use spatial grids with atomic rows close to the steps, which increase in size when approaching the center of the terrace.
The diffusion-desorption rate equations are as follows: (S12) NH3's fractional occupation of the binding site is given by . Index =0 indicates the step site, while indices 1 and − 1 indicates the terrace sites adjacent to the step. Indices 2 until − 2 are terrace sites. In the first term is the flux of molecules impinging on the surface, which is described by the molecular beam dosing function. For terrace sites, we use the sticking probabilities ( 0 T ) from Fig. S3 and for step sites, we set the sticking probability ( 0 S ) to unity. This is justified as the measured speed distributions of NH3 desorbing from Pt(332) is thermal with tr = S and exhibits a cos( ) angular distribution. In the equations above, h T/S are the hopping rate constants between adjacent atomic rows and is the diffusion constant for NH3 at terraces. The terms in parentheses are needed to account for effect of binding site occupation.
The terrace sites further away from steps are treated in non-uniform spatial elements with the width . The definition of the spatial elements is shown in Fig. S5(a). In these regions, the terrace diffusion constant is used. The resulting diffusion equations are in agreement with previously reported schemes 4 . We compare the desorption rates emerging from our nonuniform grid model to the model with only atomic rows, as suggested in Ref. 3 , using a broad range of rate parameters and choose appropriate grids to reach a relative accuracy of at least 0.04% when modelling the desorption rate. This allows especially for low step density crystals, to reduce the evaluation time to a minimum, while having accurate results suitable for the kinetic model to fit the experimental data. The diffusion constant can be derived from the hopping rate constant between adjacent atomic rows following Ref. 5 or from the well-known formula for mean-square displacement. However, we have to account for the fact that hopping between adjacent rows on a Pt(111) will have degenerate pathways, as sketched in Fig. S5(b). Thus, is given by: where ̃ is the distance between two adjacent atomic rows (see Figure S5). When the width of spatial elements, , is expressed in units of atomic rows as: =̃.
(S14) Then Eq. S11 is reduced to: (S15) The differential equations are solved following the protocol described in Ref. 6 and the corresponding desorption flux is obtained from: (S16)

S4.3. Fitting procedure and uncertainty estimation
The global optimization routine is based on least-square minimization of the residual, error-weighted difference between experiment and modelsee Ref. 1 for details. After the initial guess of the six globally optimized rate parameters the global fitting employs a cyclical two step procedure.
1. First, the desorption rate is numerically calculated (Eq. S16) based on the initial guess or the optimized rate parameter values established in the previous execution of step 2. Then, the amplitudes of the direct scattering and the trapping-desorption channel are fitted (locally) to each individual kinetic trace. For description of direct scattering we use a function proportional to the incident molecular beam B ( ) in Eq. S8 and Eq. S9. The residual over the total dataset is calculated. 2. The local fitting procedure of the direct-scattering and trapping desorption amplitudes ensures the lowest possible residual or the present set of rate parameters. The rate parameters are modified following the least squares procedure and step 1 is repeated until convergence of the residual is reached. The uncertainties of the fitted parameters were estimated by displacing the rate parameter of interest from its best fit value in small increments, fixing it and re-optimizing the remaining five parameter globally until a noticeable mismatch (deviations beyond the noise level of the data) between the re-fitted model and the data was observed. The parameter values between which an accurate fit could be obtained define the uncertainty range (see arrow with error-bar of Figure 3(b) and 3(c)). From the present analysis, we find that the strongest parameter correlation is the correlation between and a of the individual rate constants (increasing prefactor leads to increasing activation energy). We also find a strong correlation between h T and h S (which have an almost constant ratio). Correlations for other parameters are less pronounced. For terrace desorption (see Figure 3(b) and (c)) the error analysis compares well with a complementary error estimation (see Section 3.2 of the main text) reflecting its reasonability. (S17)

S5. Coverage dependence of NH3 binding energies at
Here the vibrational partition function of individual DOFs appear in an overall product form. NH3*'s possesses six internal modes-1−5 and umb that correlate to the normal modes of the gas-phase molecule. Due to adsorption to the surface, umb no longer inverts and two of rotational DOFs with rotation axes perpendicular to NH3*'s symmetry axis, become frustrated rotations-and -and translation perpendicular to the surface--correlates to the external stretch of the adsorbate. All of these modes are treated as simple quantum harmonic oscillators, where, is the harmonic frequency of the mode at the most stable binding site (see Table 2 of the main text). The rotation around NH3* symmetry axis remains free, validated from previous work 8 , and is described by the 1Drotational partition function: where is the rotational constant and is the corresponding symmetry number.
(S20) ( ; , ) and ( ; , ) contain adjustable parameters, and , that are used to reproduce the harmonic frequency of and obtained from DFT, and the experimentally derived energy barrier for site-to-site hopping. To match both of these quantities, the commonly employed 1D model potential functions 9-10, 12 are better replaced by functions of the following form: where, is or , is the length of the unit cell (2.77 Å), 0,h T is the experimentally derived hopping barrier (0.71 eV). The hindered translation potential used in uTST ( ; = 2, = 0.47) is shown in Figure S7. The configuration integral above is classical and converges to 0 at low temperatures. We ensure proper quantum behavior of the partition function-i.e. ≥ 1-at low temperature by including the term qHO cHO ⁄ . This correction term is described in Ref. 10 in detail.

S6.2. NH3* -coupled for use in the coupled TST model
Despite energy sampling of a multidimensional potential energy surface for NH3 at Pt(111) would be the most straightforward approach it is computationally demanding. Thus, we require clever approximative methods for the description of the potential energy going beyond the harmonic approximation and the assumption of normal modes.
Because the coupling of DOFs for partition functions in TST is not common and there are many ways this might be implemented, we first illustrate our idea using a 2D toy system. Consider an atom on a one-dimensional catalyst surface, which exhibits displacement over a diffusion barrier along and harmonic vibration normal to the surface along . For motion along , the harmonic frequency in the z direction varies due to the changing interaction with the catalyst. This can be expressed as an -dependent force constant ( ), where the total energy expression is then approximated by: Here, the first two terms describe the kinetic energy along and , the third term describes the potential energy along the minimum energy pathway for diffusion, while the fourth term describes the potential energy for displacement. This formulation explicitly introduces coupling of and as depends on . The classical partition function for the toy system then becomes: where is the length of the unit cell. Combining Eq. S22 and S23 leads to: Here, ( ) is the harmonic stretch frequency and the coupling between and manifests itself by a modification of the configurational integral over , which now includes the -coordinate stretch frequency in the -dependent vibrational partition function cHO ( ). Since ≪ ℎ , we replace cHO ( ) by its quantum mechanical counterpart: While this approach is certainly ad hoc, it represents a useful step beyond uTST, that can be taken without significantly increasing the effort beyond the established methods and with a more realistic accounting of the adsorbate entropy.
In the NH3 Pt system, we couple four DOFsumb , , , -to hindered translation parallel to the Pt surface. The frequencies of these DOFs change most when moving along the hindered translation coordinate from the on-top binding site to the hopping transition state (see Table 2). Then, in analogy to the toy system described above, we write ad cTST in the following form: (S26) The first six terms before the integral are identical to those used in uTST and represent vibrational partition functions that are independent of and . As a consequence of coupling, the minimum energy pathway emerging from a linear combination of multiple normal modes, cannot be described by hindered translational frequencies of the initial state, as we have done in uTST approach. The in-plane interaction potential, , ( , ) is parametrized similarly as described for uTST with the difference that here ( ; , ) are fitted to the shape of the minimum energy path for hopping obtained from CI-NEB. See Figure S7, which shows the DFT minimum energy path as crosses and the best fit ( ; = 4.5, = 0.64) as a black solid line. The configuration integral now includes -and -dependent vibrational partition functions for the four DOFs whose frequencies change most dramatically when NH3 translates from the on-top site to the diffusion transition state, according to DFT calculations.
(S28) ( , ) now describes the softening of the modes as the NH3 molecule moves away from the atop binding site. Note, that we only use DFT frequency information for two structures: i,top for the top site and ,bri for the bridge site, which is the TS for hopping. The scaling with the in-plane coordinates is done assuming radial symmetry around the top site, which is a sensible approximation due to a localized bond formation between NH3 and the top site of Pt. We have optimized such that ( , ) at the bridge site does not deviate by more than 0.1% from bri and yields good agreement with the temperature dependence of the desorption rate from Pt(111). The ( , ) functions are shown as blue solid lines in Fig. S7. The uncertainties of NH3* binding energy, that were derived from the cTST fit to the data, emerge from uncertainties associated with . We have investigated the sensitivity of derived binding energies to values that appeared reasonable in shape ( ( , ) with FWHM between 0.5 and 2.0 Å) and estimated the uncertainty based on this.

S6.3. Transition state for desorption
The desorbed gas-phase molecule is taken as the transition state for desorption in this entire work. Its partition function is given by Here, 2D tr is the partition function for free translation of a 2D ideal gas: with being the mass of NH3. is the area of the Pt(111) unit cell, in which the partition function is defined. Note that translation normal to the surface is associated with the reaction coordinate in TST and is not part of ‡ .
3D rot is the classical rotational partition function of a symmetric top: where and are the rotational constants and is the associated symmetry number. The vibrational partition function can be accurately approximated by a simple quantum harmonic oscillator partition function qHO for five DOFs: the symmetric and two asymmetric stretches as well as two asymmetric bends, conventionally referred to as 1 − 5 . qHO takes the usual form: .

(S32)
Here, is the harmonic frequency of mode . The umbrella mode is treated by explicit counting of the first ten vibrational levels, umb , to account for its complex quantum structure.
All parameters needed to evaluate the partition function can be found in Refs. [13][14][15] .

S6.4. Transition state for hopping
The transition state for hopping is described with the following partition function: The first nine partition functions describes the internal NH3* vibrations, the frustrated rotations and translation (along the surface normal) as quantum harmonic oscillator with TS frequencies from Table 2. NH3*'s rotation around its symmetry axis is assumed to be free. Since we find a rather low hindered translational frequency along the -coordinate ( is taken as the reaction coordinate) we treat it by a hindered translator. Based on DFT generated harmonic frequencies, we estimate the shape of ℎ( ) using the commonly applied model potential energy function from Ref. 9 . We also include the low temperature correction to the partition function, previously described in Section 8.6.1.
We emphasize that this partition function does not take coupling of different DOFs into account. Including coupling would require vibrational frequencies of NH3 along while would be fixed at the TS. Including coupling corrections would result in a somewhat higher value of the partition function than offered by the present approach.

S7. Extraction of First-Order rate constants from the desorption rates
We determine the first-order rate constants for NH3 desorption from Pt(332) and the fast component in Pt(111) by fitting the flux ( ) vs. residence time using a function with two contributions, ( ) = × DS( ) + × TD( , des ), comprising direct scattering (DS) and a trapping-desorption (TD) with des as the desorption constant. The DS contribution has the temporal shape of the incident molecular beam, while the TD contribution is an exponential decay convoluted with the molecular beam temporal profile. For Pt(332) analysis, the full trace is fitted with Eq. S35, while for Pt(111) we cut the slow part off for the fit. The corresponding fits to the kinetic traces are shown in Fig. S8, while the fitted rate constants are shown in Fig. 3(a) (as red circles and black crosses) of the main text.

S8. Terrace-assisted desorption model from stepped surfaces
Here we derive the formula for the terrace-assisted desorption (TAD) model given in Eq. 4 of the main text. We assume a steady state between NH3 at steps and terraces including the possibility to desorb from terraces and steps. An effective first-order desorption rate constant for NH3 desorption from a stepped surface can be derived. We assume the following processes: where reflects the step density of the crystal. We define the total desorption rate and the effective first-order desorption rate constant as follows: Since NH3 has a big energy preference for steps, after a short induction period, the total desorption rate will be characteristic of the concentration pool at steps with the corresponding effective desorption rate constant. Assuming steady-state for [NH 3,T ], which is justified as shown in Section 3.5.2 of the main text we obtain: Combining Eq. S39 and S40 we obtain eff , shown as Eq. 4 of the main text and below: We emphasize that the present model is kept simple and thus works only at stepped surfaces and adsorbates with strong energy preference. The results of the TAD model are shown in Fig. S9 as a black solid line in excellent agreement with experiment rate constants ( 2 = 2.9) obtained for NH3 desorption from Pt (332) (black crosses). We also compare the results of an earlier model from Serri, Tully and Cardillo 3 (STC model) where no direct step desorption was considered (black dashed line). Neglecting direct step desorption degrades the fit ( 2 = 16.3) somewhat underestimating the experimental rate constants. Figure S9: Comparison of experimental NH3 desorption rate constants obtained with Pt (332) (black crosses) and three models for desorption from stepped surfaces. 2 in the legend is the reduced chi-square value.

S9. Coverage dependent desorption rates and adsorption isotherms
While the elementary reaction rate constants reported in this work were obtained for the low coverage limit, they can be used as a starting point to extend our understanding of the kinetics of NH3 desorption from Pt. Here we explain how we have included the coverage dependence of the desorption rate constant from terraces and steps, based on the arguments elaborated in Section. 8.4.1 of the SI.
The coverage influences both the prefactor and the binding energy of the desorption rate constant: From DFT calculations we find that NH3 binding energy decreases linearly from = 0.06 to 0.25 ML with a slope of = −1.61 eV/ML. We use the scaling of binding energy with coverage from DFT calculations, while using the experimentally derived zero-coverage binding energy ( 0,d T = 1.13 ± 0.02 eV). The resulting coverage dependent NH3 binding energy at terraces is then the following 0 ( ) / eV = 1.13 − 1.61 × ( ML ).
We emphasize that the NH3 saturation coverage of 0.25 ML is defined with respect to the Pt(111) facet. We include the coverage dependence of the prefactor, ( ), using the calculated harmonic frequencies of the adsorbate at 0.06 and 0.25 ML. We assume that the logarithm of the prefactor, which reflects the entropy difference between the initial state and the transition state, scales linearly with coverage and is described by: log 10 ( ( , )) = (1 − 4 ) log 10 ( 0.06 ( )) + 4 log 10 ( 0.25 ( )).
We note in passing that the prefactor calculated with frequencies at 0.06 ML is assumed to represent the prefactor in the zero-coverage limit. This introduces negligible changes in the results because the prefactor has only weak sensitivity to the assumptions of its scaling and the dominant contribution to the desorption rate comes from coverage dependence of NH3's binding energy.
Using the coverage dependent desorption rate constant we simulate TPD spectra from previous work with the same experimental conditions, results of this simulation are shown in Fig. S10. We find overall good agreement with TPD spectra of the monolayer peak obtained from earlier studies.
Figure S10: Previously measured NH3 TPDs from Pt(111) from Ref. 16 (black crosses). Ref. 16 does not explicitly indicate the initial coverage in each of the three TPD spectra. The low temperature feature below 200 K is due to desorption of the multilayer only seen under conditions of the highest initial ammonia dose. At lower dosing, two TPD spectra are reported that represent desorption of the NH3 monolayer at different initial coverages. The shoulder seen at 450 K emerges most likely from desorption influenced by steps. The red lines show TPD simulations of NH3 desorption from Pt(111) based on Eq. S41 at three initial NH3 coverages (see legend). It is clear that at shape of the TPD is strongly dependent on initial coverage-this reflects the influence of adsorbate-adsorbate interactions. The simulation captures this behavior reasonably well, validating our approach.
The adsorption isotherms are simulated using the diffusion-desorption kinetic model described in Section 8.4.2 of the SI. The stationary adsorption flux (at each spatial element ) is described by the Hertz-Knudsen equation: instead of the molecular beam dosing function. For the chosen range of NH 3 and the rate equations are propagated in time until the stationary-state is established, from which the coverages at steps and terraces is derived. For simulation of the KB model 17 we use a surface with only terrace sites (as assumed by the KB model) in which the coverage independent NH3 desorption rate constant is described by: The simulation results are shown in Figure 6 of the main text.