Isotope Effects on the Vibrational Relaxation and Multidimensional

Nov 14, 2008 - Joseph Tomkins and Gabriel Hanna. The Journal of Physical ... Porscha L. McRobbie , Gabriel Hanna , Qiang Shi and Eitan Geva. Accounts ...
0 downloads 0 Views 1MB Size
Download PDF
J. Phys. Chem. B 2008, 112, 15793–15800

15793

Isotope Effects on the Vibrational Relaxation and Multidimensional Infrared Spectra of the Hydrogen Stretch in a Hydrogen-Bonded Complex Dissolved in a Polar Liquid Gabriel Hanna* and Eitan Geva* Department of Chemistry and FOCUS Center, UniVersity of Michigan, Ann Arbor, Michigan 48109-1055 ReceiVed: August 14, 2008; ReVised Manuscript ReceiVed: September 19, 2008

Isotope effects on rate processes and spectra are often used to elucidate the nature of the interactions underlying molecular structure and dynamics. In this paper, we present a computational study of the effect of substituting hydrogen by deuterium in a solvated hydrogen-bonded complex on the rates of the various processes involved in the vibrational relaxation of the hydrogen/deuterium stretch and on the corresponding 1D and 2D infrared spectra. The vibrational relaxation is simulated via the mixed quantum-classical Liouville method, where the proton/deuteron is treated quantum mechanically while the other particles are treated in a classical-like manner. We find that the vibrational relaxation pathway and the rates of the various steps in it are similar for the deuterium and hydrogen stretches. However, we also find that isotope substitution modifies the 1D and 2D spectra of the stretch in a qualitative manner. The differences between the spectra are explained in terms of the narrowing and broadening of the fundamental and overtone transition frequency ranges, respectively, and the smaller transition dipole moments in the case of the deuterium stretch. Our results demonstrate that isotope substitution may have a rather dramatic effect on the infrared spectra of a vibrational mode strongly coupled to its environment even though the rate and pathway of the underlying vibrational relaxation may not be overly sensitive to it. I. Introduction The signature of vibrational relaxation processes on multidimensional infrared (IR) spectra has received much attention over the past decade.1-4 However, relatively little attention has been given to the unique features of the IR spectra of vibrational modes which are strongly coupled to their environment. In a previous paper, we studied the vibrational relaxation pathway of an excited hydrogen stretch (H-stretch) within a model of a moderately strong H-bonded complex dissolved in an aprotic polar liquid.5 We showed that in this important example of a vibrational mode which is strongly coupled to its environment, vibrational relaxation is generally a multistep process that involves (1) solvation on the significantly different potential surfaces that correspond to the ground and excited vibrational states, (2) nonadiabatic transitions between these two surfaces, and (3) intramolecular proton transfer that leads to tautomeric equilibrium on the ground-state potential surface. In a more recent paper, we studied the signatures of the solvation and proton-transfer processes in the above-mentioned model system on the corresponding 1D and 2D IR spectra.6 We showed that the 1D spectrum consists of two high frequency peaks that correspond to the ionic and covalent tautomers, as well as a low-frequency peak which represents contributions from unstable configurations in the vicinity of the transition state. We also showed that the transition-state peak results from the strong and rapid enhancement of the transition dipole moment in this region (a non-Condon effect). On the basis of those rather unusual spectral features, we proposed that it may be possible to use 2D IR spectroscopy in order to probe and resolve such processes as the relaxation from the transition state to either the ionic or covalent wells, in real time (although it should be * To whom correspondence should be addressed.

noted that the low-frequency region is difficult to access experimentally via currently available laser technologies). In this paper, we consider the effect on the vibrational relaxation pathway and IR spectra of substituting the H by deuterium (D) in the above-mentioned model system. This study is motivated by the following questions, which, to the best of our knowledge, have not been previously explored in the context of a vibrational mode strongly coupled to its environment: •What effect does isotope substitution have on the rates of the various steps involved in the vibrational relaxation pathway? •How sensitive are the 1D and 2D IR spectra to isotope substitution? •Is the relationship between the vibrational relaxation and spectra the same for both isotopomers? It should be noted that addressing these questions is particularly relevant in light of the many recent studies aimed at learning about the structure and dynamics of H-bonded systems by measuring the multidimensional spectra of the corresponding deuterated variants of those systems.4,7-31 It should also be noted that, at least in principle, the H to D isotope substitution is expected to give rise to an appreciable isotope effect because of the large mass ratio. The remainder of this paper is organized as follows. The model is briefly outlined in section II. The theoretical and computational framework used for simulating vibrational relaxation and IR spectra is outlined in section III. The results are presented in section IV. Finally, the main conclusions are summarized and discussed in section V. II. Model In this section, we provide a short outline of the Azzouz-Borgis (AB) model as used in the present study (a more detailed description is available in refs 5 and 32-34). Within this model, it is assumed that the proton/deuteron moves along a one-

10.1021/jp8072816 CCC: $40.75  2008 American Chemical Society Published on Web 11/14/2008

15794 J. Phys. Chem. B, Vol. 112, No. 49, 2008

Hanna and Geva

dimensional axis that connects the donor and acceptor. The donor, A, and acceptor, B, are modeled as single particles and parametrized to represent phenol and trimethylamine, respectively. The functional form of the intramolecular potential surface as a function of the A-B and A-H/D distances is adopted from ref 35. Polarizability effects are accounted for by assuming that the charges on A and B are explicitly dependent on the position of the proton/deuteron. Increasing the equilibrium A-B distance modifies the potential along the deuteron coordinate from a single well at short A-B distances to a double well at large A-B distances. The double-well form of the potential gives rise to tautomeric equilibrium between a covalent form and an ionic form of the complex, AD-B h A--D+B. In the covalent form, the deuteron is attached to the phenol (A), and the complex has a small (2.5 D) dipole moment. In the ionic form, the deuteron is attached to trimethylamine (B), and the resulting charge separation gives rise to a significantly larger (10.5 D) dipole moment. In this study, we set the equilibrium A-B distance to 2.7 Å, which corresponds to a moderately strong H-bond (i.e., a double-well profile with a relatively low barrier along the proton/deuteron coordinate). The solvent consists of a liquid of methyl chloride molecules, which are modeled as rigid dipolar diatomic molecules. The complex-solvent and solvent-solvent interactions are modeled in terms of site-site Lennard-Jones (LJ) and Coulomb interactions, and the corresponding force field parameters are adopted from ref 36. Although the covalent form, AD-B, is favored in the absence of the solvent, the ionic form, A-D+B, becomes more stable when the complex is dissolved in a dipolar solvent. III. Theory and Simulation Techniques The theoretical framework used in this paper for simulating vibrational relaxation and computing spectra was described in detail in refs 5 and 6 and is only briefly outlined below for the sake of completeness. A. Mixed Quantum-Classical Liouville Dynamics. The mixed quantum-classical Hamiltonian of the overall system is given by

1, and 2 correspond to the ground state, first excited state, and second excited state, respectively. The vibrational relaxation of the H/D stretch is simulated within the framework of the mixed quantum-classical Liouville (MQCL) method.5,37-48 Starting at equilibrium on the groundstate potential surface, the initial state of the overall system after impulsive excitation to the excited state is given by

FijW(Q, P, 0) ) δ(i, 2)δ(j, 2)FW 22(Q, P, 0)

W (Q, P, 0) corresponds to equilibrium on the groundHere, F22 state potential surface in the Wigner representation, which can be approximated by its classical limit

FW 22(Q, P, 0) ≈

exp{-β[KB(P) + E0(Q)]}

∫ dQ ∫ dP exp{-β[KB(P) + E0(Q)]} (4)

The observable of interest is the excited-state population, which is represented by the operator Pˆ ) |2; Q〉〈2; Q|. Its expectation value at time t is given by

P2(t) ≡ 〈Pˆ〉(t) )

∫ dQ ∫ dPFW22(Q, P, 0)PW22(Q, P, t) (5)

W (Q, P, t) is obtained by solving the MQCL equation Here, P22 with the initial conditions PW ij (Q, P, 0) ) δ(i, 2)δ(j, 2). In the adiabatic representation, the MQCL equation for an overall system observable Bˆ is given by43,47

{

}

Fi(Q) + Fj(Q) ∂ W B (Q, P, t) ) U · ∇Q + · ∇P BW ij (Q, P, t) + ∂t ij 2 1 U · d/jk(Q) 1 + S/jk · ∇P BW iωij(Q)BW ij (Q, P, t) + ik (Q, P, t) + 2 k

∑

[

∑ U · d (Q)[1 + 21 S ik

ˆ (qˆ, pˆ, Q, P) ) KB(P) + KˆH/D(pˆ) + Vˆ(Q, qˆ) H

]

]B (Q, P, t)

ik · ∇P

k

W kj

(6)

(1)

Here, as in the rest of the paper, vectors are represented by boldface symbols (e.g., A), while operators are capped (e.g., Aˆ). P ≡ (P1, P2, ...) and Q ≡ (Q1, Q2, ...) are the coordinates and momenta of the solvent and A and B groups, collectively referred to as the bath degrees of freedom (DOF), while qˆ and pˆ are the coordinate and momentum operators of the proton/ deuteron. KB(P) is the kinetic energy of the bath DOF, KˆH/D(pˆ) is the kinetic energy operator of the proton/deuteron, and Vˆ(Q, qˆ) is the overall potential energy operator. The H/D stretch is assumed to follow the bath DOF adiabatically so that the vibrational energy levels and wave functions correspond to the eigenvalues and eigenfunctions, respectively, of the adiabatic ˆ H/D(qˆ, pˆ, Q) ≡ KˆH/D(pˆ) + Vˆ(Q, qˆ) Hamiltonian, H

ˆ H/D |j;Q〉 ) Ej(Q)|j;Q〉 H

(3)

(2)

Here, {|j; Q〉} are the adiabatic vibrational states, and {Ej(Q)} are the corresponding vibrational energy levels, which also correspond to the potential energy surfaces on which the classical DOF evolve. We will adopt the convention that j ) 0,

Here, U ) (U1, U2, ...) ) (P1/M1, P2/M2, ...) are the velocities of the bath DOF, where (M1, M2, ...) are the corresponding masses, ωij(Q) ) [Ei(Q) - Ej(Q)]/p is the angular transition frequency, Fj(Q) ) -〈j; Q|∇QV(Q, qˆ)|j; Q〉 are the HellmannFeynman forces, dmn(Q) ) 〈m; Q|∇Q|n; Q〉 are the nonadiabatic coupling coefficients, and Smn ) pωmn(Q)[U · dmn(Q)]-1dmn(Q). The first term on the right-hand side of eq 6 corresponds to classical propagation of BW ij (Q, P, t) on the potential surface [Ei(Q) + Ej(Q)]/2. Thus, the “diagonal” variables, BW ii (Q, P, t), are propagated on the corresponding adiabatic potential surfaces, Ei(Q), while the “off-diagonal” variables, BW ij (Q, P, t) (with i * j), are propagated on potential surfaces that correspond to the arithmetic average of Ei(Q) and Ej(Q). The second term on the right-hand side of eq 6 leads to the creation of a cumulative phase factor, exp[i∫dτωij(τ)], during the time evolution of the off-diagonal variables (i.e., when i * j). The last two terms account for nonadiabatic transitions via the coupling of W W BW ij (Q, P, t) to Bik (Q, P, t) and Bkj (Q, P, t). B. 1D and 2D IR Spectra. The 1D IR spectrum is given by49

I(ω) ) Re

∫0∞ dtJ(t)eiωt

(7)

Isotope Effects on the Hydrogen Stretch

J. Phys. Chem. B, Vol. 112, No. 49, 2008 15795

where J(t) is the linear optical response function (ORF), which is assumed to be given by the following commonly used semiclassical expression25,49-52

J(t) )

∫ dQ0 ∫ dP0

exp[-βH0(Q0, P0)] V01(t)V10(0) × Z0

∫0t dτω10(τ)] ≡ 〈V01(t)V10(0) × t exp[-i ∫0 dτω10(τ)]〉 exp[-i

(8)

Here, H0(Q0, P0) ) KB(P0) + E0(Q0) and Z0 ) ∫dQ0∫dP0 exp[-βH0(Q0, P0)] are the classical bath Hamiltonian and partition function that correspond to the ground vibrational state; Vij(t) ) 〈i; Qt|qˆ|j; Qt〉 is the transition dipole moment; pωij(t) ) Ei(Qt) - Ej(Qt) is the energy gap between vibrational states i and j; and Qt represents the trajectory of the bath DOF obtained via classical propagation on the ground state potential surface, E0(Q), with the initial conditions (Q0, P0). The Condon approximation amounts to neglecting the fluctuations in the transition dipole moments due to the dependence of the wave functions on the bath configuration. It corresponds to assuming that the transition dipole moments in eq 8 can be substituted by their average values. The second-order cumulant (2OC) approximation corresponds to truncating the cumulant expansion of 〈exp(i∫0t dτω10(Qτ))〉 at second order. Explicit expressions for the linear ORF under the Condon and 2OC approximations can be found in ref 6. 2,6,25,53 The 2D IR spectra are defined by

∫0∞ dt1 ∫0∞ dt3{ei(ω t +ω t )Rnr(t3, t2, t1) +

I(ω1, t2, ω3) ≡ Re

11

33

ei(-ω1t1+ω3t3)Rr(t3, t2, t1)}

(9)

The so-called nonrephasing (Rnr) and rephasing (Rr) signals are assumed to be given by the following commonly used semiclassical expressions Rnr(t3, t2, t1) ) 2〈V10(0)V01(t1)V10(t1 + t2)V01(t1 + t2 + t3) × exp[-i

∫

t1

0

dτω10(τ) - i

∫

t1+t2+t3

t1+t2

dτω10(τ)]〉 -

〈V10(0)V01(t1)V21(t1 + t2)V12(t1 + t2 + t3) × exp[-i

∫

t1

0

dτω10(τ) - i

∫

t1+t2+t3

t1+t2

dτω21(τ)]〉

(10)

and Rr(t3, t2, t1) ) 2〈V01(0)V10(t1)V10(t1 + t2)V01(t1 + t2 + t3) × exp[i

∫

t1

0

dτω10(τ) - i

∫

t1+t2+t3

t1+t2

dτω10(τ)]〉 -

〈V01(0)V10(t1)V21(t1 + t2)V12(t1 + t2 + t3) × exp[i

∫

t1

0

dτω10(τ) - i

∫

t1+t2+t3

t1+t2

dτω21(τ)]〉

(11)

respectively. A more detailed discussion on the measurement of these signals as well as the corresponding expressions under the Condon and 2OC approximations can be found in ref 6 and references therein. C. Simulation Techniques. The simulation techniques used in this paper are similar to those described in refs 5, 6, and 33 and will therefore be only briefly outlined below. Simulations

Figure 1. The ground (0), first excited (1), and second excited (2) free-energy surfaces as a function of the solvent polarization for ADB (red) and AHB (blue).

were performed with a single AHB/ADB complex and 255 methyl chloride molecules in a cubic simulation box with periodic boundary conditions at a temperature and density of 250 K and 0.012 Å-3, respectively. The LJ potentials were spherically truncated at Rc ) 13.8 Å and shifted accordingly. The Coulomb potentials were smoothly truncated to zero at Rc ) 13.8 Å. The adiabatic surfaces were computed “on-the-fly” by diagonalizing the protonic/deuteronic Hamiltonian for different instantaneous bath configurations using the LAPACK DSYGV routine. Adiabatic dynamics was simulated by solving the classical equations of motion via the velocity Verlet algorithm with a time step of 1 fs. Initial bath configurations were generated via a constant-temperature MD simulation using a Nose´-Hoover thermostat.54 Simulations of nonadiabatic dynamics were based on numerically solving the MQCL equation via the sequential short-time propagation algorithm and employed the momentum jump approximation.5,33,55 The 1D and 2D Fourier transforms required for computing the 1D and 2D IR spectra, respectively, were carried out numerically via the FFT routine using 2048 grid points and a time interval of 1 fs. IV. Results A. Free-Energy Surfaces and Transition Frequencies. Figure 1 shows the ground, first excited, and second excited free-energy surfaces (FESs) for ADB (red) and AHB (blue), as a function of the solvent polarization, ∆E, which is defined as the difference between the solvent electrical potentials at the minima of the covalent and ionic wells and serves as a convenient collective solvent coordinate.5,6,56,57 The three FESs are very different from each other for both ADB and AHB, which is a manifestation of the strong coupling between the vibrational mode and the bath. The ground-state FES has a very similar double-well form in both cases, which reflects coexistence between the ionic and covalent tautomers. The equilibrium compositions of the two species are also similar in both cases, namely, 69% ionic and 31% covalent in ADB versus 65% ionic and 35% covalent in AHB. The free-energy barrier in the case of ADB is observed to be somewhat higher than that in the case of AHB, thereby giving rise to a kinetic isotope effect on the proton/deuteron transfer rate constant. More specifically, the adiabatic deuterontransfer rate constant is 0.12 ps-1, which is somewhat slower than the adiabatic proton-transfer rate constant, which is 0.16 ps-1.34 The first excited FESs for ADB and AHB are observed to be rather similar at intermediate polarizations. However, the first excited FES of ADB drops below that of AHB at the high and low polarization regions, which correspond to the regions of

15796 J. Phys. Chem. B, Vol. 112, No. 49, 2008

Figure 2. Upper panel: The distributions of transition frequencies ω10 and ω21 for ADB and AHB, as obtained from equilibrium MD simulations on the ground-state potential surface. Lower panel: The averaged transition dipole moments V10 and V21 for ADB and AHB as a function of the corresponding transition frequencies (calculated from equilibrium MD simulations on the ground-state potential surface). Note that the blue curves are plotted against the upper x-axis, while the red curves are plotted against the lower x-axis; not also that the lower x-axis is scaled by a factor of 1/(2)1/2 relative to the upper x-axis.

the ionic and covalent wells on the ground-state FES, respectively. The second excited FES of ADB is significantly lower than that of AHB and is observed to have a qualitatively different shape within a wide range of solvent polarizations. More specifically, while the second excited FES of AHB has the shape of a very broad single well, that of ADB has the shape of a double well with the ionic well slightly more stable and separated by a discernible free-energy barrier from the covalent well. We next consider the fundamental transition frequency between the ground and first excited state, ω10, and the overtone transition frequency between the first and second excited states, ω21. The distribution of those transition frequencies as obtained from equilibrium simulations on the ground-state potential surface are shown in the upper panel of Figure 2. The peak of the ω10 distribution in the case of ADB is seen to be red-shifted by a factor of 1/(2)1/2 relative to that in the case of AHB, which is consistent with the expected shift in the transition frequency of a harmonic stretch associated with doubling the mass. However, overall rescaling of the ω10 distribution in the case of AHB by a factor of 1/(2)1/2 along the ω10 axis yields an ω10 distribution which is significantly wider than the actual ω10 distribution in the case of ADB. Thus, the ω10 distribution in the case of ADB is observed to be both red-shifted and effectively narrower relative to that in the case of AHB. Although the ω21 distributions have a less pronounced peak, one can still observe that the peak in the case of ADB is also red-shifted relative to that in the case of AHB by a factor of ∼1/(2)1/2. However, here too, the isotope effect cannot be described solely in terms of rescaling the ω21 axis. More specifically, rescaling the ω21 distribution in the case of AHB by a factor of 1/(2)1/2 along the ω21 axis yields an ω21 distribution which is effectively narrower than the actual ω21 distribution in the case of ADB. Thus, substituting H for D is observed to effectively narrow the ω10 distribution and broaden the ω21 distribution. As we will show below, this observation is key for understanding the effect of isotope substitution on the corresponding IR spectra. It should be noted that the ground-state FESs of ADB and AHB almost coincide (see Figure 1), which implies that the

Hanna and Geva

Figure 3. (A) A schematic view of the four steps involved in the vibrational relaxation of the H/D stretch, following impulsive excitation from the ground state to the first excited state. (B) The vibrational energy of the first excited state as a function of time during solvation on the excited-state surface. (C) The nonadiabatic relaxation of the excited-state population starting at equilibrium on the first excited potential surface. (D) The vibrational energy of the ground state as a function of time during the solvation on the ground-state surface that follows the nonadiabatic transition from the first excited state to the ground state.

transition frequency distributions shown in the upper panel of Figure 2 are based on sampling similar subsets of bath configurations. Thus, the difference in the width of the distributions must reflect the fact that isotope substitution affects the way in which the transition frequencies respond to changes in the solvent configuration. Indeed, the above-mentioned trends are consistent with the fact that for bath configurations most likely to be sampled at equilibrium on the ground-state potential, the first excited FES changes more slowly as a function of polarization in the case of ADB, while the opposite is true in the case of the second excited FES (see Figure 1). Hence, the corresponding ω10 and ω21 distributions are effectively narrower and broader, respectively, in comparison to AHB. B. Vibrational Relaxation. Adopting the same scheme as that in ref 5, we consider the pathway and rates of vibrational relaxation following impulsive photoexcitation of the system at equilibrium on the ground-state potential surface to the first excited state. As in the case of AHB, we have found that vibrational relaxation in the case of ADB follows a pathway that involves the following four consecutive steps (see Figure 3A): (1) Solvation on the first excited-state potential surface. (2) A nonadiabatic transition between the first excited potential surface to the ground potential surface. (3) Solvation on the ground potential surface. (4) Deuteron-/proton-transfer reaction to obtain equilibrium ionic and covalent subpopulations. The dynamics during steps (1)-(3) as obtained for ADB and AHB is shown in Figure 3B-D, respectively. The main conclusion based on Figure 3 is that isotope substitution has a rather small effect on the time scales of steps (1)-(3) of the vibrational relaxation pathway. The ionic/covalent branching ratio at the end of step (3) is found to be equal to 0.85 in the case of ADB, which is about a factor of 2.5 smaller than the corresponding equilibrium branching ratio of 2.2.34 Thus, at the end of step (3), the system is in a nonequilibrium state. Relaxation to the equilibrium ionic/covalent branching ratio is

Isotope Effects on the Hydrogen Stretch

Figure 4. 1D IR spectra of the D (red) and H (blue) stretches as obtained within the non-Condon (upper panel), Condon (middle panel), and 2OC (lower panel) treatments. Note that the blue curves are plotted against the upper x-axis, while the red curves are plotted against the lower x-axis; note also that the lower x-axis is scaled by a factor of 1/(2)1/2 relative to the upper x-axis.

achieved via proton/deuteron transfer (step (4)). As was shown in ref 34, isotope substitution also has a rather small effect on the time scale of step (4), namely, the deuteron-transfer rate constant is found to be 25% smaller than the proton-transfer rate constant. It should be noted that the vibrational energy exhibits a damped oscillatory behavior during the solvation steps (steps (1) and (3)). The oscillations can be attributed to intramolecular energy exchange between the H/D stretch and the A–B stretch, while the damping corresponds to energy loss to the solvent. The fact that the oscillations in the excited-state solvation (step (1)) are smaller in the case of the D stretch can be attributed to its larger mass, which makes it less sensitive to motion along the A–B stretch. In light of this, the fact that the opposite trend is observed in the ground-state solvation may appear surprising at first sight. However, it should be noted that in this case, the initial states of the bath are obtained immediately after the nonadiabatic transition has taken place, which do not correspond to thermal equilibrium and are isotope-dependent. Indeed, we have found that the amplitude of the A–B stretch right after the nonadiabatic transition is significantly larger than its typical values at thermal equilibrium and even more so in the case of the D stretch than in the case of the H stretch. Thus, the larger amplitude of the oscillation in the case of the ground-state solvation of the D stretch can be attributed to the fact that the A–B stretch starts out in a more excited state. C. One-Dimensional Spectra. The absorption spectra of the D and H stretches as obtained within the non-Condon, Condon, and 2OC treatments are compared in Figure 4. We start by considering the non-Condon 1D spectrum of the D stretch, which is seen to be qualitatively different from that of the H stretch (see upper panel in Figure 4). More specifically, the spectrum of the D stretch consists of a single relatively narrow band which is centered at ∼1800 cm-1 and can be attributed to contributions from the covalent and ionic subpopulations (although their individual contributions cannot be resolved). In

J. Phys. Chem. B, Vol. 112, No. 49, 2008 15797 contrast, the non-Condon spectrum of the H stretch consists of three similarly narrow bands centered at ∼200, ∼2250, and ∼2500 cm-1. The two, partially overlapping, high-frequency bands can be assigned to the covalent (∼2250 cm-1) and ionic (∼2500 cm-1) subpopulations, which can be resolved because of motional narrowing. The low-frequency band results from a strong enhancement of the transition dipole moment in the lowfrequency region, which is attributed to bath configurations in the transition-state region.6 The origin of the transition-state peak in the case of the H stretch can be traced back to the strong and rapid enhancement of the transition dipole in the low-frequency region, which compensates for the low probability density of bath configurations in this region. Thus, the absence of a noticeable transitionstate peak in the 1D spectrum of the D stretch can be traced back to the fact that the transition dipole moment V10 in the vicinity of the transition state is approximately nine times smaller than that of the H stretch (see lower panel of Figure 2). We next turn to the 1D spectra as obtained within the framework of the Condon approximation (see middle panel of Figure 4). For both H and D stretches, the Condon approximation gives rise to a very broad shoulder which extends from the ionic and covalent high-frequency bands down to very low transition frequencies which correspond to transition-state configurations. It should be noted that the disappearance of this broad shoulder in the non-Condon spectrum of the H/D stretch represents a very pronounced non-Condon effect. More specifically, the fact that the transition dipole moment is weak in the intermediate transition frequency region leads to a solventinduced selection rule that renders transitions in this region forbidden. Furthermore, the fact that the transition dipole moment V10 in this region is approximately four times smaller for the D stretch than that for the H stretch (see lower panel of Figure 2) explains why this shoulder is less pronounced in the 1D spectrum of the D stretch. Insight into the origin of the non-Condon effects can be obtained by examining the wave functions, ψi ) 〈q|i; Q〉, as a function of the quantum coordinate, q, at representative bath configurations that correspond to different values of ∆E (see Figure 5). To this end, we consider two bath configurations, one which corresponds to the vicinity of the transition state (see top right panel of Figure 5) and the other which corresponds to the intermediate frequency region (see top left panel of Figure 5). For both configurations, the ground- and first excited-state wave functions are more localized in the case of the D stretch than in the case of the H stretch. As a result, they overlap more poorly such that the magnitude of the integrand ψ1qψ0 in the case of the D stretch is lower than that in the case of the H stretch (see bottom panels of Figure 5). In turn, the transition dipole moments V10 are smaller for the D stretch than those for the H stretch. Thus, the absence of the low-frequency band and the lack of intensity in the intermediate frequency region in the non-Condon 1D spectrum of the D stretch can be traced back to the localization of the wave functions for the heavier isotope. Finally, we examine the absorption spectra of the D and H stretches as obtained within the 2OC treatment (see bottom panel of Figure 4). In both cases, the 2OC spectrum is seen to be inconsistent in almost every respect when compared to the Condon and non-Condon spectra. However, the differences appear to be somewhat smaller in the case of the D stretch because of the smaller range of transition frequencies covered and the lack of such spectral features as a pronounced shoulder in the Condon spectrum and a multiband structure in the nonCondon spectrum.

15798 J. Phys. Chem. B, Vol. 112, No. 49, 2008

Hanna and Geva equilibrium distributions of the fundamental (ω10) and overtone (ω21) transition frequencies (see upper panel of Figure 2), alongside the dependence of the transition dipole moments V10 and V21 on the transition frequencies ω10 and ω21, respectively (see lower panel of Figure 2). More specifically, the reason for the weaker negative feature in the case of the D stretch can be explained by the facts that 〈V21〉 is somewhat smaller than 〈V10〉 and that the intensity of the peak in the ω10 distribution is much higher than that in the ω21 distribution. In contrast, the reason for the stronger negative feature in the case of the H stretch can be explained by the facts that 〈V21〉 is significantly larger than 〈V10〉 and that the intensity of the peak in the ω21 distribution is slightly higher (and broader) than that in the ω10 distribution. It should be noted that the inability to resolve the individual contributions of the ionic and covalent species in the case of the D stretch implies that one cannot use 2D IR spectroscopy in order to probe the deuteron-transfer process. This should be contrasted with the case of the H stretch, where the contributions of the ionic and covalent species can be resolved in the 1D spectrum such that proton transfer can be probed via 2D IR spectroscopy.6 Indeed, we have observed that the 2D spectra of the D stretch do not change significantly at t2 > 250 fs, which should be contrasted with the case of the H stretch where we previously observed the emergence of off-diagonal peaks due to proton transfer. V. Discussion and Conclusions

Figure 5. Upper panels: The vibrational wave functions that correspond to the ground (solid) and first excited (dashed) vibrational energy levels of the H (blue) and D (red) stretches at two representative bath configurations. The bath configuration in the left panels (∆E ∼ 0.01 eC/Å) corresponds to the region between the covalent well and transition state, and that in the right panels (∆E ∼ 0.0138 eC/Å) corresponds to the transition-state region. Lower panels: The integrand of the transition dipole moment V10 for the H (blue) and D (red) systems at these two representative bath configurations.

D. Two-Dimensional Spectra. The 2D IR spectra of the D stretch as obtained within the Condon and non-Condon treatments are shown in Figures 6 and 7, respectively, at t2 ) 0 and 250 fs (the 2OC treatment was found to be very inaccurate for this system and is therefore not shown). A comparison between the 2D spectrum of the D stretch and the corresponding 2D spectrum of the H stretch6 gives rise to the following observations: • The peak intensity of the spectral signal in the case of the D stretch is stronger than that in the case of the H stretch. • The 2D Condon and non-Condon spectra of the D stretch are rather similar. This should be contrasted with the very pronounced non-Condon effects in the case of the H stretch that lead to qualitatively different 2D spectra. • The negative spectral feature, which is attributed to the Liouville pathway that involves the overtone transition, is found to be much weaker than the positive feature, which is attributed to Liouville pathways that only involve the fundamental transition. This should be contrasted with the 2D spectrum of the H stretch, where the negative feature was found to be dominant. The lower peak intensity in the case of the H stretch 2D spectra (relative to the D stretch 2D spectra) can be explained by the smearing of the signal over a larger range of frequencies due to the pronounced shoulder in the ω10 distribution. The similarity between the Condon and non-Condon spectra in the case of the D stretch can be explained by the lack of a pronounced shoulder in the ω10 distribution. In order to understand the last observation, it is useful to consider the

In this paper, we have considered the effect of isotope substitution on the pathway and rates of vibrational relaxation and IR spectra in the case of a vibrational mode which is strongly coupled to its environment. To this end, we have performed a systematic analysis of the isotope effect in the case of a hydrogen-bonded complex in a polar liquid where the hydrogen is substituted by deuterium. Our main findings are summarized and discussed below. We have found that isotope substitution has a rather small effect on the vibrational relaxation pathway, as well as on the rates of the individual steps which constitute it. More specifically, the rates of solvation on the excited and ground potential surfaces and of the nonadiabatic transition between them are not significantly affected by isotope substitution, and the protontransfer rate only slows down by about 25%. The insensitivity of the solvation rates to isotope substitution can be attributed to the similarity between the relevant regions of the ground and first excited potential energy surfaces that correspond to the two isotopomers. It should also be noted that in the case of this strongly coupled system, solvation processes allow the bath to adjust so as to minimize the energy mismatch between the ground and excited states. As a result, the two states become strongly coupled, and the probability that the bath will induce a transition between these states becomes high. The actual transition rate depends on the nonadiabatic coupling strength

|

| |

U · d10(Q) ) U ·

〈1;Q|∇QV(Q, qˆ)|0;Q〉 pω10(Q)

|

(12)

The distribution of transition frequencies ω10 at equilibrium on the excited state in the case of the D stretch is peaked at a lower value than that in the case of the H stretch. As a result, the denominator of eq 12 is smaller for the D stretch, thereby favoring the heavier isotope. However, the numerator of eq 12 is actually higher for the H stretch than that for the D stretch,

Isotope Effects on the Hydrogen Stretch

J. Phys. Chem. B, Vol. 112, No. 49, 2008 15799

Figure 6. Left panels: 2D IR spectra of the D stretch at t2 ) 0 and 250 fs, as obtained via the Condon treatment (top). Also shown are the individual contributions from Liouville pathways that only involve the fundamental transitions (middle) and from the Liouville pathway that also involves the overtone transition (bottom). Right panels: Same as on the left but for the H stretch.

Figure 7. Same as Figure 6 within the non-Condon treatment.

thereby favoring the lighter isotope. Thus, the insensitivity of the population relaxation rate to isotope substitution arises from a cancelation of these two opposing tendencies in the nonadiabatic coupling strength. It is interesting to contrast this result with the case of a weakly coupled system, such as H2/D2 in liquid argon, where the vibrational energy relaxation rate was also found to be rather insensitive to isotope substitution.58 In such a weakly coupled system, vibrational relaxation does not involve solvation processes, and the vibrational energy relaxation rate constant can be described within the framework of Fermi’s golden rule.

In this case, the ratio of vibrational energy relaxation rate constants for the two isotopomers is given by58

k0r1(H2) ) k0r1(D2)



µ(D2) C˜H2(ωH2) × µ(H2) C˜D (ωD ) 2 2

(13)

where µ(A2) is the reduced mass, ω(A2) is the vibrational ˜ A (ω) is the spectral density associated with the frequency and C 2 autocorrelation function of the force induced by the bath along

15800 J. Phys. Chem. B, Vol. 112, No. 49, 2008 the vibrational mode coordinate. The factor (µ(D2)/µ(H2))1/2 ) (2)1/2 weakly favors the heavier isotopomer. However, a potentially stronger source for an isotope effect is provided by the factor C˜H2(ωH2)/C˜D2(ωD2). On the one hand, the exponential gap law behavior of C˜A2(ω) favors the heavier isotope because it corresponds to a lower vibrational frequency. On the other hand, at a given frequency, C˜A2(ω) decreases with increasing isotope mass, thereby favoring the lighter isotope.58 In the case of H2/D2 in liquid argon, these opposing tendencies essentially cancel each other out, thereby giving rise to a rather negligible isotope effect on the rate of vibrational energy relaxation. On the basis of this argument, one may draw an analogy between the roles of the exponential gap law and the transition frequencies in favoring the heavier isotope for the weakly coupled system and the strongly coupled system studied in this paper, respectively. Also, the fact that the value of C˜A2(ω) is larger for the lighter isotope is analogous to the fact that the numerator of eq 12 is also larger for the lighter isotope. In contrast to the vibrational relaxation pathway and rates, we have found that both 1D and 2D IR spectra are highly sensitive to isotope substitution. It is important to note that the isotope effect on the spectra cannot be explained in terms of simply rescaling the frequencies by the factor (mD/mH)1/2 ) (2)1/2, as in the case of weakly coupled systems. More specifically, we have found that the spectra can be used for probing the rather pronounced effect that isotope substitution has on the distributions of the fundamental and overtone transition frequencies and dipole moments, which in turn reflect variations in the adiabatic potential energy surfaces and wave functions. Our results also suggest that the relationship between vibrational relaxation and the spectra can depend strongly on the isotopomer. Thus, for the model system considered here, a very similar vibrational relaxation behavior has been observed to lead to very different 1D and 2D spectra. Furthermore, we have shown that one’s ability to use spectroscopy in order to probe vibrational relaxation is strongly dependent on the choice of isotopomer. For example, for the model system considered here, the spectra of the H stretch are much richer than those of the D stretch since they include a pronounced transition-state peak and allow one to resolve the individual contributions from the ionic and covalent subpopulations. Finally, we would like to emphasize that although the analysis presented in this paper was performed in the context of a specific model, we believe that it provides useful general insight for the interpretation of isotope effects in strongly coupled systems. This is because many of our conclusions rely on general features of the potential energy surfaces that are expected to be relevant in other hydrogen-bonded systems. Extending the study of isotope effects on vibrational relaxation and spectra to other hydrogen-bonded systems as well as other strongly coupled systems is obviously highly desirable and will be the subject of future work. Acknowledgment. This project was supported by the Petroleum Research Fund and the National Science Foundation through Grants PHY-0114336 and CHE-0809506. References and Notes (1) Mukamel, S. Annu. ReV. Phys. Chem. 2000, 51, 691. (2) Khalil, M.; Demirodoven, N.; Tokmakoff, A. J. Phys. Chem. A 2003, 107, 5258. (3) Jonas, D. M. Annu. ReV. Phys. Chem. 2003, 54, 425. (4) Cho, M. Chem. ReV. 2008, 108, 1331. (5) Hanna, G.; Geva, E. J. Phys. Chem. B 2008, 112, 4048. (6) Hanna, G.; Geva, E. J. Phys. Chem. B 2008, 112, 12991.

Hanna and Geva (7) Loparo, J. J.; Roberts, S. T.; Tokmakoff, A. J. Chem. Phys. 2006, 125, 194521. (8) Loparo, J. J.; Roberts, S. T.; Tokmakoff, A. J. Chem. Phys. 2006, 125, 194522. (9) Gundogdu, K.; Bandaria, J.; Nydegger, M.; Rock, W.; Cheatum, C. M. J. Chem. Phys. 2007, 127, 044501. (10) Fecko, C. J.; Eaves, J. D.; Loparo, J. J.; Tokmakoff, A.; Geissler, P. L. Science 2003, 301, 1698. (11) Heyne, K.; Huse, N.; Nibbering, E. T. J.; Elsaesser, T. Chem. Phys. Lett. 2003, 382, 19. (12) Heyne, K.; Huse, N.; Dreyer, J.; Nibbering, E. T. J.; Elsaesser, T.; Mukamel, S. J. Chem. Phys. 2004, 121, 902. (13) Cowan, M. L.; Bruner, B. D.; Huse, N.; Dwyer, J. R.; Chugh, B.; Nibbering, E. T. J.; Elsaesser, T.; Miller, R. J. D. Nature 2005, 434, 199. (14) Asbury, J. B.; Steinel, T.; Stromberg, C.; Gaffney, K. J.; Piletic, I. R.; Goun, A.; Fayer, M. D. Phys. ReV. Lett. 2003, 91, 23742. (15) Asbury, J. B.; Steinel, T.; Fayer, M. D. J. Phys. Chem. B 2004, 108, 6544. (16) Steinel, T.; Asbury, J. B.; Corcelli, S. A.; Lawrence, C. P.; Skinner, J. L.; Fayer, M. D. Chem. Phys. Lett. 2004, 386, 295. (17) Kolano, C.; Helbing, J.; Kozinski, M.; Sander, W.; Hamm, P. Nature 2006, 444, 469. (18) Woutersen, S.; Mu, Y.; Stock, G.; Hamm, P. Chem. Phys. 2004, 266, 137. (19) Huse, N.; Bruner, B. D.; Cowan, M. L.; Dreyer, J.; Nibbering, E. T. J.; Miller, R. J. D.; Elsaesser, T. Phys. ReV. Lett. 2005, 95, 147402. (20) Park, S.; Fayer, M. D. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 16731. (21) Auer, B.; Kumar, R.; Schmidt, J. R.; Skinner, J. L. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 14215. (22) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2002, 117, 5827. (23) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 118, 264. (24) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 119, 3840. (25) Schmidt, J. R.; Corcelli, S. A.; Skinner, J. L. J. Chem. Phys. 2005, 123, 044513. (26) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2002, 117, 8847. (27) Piryatinski, A.; Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 118, 9664. (28) Piryatinski, A.; Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 118, 9672. (29) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 119, 1623. (30) Corcelli, S. A.; Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2004, 120, 8107. (31) Schmidt, J. R.; Corcelli, S. A.; Skinner, J. L. J. Chem. Phys. 2004, 121, 887. (32) Azzouz, H.; Borgis, D. J. Chem. Phys. 1993, 98, 7361. (33) Hanna, G.; Kapral, R. J. Chem. Phys. 2005, 122, 244505. (34) Hanna, G.; Kapral, R. J. Chem. Phys. 2008, 128, 164520. (35) Lippincott, E. R.; Schroeder, R. J. Chem. Phys. 1955, 23, 1099. (36) Hammes-Schiffer, S.; Tully, J. C. J. Chem. Phys. 1994, 101, 4657. (37) Martens, C. C.; Fang, J. J. Chem. Phys. 1997, 106, 4918. (38) Donoso, A.; Martens, C. C. J. Chem. Phys. 1998, 102, 4291. (39) Donoso, A.; Martens, C. C. J. Chem. Phys. 2000, 112, 3980. (40) Donoso, A.; Kohen, D.; Martens, C. C. J. Chem. Phys. 2000, 112, 7345. (41) Schu¨tte, C. Konard-Zuse-Zentrum fu¨r informationstechnik; Berlin, Germany, June 1999; Preprint pp SC 99-10. (42) Santer, M.; Manthe, U.; Stock, G. J. Chem. Phys. 2001, 114, 2001. (43) Kapral, R.; Ciccotti, G. J. Chem. Phys. 1999, 110, 8919. (44) Nielsen, S.; Kapral, R.; Ciccotti, G. J. Chem. Phys. 2000, 112, 6543. (45) Wan, C.; Schofield, J. J. Chem. Phys. 2000, 113, 7047. (46) Wan, C.; Schofield, J. J. Chem. Phys. 2000, 112, 4447. (47) Kapral, R. Annu. ReV. Phys. Chem. 2006, 57, 129. (48) Shi, Q.; Geva, E. J. Chem. Phys. 2004, 121, 3393. (49) Mukamel, S. Principles of Nonlinear Optical Spectroscopy; Oxford: New York, 1995. (50) Saven, J. G.; Skinner, J. L. J. Chem. Phys. 1993, 99, 4391. (51) Ka, B. J.; Geva, E. J. Chem. Phys. 2006, 125, 214501. (52) Stephens, M. D.; Saven, J. G.; Skinner, J. L. J. Chem. Phys. 1997, 106, 2129. (53) Khalil, M.; Demirodoven, N.; Tokmakoff, A. Phys. ReV. Lett. 2003, 90, 47401. (54) Evans, D. J.; Holian, B. L. J. Chem. Phys. 1985, 83, 4069. (55) MacKernan, D.; Kapral, R.; Ciccotti, G. J. Phys: Condens. Matter 2002, 14, 9069. (56) Marcus, R. A.; Sutin, N. Biochim. Biophys. Acta 1985, 811, 265– 322. (57) Warshel, A. J. Phys. Chem. 1982, 86, 2218–2224. (58) Navrotskaya, I.; Geva, E. J. Phys. Chem. A 2007, 111, 460.

JP8072816