ARTICLE pubs.acs.org/JPCC
Theoretical Study of the Dynamics of Collisions Between HCl and ω-Hydroxylated Alkanethiol Self-Assembled Monolayers William A. Alexander and Diego Troya* Department of Chemistry, Virginia Tech, Blacksburg, Virginia 24061-0212, United States
bS Supporting Information ABSTRACT: We present a classical-trajectory study of collisions of HCl with hydroxylated alkanethiol self-assembled monolayers. The potential-energy surface used in the trajectory propagation is a combination of the standard OPLS force field to describe the surface and an analytical potential for the gas/surface interaction developed in this work. The gas/surface potential has been derived based on high-quality electronic-structure calculations of model HCl-alcohol systems in the gas phase and includes a flexible Buckingham term and a Coulombic term. The results of the trajectories calculations are in good agreement with recent molecular-beam experiments on the same system, thereby lending support to the accuracy of the calculations. The collision dynamics differ vastly from prior scattering studies involving rare gases and CO, primarily because the gas/surface attraction governed by hydrogen bonding dramatically increases the ability of the gas molecules to trap on the surface for extended times. The properties of the desorbing HCl molecules are largely insensitive to the initial collision energy and are only mildly affected by the incident angle. An analysis of the reaction mechanism reveals the distinct dynamics of trajectories that either recoil from the surface directly or undergo multiple collisions with the surface and result in thermalization.
’ INTRODUCTION Studying energy transfer in collisions of gases with surfaces yields information important to a fundamental understanding of gas/surface chemical dynamics, which is relevant in a variety of fields. For instance, the accommodation and transport of gases on hydrocarbon surfaces plays an important role in atmospheric chemistry. Surfactant-coated aerosols are abundant in the troposphere, and these coatings represent a barrier for gaseous species to enter and react with species interior to the droplet.1 The oxidation of these organic surfaces by various reactive pathways typically confers hydrophilic character to the surface, which may aid in dissolution of polar atmospheric molecules into the aerosol. Once in the aerosol, the gas species may undergo subsequent reaction. While HCl is a sink species for catalytic, ozone-degrading Cl radicals, HCl dissolved into aerosol particles, or other particles such as water droplets and ice crystals, may react to recover Cl2, which upon photolysis restarts the ozone-destroying catalytic process. Even before dissolution into the particle, hydroxyl sites on the surfaces can bind the gas-phase HCl through hydrogen-bonding interactions, initiate reaction by dissociating HCl, or serve as a surface-bound reagent for other reactions.2 Each of these scenarios is mediated by the initial gas/ surface collision, where the interactions of HCl with the organic surface determine whether the molecule will scatter away directly or accommodate on the surface. The dynamics of this initial collision between HCl and organic surfaces has recently been studied using time-of-flight molecular-beam r 2011 American Chemical Society
scattering techniques by Lohr et al.3,4 Hydrophobic and hydrophilic surfaces were modeled with straight-chain CH3-(CH2)15-SH and terminally hydroxylated OH-(CH2)16-SH alkanethiol self-assembled monolayer surfaces on gold, respectively. In those experimental studies, it was found that Ar and HCl had similar accommodation fractions and energy-transfer dynamics when colliding with the nonpolar hydrocarbon surface. However, these gases showed remarkably different dynamics when scattering from the hydroxylated surface. Using the traditional separation of the products’ energy distribution into an impulsive scattering (IS) and thermal desorption (TD) channel,5 the experiments indicate that while 43% of the Ar atoms recoil from the OH-SAM with a thermal energy distribution at the surface temperature when colliding at 85 kJ 3 mol-1 and 30, 73% of the HCl thermalize at the same conditions, showing the effect of gas/surface hydrogen bonding on the dynamics of interfacial scattering. The results that Ar and HCl exhibit similar scattering dynamics in collisions with CH3-SAMs imply that kinematics, rather than the weak gas/surface chemical forces in that system, determine the scattering outcome, such that gases with similar size, mass, and incident energy follow similar scattering pathways regardless of their electronic structure. This result is consistent with recent studies of Ar and DCl scattering from squalane (a liquid alkane at room Received: November 7, 2010 Revised: December 11, 2010 Published: January 18, 2011 2273
dx.doi.org/10.1021/jp1106422 | J. Phys. Chem. C 2011, 115, 2273–2283
The Journal of Physical Chemistry C temperature)6 and with earlier experimental work on the scattering of Ne, H2O, NH3, and CH4, also from squalane.7 The trend that molecules with similar mass display similar behavior when scattering from the nonpolar surfaces has also been corroborated by studies using alkanethiol self-assembled monolayers.8 Another revealing result of Lohr and co-workers was that the OH-SAM is a more rigid collision partner than a regular CH3-SAM surface due to the extensive hydrogen-bonding network formed by the terminal OH groups. The larger rigidity of the OH-SAM can be appreciated by the percentage of collisions in which Ar atoms recoil directly from the surface (57%), which is notably larger than for a less rigid regular CH3-SAM (39%). The extensive hydrogen-bonding network in the OH-SAM has been recently thoroughly studied via molecular-dynamics simulations by the Hase group.9 The simulations reveal that, in comparison to a normal alkanethiol SAM, the OH-SAM is disordered, with extensive, dynamic hydrogen bonding between the OH groups, giving rise to local areas on the surface where clusters and chains of hydrogen-bonded hydroxyl groups form a rigid “glassy” surface. To form these networks, the hydroxyl groups can approximate to an O-O distance of ∼3 Å, in comparison to the 5 Å nearest-neighbor spacing of the sulfur anchor groups. This clustering results in small, unordered domains in which the underlying hydrocarbon regions of the alkyl chain are exposed. Trajectory calculations of Ar scattering from the CH3- and OH-SAM surfaces agree with the experimental findings that Ar recoils from the OH-SAM with higher average final energies than from the CH3-SAM and that penetration into the monolayer is decreased on the hydroxyl surface.9 Although ab initio calculations of the potential-energy surface show the Ar/OH-SAM potential to be more attractive than the Ar/CH3-SAM potential, the rigidity of the hydrogenbonding network results in Ar scattering with slightly higher final energies in direct scattering events. A number of other dynamics investigations of HCl-surface collisions involving hydroxylated systems have been performed recently using both experiment and theory. The Nathanson group provided extensive insight into the scattering, trapping, and dissociation channels when HCl collides with sulfuric acid10,11 and glycerol.6,12-14 The liquid scattering experiments reveal high trapping probabilities, in qualitative agreement with the work of Lohr et al. on SAM surfaces. Pettersson and coworkers performed molecular-dynamics and molecular-beam studies of HCl and Ar in collision with water ice surfaces15-18 and demonstrated that these species transfer energy efficiently to the surface. Additionally, they observed distinct scattering channels for the recoiling HCl molecule: direct inelastic scattering, trapping followed by prompt desorption, and long-time uptake.18 The first channel is analogous to the aforementioned impulsive scattering (IS) channel observed in rare-gas scattering, while the second and third are due to molecules which thermalize on the surface, with the third being characterized by the formation of realtively strong hydrogen bonds between the gas and the surface. The studies of Lohr and co-workers of HCl scattering from OH-SAM surfaces also showed evidence for this third, long-lived trapping channel.3 With this study, we aim to investigate the mechanistic details of HCl’s initial collision with a model hydroxylated organic surface via use of classical trajectories. We propagate the trajectories using our own potential-energy surface (PES) derived from high-accuracy ab initio calculations, the implementation of which is extensively detailed within. Our general approach is
ARTICLE
analogous to work we and others have done previously in studies of rare-gas9,19-25 and small-molecule26-30 collisions with model organic surfaces. In all of those studies, the organic surface was designed to mimic the self-assembled monolayers (SAMs) used extensively in many experimental scattering studies.9,21,31-38 However, the increased complexity of describing the hydrogen-bonding attractions between the gas and the surface makes the HCl/ OH-SAM system substantially more challenging to model than the previous work involving straight-chain alkanethiol surfaces and rare gases. An accurate description of a hydrogen-bonding system such as HCl/OH-SAM requires that the long-range character of the gas/surface interactions must be captured within the PES. This problem is not as apparent in simpler systems, such as rare gases scattering from alkanethiol SAMs, where short-range repulsions and kinematic effects dominate the scattering outcome19,21 and the attractive dispersion interactions decrease rapidly with the separation between the gas and the surface. The remainder of the paper is divided into two parts. In the first, we describe our methods for deriving an analytical two-body pairwise potential for the HCl/OH-SAM system derived by fitting high-accuracy electronic-structure calculations to model gas-phase systems. In the second section, we enlist our gas/ surface model to perform classical trajectory simulations of HCl scattering from hydroxyl-terminated self-assembled monolayers. We validate our model with comparison to existing molecularbeam scattering experiments, subsequently investigate the influence of collision energy and incidence angle on the scattering dynamics, and provide details of the microscopic collision mechanism.
’ COMPUTATIONAL DETAILS Potential-Energy Surfaces. The potential-energy surfaces employed to evolve the HCl/SAM trajectories are analgous to those described in detail in our prior work on gas/SAM collisions.19-21,26,27 Briefly, we divide the global potential into three terms: the potential describing the organic monolayer (surface potential), the potential describing the HCl molecule (gas potential), and the potential for the HCl/SAM interactions (gas/surface potential). We employ an explicit-atom model using the OPLSAA force field39 for the surface potential, as this standard force field bears out the experimental structure of regular SAMs,40 with some modifications, as listed in Table S1 of the Supporting Information. With this parameter set, the OH-SAM behaves similarly to the behavior described by Hase and co-workers,9 with localized hydrogen-bonding clusters and chains which break and reform over time. While the nature of the hydrogen-bonding network is dynamic, the overall extent of hydrogen bonding is consistent throughout the simulation. To form hydrogen bonds, the OH groups move closer together (∼3 Å O-O distance) than the normal SAM nearest neighbor spacing of 5 Å for the sulfur head groups. As such, the normal ordered and periodic nature of the SAM surface is disrupted, creating domains of glassy hydrogenbonded islands alternating with crevices in which the interior hydrocarbon backbone is exposed. These two different domain types represent quite different collision partners and should give rise to substantially different scattering dynamics when HCl impinges on each. In the experiment these effects will be averaged, but in the current study, we can differentiate between these two impact domains to gain further insight into the collision dynamics. 2274
dx.doi.org/10.1021/jp1106422 |J. Phys. Chem. C 2011, 115, 2273–2283
The Journal of Physical Chemistry C
ARTICLE
Figure 1. Schematic of the HCl-CH3OH orientations investigated with ab initio calculations in this work.
A standard Morse function is used to describe the gas potential h i2 ð1Þ VHCl ¼ D 1 - e - βðr - re Þ where D, β, and re are the bond dissociation energy, curvature of the potential, and equilibrium bond length, respectively. To appropriately describe the HCl molecule, D and re are taken to be 103.15 kJ 3 mol-1 41 and 1.2750 Å.42 The value of β = 1.9005 Å-1 was determined by fitting to established experimental values of HCl’s spectroscopic constants.42 The resulting fundamental frequency and rotational constant obtained with our Morse potential are ωe = 2990.9 cm-1 and Be = 10.582 cm-1, which are in excellent agreement with the experimental values of ωe = 2990.9 cm-1 and Be = 10.596 cm-1.42 Regarding the gas/surface potential, earlier work demonstrated that popular force fields commonly fail to model accurately intermolecular interactions involving gases with SAMs.19 Instead, analytic potentials must be derived from high-quality ab initio calculations26,43 for adequate accuracy, and we followed such strategy in this work for the HCl/OH-SAM system. In the following, we first describe the ab initio calculations used to capture the intermolecular interactions of HCl with OH-SAMs and then detail how we use those electronic-structure calculations to develop pairwise analytical potentials that can be conveniently used to propagate classical trajectories. Ab Initio Calculations. To model the interactions between HCl and the OH-SAMs, we calculated intermolecular potentialenergy curves for the HCl-CH3OH system by scanning the HCl-CH3OH center-of-mass coordinate from the asymptote to repulsive energies of about 400 kJ 3 mol-1 using electronic-structure methods. The separation between the points of each scan is 0.1 Å. The equilibrium geometry of the CH3OH molecule as well as the equilibrium HCl bond length have been held fixed throughout the scans. For adequate coverage of the potential-energy surface, we investigated a total of eight different approaches of HCl to the CH3OH molecule, which are depicted in Figure 1. In the “facial” approaches, f(H-C) and f(Cl-C), the HCl approaches the center of the triangle formed by the three methyl hydrogens, with the scan axis collinear to the C-O bond. In “bisect” approaches, b(H-O) and b(Cl-O), HCl approaches the hydroxyl oxygen collinearly to the scan axis that bisects the methanol C-O-H angle ( — (CO-H) = 108.8, so the angle formed by the H-Cl and C-O bond axes is 125.6). In the configurations denoted h(H-H) and h(Cl-H), HCl approaches the hydroxyl hydrogen collinear to the O-H bond. Finally, we scanned the HCl molecule collinear to the methyl hydrogen located anti to the hydroxyl, denoted a(H-H) and a(Cl-H) approaches.
Figure 2. Calculated intermolecular potential energy for the HCl-O(H)CH3 system as a function of the H(Cl)-O distance for the b(H-O) approach orientation (see Figure 1) with various basis sets.
The electronic Schr€odinger equation has been solved at each step of the scans using second-order M€oller-Plesset perturbation theory in combination with the double-, triple-, and quadruple-ζ family of correlation-consistent basis sets of Dunning,44-46 augmented with diffuse functions (aug-cc-pVDZ, aug-cc-pVTZ, and aug-cc-pVQZ, respectively). Coupled cluster calculations with explicit single and double excitations and perturbative treatment of triple excitations (CCSD(T)) have also been carried out with the aug-cc-pVDZ basis set. The focal-point approach of Allen and coworkers47,48 has been used to estimate CCSD(T) energies with the aug-cc-pVTZ and aug-cc-pVQZ basis sets from MP2 calculations with those basis sets. (Hereafter, focal-point CCSD(T) energies will be referred to as fp-CCSD(T).) The focal-point approach is based on the observation that the difference between MP2 and CCSD(T) energies is essentially independent of basis set for highquality basis sets. By calculating the differences between MP2 and CCSD(T) energies with an affordable basis set (e.g., augcc-pVDZ), for a variety of intermolecular geometries, CCSD(T) energies with larger basis sets (e.g., aug-cc-pVTZ, aug-cc-pVQZ) can be estimated from MP2 calculations at those geometries (fpCCSD(T) data). The legitimacy of this approach has been demonstrated by our group and others in the description of various pairs of interacting species.43,49,50 Complete basis set (CBS) estimates are obtained for both MP2 and fp-CCSD(T) calculations using the two-point extrapolation procedure of Halkier et al.51 We removed the basis-set superposition error using the standard counterpoise method52 in all of the points of the calculated potential-energy surfaces. The electronic-structure calculations have been carried out with the Gaussian03 suite of programs.53 Figure 2 shows the intermolecular potential-energy curve of the HCl-CH3OH approach, b(H-O), that yields the deepest well as a function of the separation between the closest interatomic distance of the interacting molecules as predicted by various ab initio methods. The intermolecular potential-energy curve shows the expected features of a steep repulsive wall at short distances and a relatively deep attractive well at longer separations due to stabilizing intermolecular interactions dominated by the formation of a hydrogen bond. An examination of the dependence of the intermolecular potential on the basis set reveals that, as expected, for MP2 calculations, an increase in the size of the basis set results in lower intermolecular energies. In addition, the location of the potential 2275
dx.doi.org/10.1021/jp1106422 |J. Phys. Chem. C 2011, 115, 2273–2283
The Journal of Physical Chemistry C
ARTICLE
wells occurs at shorter approach distances with larger basis sets. The MP2 method converges quickly with increasing basis set size, so that the difference between aug-cc-pVDZ and aug-cc-pVTZ energies is substantially larger than that between the aug-cc-pVTZ and augcc-pVQZ data. Quantitatively, the root-mean-square deviation (rmsd) between MP2/aug-cc-pVDZ and aug-cc-pVTZ energies is 6.80 kJ 3 mol-1 for the overall 315 points calculated for all eight approaches, while the rmsd between aug-cc-pVTZ and augcc-pVQZ data decreases to 1.93 kJ 3 mol-1. Comparison of the CCSD(T) and MP2 data with the augcc-pVDZ basis set reveals how treating electron correlation influences the characteristics of the potential-energy surface. We find the differences between MP2 and CCSD(T) data to be sizable, with a rmsd of 3.23 kJ 3 mol-1 over all approaches, which is nearly twice the influence of increasing the basis from the triple to the quadruple basis set. This indicates that for this system, inclusion of CCSD(T)-level data will greatly increase the accuracy of our potential. Analytic HCl/OH-SAM Potentials. Using the highest level ab initio information for the HCl-CH3OH system (fp-CCSD(T)/ CBS), we constructed two-body analytic potentials that can be used to describe the HCl/OH-SAM gas/surface potential-energy surface. We considered two different potential forms constructed as sums of two-body functions. The first functional form is the familiar Buckingham potential, where each two-body term is expressed in the form Vij ¼ Aij e - Bij rij þ
Cij n rij ij
ð2Þ
where rij is the internuclear distance between the atoms of each pair and Aij, Bij, Cij, and nij are adjustable parameters specific to each pair of atoms. We will refer to this function hereafter as “Fit A.” We used a nonlinear least-squares procedure to obtain the values of the Aij, Bij, Cij, and nij parameters that minimize the relative differences between analytic energies obtained with the Buckingham potentials and the fp-CCSD(T)/CBS data for all of the points contained in the potential-energy-surface scans described before. While the Buckingham function normally includes an inverse r6 term, designed to describe attraction due to dispersion, by allowing this exponent to vary, the generalized exponential function is in principle flexible enough to incorporate longer-ranged interactions which arrise from partial localization of charge within the molecules. We tested this idea with a second fit which includes addition of an explicit Coulombic term, making the function of the form Vij ¼ Aij e - Bij rij þ
Cij qi qj nij þ rij 4πε0 rij
ð3Þ
where the additional terms qi and qj are the partial charges of the interacting atomic pair and ε0 is the permittivity of free space. We will denote this function “Fit B.” Atomic partial charges for the methanol molecule were taken from the OPLSAA force field and are -0.683, þ0.145, þ0.418, and þ0.060 e for the O, C, H(O), and H(C) atoms, respectively.39 Charges on the HCl molecules were taken to be þ0.176 and -0.176 e for the H and Cl atoms, respectively.41 We found that with the inclusion of the Coulombic term, allowing the nij exponent to vary did not result in significant improvement of the overall fit and additionally resulted in no significant differences in the molecular dynamics results of test runs. We did find, however, slight improvement if we allowed nij to vary only for the Cl-H(C) and Cl-H(O) pairs. For both functional forms, we found the H-H(O) and H-H(C) terms to be almost entirely repulsive,
Table 1. Parameters of the Analytic Potentials Describing the HCl/OH-SAM Interactions for Fit Aa pair
Aij
Bij
Cij
nij
Cl-C
322 819.9
3.747
-30 770.076
6.379
Cl-O Cl-H(C)
218 063.1 4205.6
4.294 3.113
-20 791.229 -29.980
10.070 7.22
Cl-H(O)
4863.8
3.253
-52.781
6.857
H-C
4912.4
3.326
-118.257
9.569
H-O
4743.5
3.516
-111.921
1.549
H-H(C)b
0.0
n/a
24.112
16.046
H-H(O)b
0.0
n/a
30.932
6.012
a
Units are such that if internuclear distances are given in Angstroms, then the potential energy is in kcal 3 mol-1. b The H-H(C) and H-H(O) pairs were found to be almost wholly repulsive; so, we simplified this term by setting Aij to zero; when Cij is positive the inverse rnij term in the gas/surface potential becomes repulsive.
so for those terms, either the Aij (for Fit A) or the Cij (for Fit B) parameter has been set to zero to simplify the function. The two-body Buckingham potentials provide an adequate representation of all eight approaches of HCl to the CH3OH molecules. Each of the fits includes a total of about 315 points distributed roughly evenly between the eight approaches and covers regions of the potential-energy surface from the asymptote up to ∼400 kJ 3 mol-1. Because we would like to implement these potentials in molecular dynamics simulations at superthermal (∼1 eV) scattering energies, particular care was taken to maintain an accurate description of the repulsive walls while ensuring accuracy in the region of the attractive well. To achieve this, different weighting schemes were implemented for each fit. In Fit A, we weighed points whose energies were >350 kJ 3 mol-1 and
