Vibrational Energy Relaxation of a Hydrogen-Bonded Complex

4048

J. Phys. Chem. B 2008, 112, 4048-4058

Vibrational Energy Relaxation of a Hydrogen-Bonded Complex Dissolved in a Polar Liquid via the Mixed Quantum-Classical Liouville Method Gabriel Hanna and Eitan Geva* Department of Chemistry and FOCUS center, UniVersity of Michigan, Ann Arbor, Michigan 48109-1055 ReceiVed: August 1, 2007; In Final Form: January 8, 2008

The vibrational energy relaxation (VER) of the hydrogen stretch in a linear hydrogen-bonded complex dissolved in a polar solvent is studied. The study is based on the Azzouz-Borgis model [Azzouz, H.; Borgis, D. J. Chem. Phys. 1993, 98, 7361], which is known to account for many important features of real hydrogenbonded systems, including ionic-to-covalent tautomerism and a broad distribution of hydrogen stretch frequencies. A description of VER in this strongly coupled system is considered, which consists of the following three consecutive steps: (1) solvation on the adiabatic excited vibrational surface; (2) nonadiabatic transition from the excited to the ground adiabatic vibrational surface; and (3) solvation on the adiabatic ground vibrational surface. The relaxation dynamics during those three steps were simulated via the mixed quantum-classical LiouVille method, where the hydrogen is treated quantum-mechanically, while the other particles are treated in a classical-like manner. The solvation on the excited-state surface was found to occur rapidly (∼0.5 ps) and to involve energy exchange with both the intramolecular and intermolecular degrees of freedom. It was also found that, while energy is released to the solvent during the solvation of the covalent tautomer, the solvation of the ionic tautomer involves absorption of energy from the solvent. The decrease in the transition frequency during the solvation process also facilitates the nonadiabatic transitions, which occur rapidly (∼0.8 ps) thereafter. The nonadiabatic transitions were shown to be induced by interactions with a large number of solvent molecules and to be relatively insensitive to their location relative to the complex. Finally, solvation on the ground-state surface was seen to occur on a time scale of ∼1.0 ps and leads to nonequilibrium ionic and covalent subpopulations. Equilibration on the ground-state surface occurs on a significantly slower time scale (∼7.6 ps). Our results shed new light on the problem of VER in strongly coupled condensed phase systems that lie outside the range of validity of the Landau-Teller formula.

I. Introduction Vibrational energy relaxation (VER) is the process by which an excited vibrational mode releases its excess energy to a bath of intermolecular and/or intramolecular accepting modes. Virtually all chemical phenomena involve such VER processes, and much attention has been given over the last several decades to developing experimental techniques that can probe them on a wide range of time scales.1-39 The wealth of detailed experimental information on VER has also motivated many theoretical studies that attempted to provide a molecular interpretation for the observed time scales and pathways.40-74 Many of those theoretical studies have been based on the assumption that the relaxing mode is weakly coupled to a bath of accepting modes. Under those conditions, VER can be described in terms of a Markovian Master equation for the populations of the vibrational states, where the state-to-state rate constants are given in terms of the Fourier transform, at the vibrational transition frequency, of certain bath correlation functions. The validity of this approach, which will be referred to in the remainder of this paper as corresponding to the Landau-Teller (LT) regime,13,75,76 rests on several assumptions. First, the fact that the rate constants are given in terms of free-bath correlation functions implies that the bath Hamiltonian is not affected by the vibrational state of the relaxing mode. Second, the existence of a well-defined vibrational transition frequency and the fact that population * Corresponding author. E-mail: [email protected].

relaxation is not coupled to phase relaxation (vibrational dephasing) require that the coupling to the bath does not give rise to large fluctuations of the vibrational energy levels. Third, the fact that VER can be described in terms of a Markovian Master equation implies a separation of time scales such that VER occurs on a time scale that is significantly slower than that of the bath fluctuations inducing it. The validity of the assumptions underlying the LT regime becomes questionable in cases that involve strong coupling between the relaxing mode and the bath modes. Indeed, several recent studies have explored various routes for going beyond the LT regime by using time-dependent perturbation theory56,57,65,77-79 or generalizing the LT formula to cases involving extensive energy level fluctuations.57,78 However, those attempts were still based on a perturbative treatment of the coupling between the vibrational mode and the bath. It is therefore highly desirable to develop alternative nonperturbatiVe methods for simulating VER in such systems. Such methods can be based on describing VER in terms of nonadiabatic transitions between adiabatic vibrational states. This in turn suggests the use of a mixed quantum-classical treatment, where the relaxing mode is treated quantum-mechanically, while the bath modes are treated in a classical-like, trajectory-based manner. Indeed, mixed quantum-classical methods have been employed for studying VER in a few previous studies, including the use of the surface-hopping and mean-field methods for studying VER of the hydrogen stretch in a hydrogen-bonded

10.1021/jp076155b CCC: $40.75 © 2008 American Chemical Society Published on Web 03/11/2008

VER of Hydrogen-Bonded Complex in a Polar Solvent complex in polar liquid solution80,81 and of dihalogens in noble gas liquids63,64,66,82 and the use of the mean-field method for studying VER of CN- and NO in water.65,67,68 In this paper, we study the VER of the hydrogen stretch in a moderately strong hydrogen-bonded complex dissolved in a polar liquid. This model system represents an important class of systems for which the assumptions underlying the LT regime are not valid. Our approach to studying VER in this system is based on the recently proposed mixed quantum-classical LiouVille (MQCL) method.83-96 Within this approach, VER is described as a nonequilibrium process involving solvation on different vibrational adiabatic potential energy surfaces and nonadiabatic transitions between them (see Section III). Importantly, this approach does not rely on a perturbative treatment of the interaction between the relaxing mode and the bath modes and allows for a direct analysis of the molecular mechanism underlying VER. The remainder of this paper is organized as follows: The model of a hydrogen-bonded complex in a dipolar solution, which was developed by Azzouz and Borgis,97 is described in Section II. The theoretical framework for treating VER as a nonadiabatic process is outlined in Section III. The MQCL method and its application to VER are described in Section IV. The simulation techniques are described in Section V. The results are presented and discussed in Section VI. Finally, the main conclusions are summarized in Section VII. II. Model The study presented in this paper is based on the AzzouzBorgis (AB) model of a hydrogen-bonded phenol-trimethylamine complex dissolved in liquid methyl chloride.97 In this model, the hydrogen is assumed to move along a 1D donoracceptor axis, where the proton donor, A, and acceptor, B, are modeled as single particles and parametrized to represent phenol and trimethylamine, respectively.97 The functional form of the intramolecular potential surface as a function of the A-B and A-H distances is adopted from ref 98. Polarizability effects are accounted for by assuming that the charges on A and B are explicitly dependent on the position of the proton. Within this model, the strength of the hydrogen bond increases as the equilibrium A-B distance decreases. More specifically, increasing the equilibrium A-B distance modifies the potential along the A-H coordinate from a single well at short A-B distances to a double-well at large A-B distances. The doublewell form of the potential gives rise to tautomeric equilibrium between a covalent form and an ionic form of the complex: AH - B h A- - H+B. In the covalent form, the proton is attached to the phenol (A) and the complex has a small (2.5 D) dipole moment due to the polarity of the A-H bond. In the ionic form, the proton is attached to trimethylamine (B), and the resulting charge separation gives rise to a significantly larger (10.5 D) dipole moment. Finally, it should be noted that the AB model provides a somewhat oversimplified description of hydrogen-bonded complexes (see Section VII). In our case, the equilibrium A-B distance was set to 2.7 Å, which corresponds to a double-well profile along the A-H coordinate with a relatively low barrier (i.e., a moderately strong hydrogen bond). 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 were adopted from ref 99. The covalent form, AHB, is favored in the absence of the solvent. However, the ionic

J. Phys. Chem. B, Vol. 112, No. 13, 2008 4049 form, A--H+B, becomes somewhat more stable when the complex is dissolved in a polar solvent. The AHB complex is analogous to a (rather floppy) asymmetric linear triatomic molecule. The large mass disparity between the hydrogen and the A and B groups implies that the “symmetric” normal mode can be viewed as a vibration of the A-B distance, while the “asymmetric” normal mode can be described in terms of mostly the motion of the proton relative to the essentially fixed A and B groups. A unique feature of VER in polyatomic molecules is that it can occur via either an intramolecular pathway (i.e., by energy transfer between the normal modes) or an intermolecular pathway (i.e., by releasing the excess energy to the solvent).40,43,44,51,52,58,72,73 This should be contrasted with the case of diatomic molecules where VER must follow an intermolecular pathway. The present study therefore also provides an opportunity for exploring the VER pathway in a strongly coupled polyatomic molecule. Previous studies of the AB model have been focused on calculating the rate constant for proton transfer.97,99-107 Studies that explored the VER of the hydrogen stretch in the AB model have been previously reported in refs 48, 80, and 81. The study in ref 48 was based on a classical treatment of VER within the framework of the LT regime, and assumed that the hydrogen is fixed at its equilibrium position in the covalent form and that the vibrational transition between the ground and first excited vibrational states corresponds to a sharp resonance. This study also did not account for the effect of intramolecular energy exchange between the hydrogen stretch and A-B stretch. Those assumptions are rather restrictive, as noted by the authors of ref 48. More specifically, the hydrogen position is rather delocalized, and the hydrogen stretch transition frequency fluctuates rapidly over a very broad range of frequencies, which actually makes the AB model suitable for studying VER in a system that is outside the LT regime. Furthermore, we will show that the VER of the hydrogen stretch involves energy exchange with the A-B stretch. The studies in refs 80 and 81 were based on treating the hydrogen stretch quantum-mechanically. In ref 80, VER was treated as a nonadiabatic nonequilibrium process within the framework of the surface-hopping method, which is distinctly different from the MQCL method employed here. The analysis in ref 81 included energy transfer to a quantized A-B stretch. However, the calculation of the VER rate was based on assuming that the energy gap between the vibrational states is small and employed a Fermi’s golden rule-type expression, which assumes weak coupling between the relaxing mode and accepting modes. III. Vibrational Energy Relaxation The light mass of the hydrogen relative to the other particles suggests that a mixed quantum-classical treatment, where the hydrogen is treated quantum-mechanically while the solvent and the A and B groups are treated in a classical-like manner, may be appropriate. The mixed quantum-classical Hamiltonian of the overall system is given by

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

(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) from now on), qˆ and pˆ are the coordinate and momentum operators of the proton, KB(P) is the kinetic energy of the bath DOF, K ˆ P(pˆ ) is

4050 J. Phys. Chem. B, Vol. 112, No. 13, 2008

Hanna and Geva

Figure 1. The ground and first excited FESs, in units of kBT, as a function of the solvent polarization, ∆E (see eq 3), for different constrained values of the A-B distance: 2.60 Å (red), 2.63 Å (green), 2.70 Å (blue), 2.76 Å (magenta). Also shown is the FES for an unconstrained A-B distance (black).

the kinetic energy operator of the proton, and Vˆ (Q,qˆ ) is the overall potential energy operator. The mass scale separation implies that one may assume that the vibrational energy levels of the hydrogen follow the bath DOF in an adiabatic manner. It is therefore reasonable to identify the vibrational energy levels of the hydrogen stretch with the eigenvalues of the protonic Hamiltonian, H ˆ P(qˆ , pˆ ,Q) ≡K ˆ P(pˆ ) + Vˆ (Q,qˆ ), where

H ˆ P|j;Q〉 ) Ej(Q)|j;Q〉

(2)

Here, {|j;Q〉} are the adiabatic vibrational states, and {Ej(Q)} are the corresponding vibrational energy levels (both explicitly dependent on the coordinates of the bath DOF). It should also be noted that Ej(Q) corresponds to the multidimensional potential energy surface, which dictates the adiabatic dynamics of the bath DOF in the jth vibrational state. The dependence of Ej(Q) on Q implies that the vibrational energy levels will fluctuate in response to variations in the positions of the bath DOF. The extent of those fluctuations can be inferred from Figure 1, which shows the ground and first excited free energy surfaces (FESs) as a function of the solvent polarization, for different constrained values of the A-B distance. The solvent polarization is defined by108,109

∆E(Q) )

zi,ae ∑ i,a

(

|Qai

1 - s|

1

|Qai

- s′|

)

(3)

and corresponds to the difference between the solvent electrical potentials at the minima of the covalent and ionic wells. More specifically, zi,ae (e ) 1.602 × 10-19 C) is the charge on atom a in solvent molecule i, and s and s′ correspond to the center of mass of the complex, and a point displaced by 0.56 Å relative to it toward the B group, respectively.100 Thus, the stability of the covalent form relative to the ionic form increases with decreasing ∆E(Q). The two FESs shown in Figure 1 are clearly very different, with the ground state FES having a double-well form, while the excited state FES has a single-well form. This implies that the configurations sampled by the solvent in the ground state will differ considerably from those sampled in the excited state. While the proton is in the ground state, the solvent will favor configurations that correspond to either low polarizations or high polarizations, with a low likelihood for configurations with intermediate values of the polarization. This behavior is a signature of the above-mentioned tautomeric equilibrium be-

Figure 2. The distributions of transition frequencies between the ground and first excited adiabatic vibrational states, as obtained by equilibrium sampling on the ground (blue line) and first excited (red line) potential energy surfaces.

tween the covalent form of the complex, which favors lowpolarization configurations, and the ionic form of the complex, which favors high-polarization configurations. It should be noted that significant subpopulations of both forms coexist at equilibrium on the ground vibrational state potential energy surface (the equilibrium composition is 65% ionic and 35% covalent). It has also been shown that proton transfer, which converts the complex between the covalent and ionic forms, is driven by the very same fluctuations in the solvent polarization (for this model, the covalent-to-ionic and ionic-to-covalent free energy barrier heights are given by 2.09 kBT and 2.88 kBT, respectively).100,110 The situation is rather different when the proton occupies the excited state. In this case, the very same configurations that were the least favorable on the ground state, that is, configurations that correspond to intermediate values of the polarization, become the most favorable ones. It should be noted that this observation implies that the potential energy surface for the bath DOF is strongly dependent on the vibrational state of the complex, which is clearly inconsistent with one of the main assumptions underlying the LT regime, namely, that the bath Hamiltonian is independent of the vibrational state. The fluctuations in the vibrational energy levels also modulate the vibrational transition frequency. In Figure 2, we show the distributions of transition frequencies between the ground and first excited adiabatic vibrational states, as obtained by equilibrium sampling on the ground and first excited potential energy surfaces. The distributions are seen to be very different from each other, which reflects the fact that the configurations sampled by the solvent are strongly dependent on the vibrational state. The two distributions are also seen to be very broad (∼2500 cm-1 on the ground state potential surface and ∼1000 cm-1 on the excited-state potential surface), thereby making it difficult to assign a single value to the transition frequency. The latter observation is inconsistent with another major assumption underlying the LT regime, namely, that one can assign a single value to the transition frequency (which is used as input to the LT formula). Within the LT regime, VER can be described in terms of a Markovian Master equation involving rate constants that are independent of the initial vibrational states. This may not necessarily be the case beyond the LT regime, where the time scale and pathway of VER may depend on the initial state. In this paper, we consider a scenario, originally proposed by Kim et al. in ref 81, which is likely to be realized experimentally in a pump-probe experiment. We assume that the system is initially in the ground vibrational state, so that the corresponding ensemble of solvent configurations consists of two subsets of

VER of Hydrogen-Bonded Complex in a Polar Solvent

J. Phys. Chem. B, Vol. 112, No. 13, 2008 4051

Bˆ W(Q,P,t) )

∫ dZ exp[iZ‚P/p]〈Q - Z/2|Bˆ (t)|Q + Z/2〉

(4)

where Z ) (Z1, Z2, ...). The operator Bˆ W(Q,P,t) can be conveniently represented by a matrix in terms of the adiabatic basis set, {|j;Q〉}, where the corresponding matrix elements are defined by BW ˆ W(Q,P,t)|j;Q〉. The matrix eleij (Q,P,t) ) 〈i;Q|B ments of the partially Wigner transformed density operator are ˆ W(Q,P,t)|j;Q〉, where similarly given by FW ij (Q,P,t) ) 〈i;Q|F Figure 3. A schematic view of VER in the AB model.

configurations, which are characterized by either high or low polarizations. We assume that, at time t ) 0, the complex is photoexcited to the first excited vibrational state by a broad band pulse. The VER that follows photoexcitation then involves three distinctly different relaxation processes (see Figure 3):81 • Excited-state solVation, that is, rearrangement of the solvent configuration toward a configuration that is favored by the excited vibrational state. This nonequilibrium process has no analogue in the LT regime, where it is assumed that the vibrational state does not affect the ensemble of configurations sampled by the bath. The solvation process can be probed experimentally, for example, by measuring the corresponding dynamical Stokes shift.60,111-116 Interestingly, although the solvation process involves energy exchange between the proton and the bath DOF, it does not involve population transfer, since the proton remains on the adiabatic excited-state potential surface throughout the process. • Vibrational population relaxation, that is, the transition of the system between the excited and ground adiabatic vibrational states. The process of vibrational population relaxation is driven by the nonadiabatic coupling terms between those two states. It should be noted that, within the LT regime, one calculates the rate constant for this nonadiabatic transition via Fermi’s golden rule in the limit where the ground and excited adiabatic surfaces have the exact same form. • Ground state solVation, that is, rearrangement of the solvent configuration toward a configuration that is favored by the ground vibrational state. As for the excited-state solvation, this nonequilibrium process lacks an analogue in the LT regime and involves energy exchange without population relaxation. IV. Mixed Quantum-Classical Liouville Dynamics As discussed in Section III, VER in the AB model involves nonadiabatic transitions. In this paper, we will describe the underlying nonadiabatic dynamics within the framework of the MQCL method.83-92,94,95 This approach is based on formulating the full quantum dynamics within the partial Wigner representation (i.e., Wigner transforming with respect to the heavy DOF) and taking the limit in which the mass ratio between the light and heavy particles becomes small. The resulting equation of motion describes the time evolution of a quantum subsystem of light particles coupled to a classical subsystem of heavy particles. More specifically, consider an overall system observable that is represented by the operator Bˆ . Its time evolution in the Heisenberg picture is given by Bˆ (t) ) exp(iH ˆ t/p)Bˆ exp(-iH ˆ t/ p), where H ˆ is the overall Hamiltonian. The partial Wigner transform of Bˆ (t) over the heavy DOF is then given by

Fˆ W(Q,P,t) ) 1 (2πp)N

∫ dZ exp[iZ‚P/p]〈Q - Z/2|Fˆ (t)|Q + Z/2〉

(5)

Thus, the expectation value of the observable at time t is given by

〈Bˆ (t)〉 ) Tr{Fˆ (0)Bˆ (t)} )

∑ij ∫ dQ ∫ dPFWij (Q,P,0)BWji (Q,P,t)

(6)

The time evolution of BW ij (Q, P, t) as dictated by the MQCL equation in the adiabatic representation is given by85,88-92,94 ∂ ∂t

{

BW ij (Q,P,t) )

U‚∇Q +

}

Fi(Q) + Fj(Q) W ‚∇P BW ij (Q,P,t) + iωij(Q)Bij (Q,P,t) + 2

[

1

∑ U‚d (Q) 1 + 2 S ‚∇ / jk

k

/ jk

[

P

]

BW ik (Q,P,t) + 1

∑ U‚d (Q) 1 + 2 S ‚∇ ik

k

ik

P

]

BW kj (Q,P,t) (7)

Here, U ) (U1, U2, ...) ) (P1/M1, P2/M2, ...) are the velocities of the heavy 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). It should be noted that, within the MQCL method, one describes the time evolution of an observable in terms of the time evolution of n2 phase-space variables, {BW ij (Q,P,t)|i, j ) 1, ..., n}, where n is the dimension of the quantum subsystem’s Hilbert space. The first term on the right-hand side of eq 7 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 7 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). The MQCL equation has several features that make its application to vibrational relaxation particularly attractive, namely, (1) it allows for energy exchange between the “quantum” and “classical” DOF; (2) it satisfies energy conservation in the sense that the expectation value of the Hamiltonian in eq


Hanna and Geva

1 is a constant of the motion; (3) its stationary solution coincides with the exact quantum mechanical equilibrium state up to first order in p;117 and (4) it treats the vibrational populations and coherences on a similar footing, thereby automatically accounting for the coupling between them. The model system considered here can be described in terms of the ground and first excited vibrational states, {|1;Q〉, |2;Q〉}. The dynamics can therefore be described in terms of four phase W W space variables, namely BW 11(Q,P,t), B12(Q,P,t), B21(Q,P,t) and W B22(Q,P,t). Our particular choice of initial state implies that W FW ij (Q,P,0) ) δ(i, 2)δ(j, 2)F22(Q,P,0)

(8)

with FW 22(Q,P,0) corresponding to equilibrium on the ground state potential surface. We will also assume that FW 22(Q,P,0) can be approximated by its classical limit:

FW 22(Q,P,0) ≈

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

∫ dQ ∫ dQ exp{-β[KB(P) + E1(Q)]}

(9)

where β ) 1/kBT. The observable of interest is the excitedstate population, which is represented by the operator Pˆ ) |2;Q〉〈2;Q|. Its expectation value at time t is obtained by substituting FW ij (Q,P,0) from eq 8 into eq 6:

P2(t) ≡ 〈Pˆ 〉(t) )

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

(10)

Here, PW 22(Q,P,t) is obtained by solving eq 7 with the initial conditions PW ij (Q,P,0) ) δ(i,2)δ(j,2).

trajectories of lengths 0.5-3.0 ps (depending on the time scale of the process one wishes to simulate). Simulations of nonadiabatic dynamics were based on numerically solving the MQCL equation for the excited-state population. The numerical solution of eq 7 was carried out via the sequential short-time propagation algorithm and employed the momentum jump approximation.100,121 This algorithm is based on averaging over a large number of “nonadiabatic trajectories” (the overall number is dictated by convergence). The following quantities are monitored during each MQCL trajectory: RjR′j, Qj, Pj, and Fj. They correspond to the potential energy surface (i.e., [ERj(Q) + ER ′j (Q)]/2, values of the coordinates and momenta of the bath DOF and a cumulative factor Fj (see below), respectively, at the jth time step. At the initial time, R0R′0 ) 22, {Q0, P0} are sampled from the classical ground state equilibrium ensemble, and F0 ) 1.0. Propagation by one time step ∆t from {RjR′j, Qj, Pj, Fj} to {Rj+1R′j+1, Qj+1, Pj+1, Fj+1} involves the following steps: • Propagate the system classically for one time step ∆t on the potential energy surface [ERj(Q) + ER′j(Q)]/2, starting at Qj, Pj and ending at Qj+1, Pj+1. • Calculate the corresponding phase factor WRjR′j ) exp[ iωRjRj′∆t]. • Choose with a probability of 1/2 whether the attempted nonadiabatic transition would involve changing Rj into β or R′j into β. For example, if RjR′j ) 12, then we can either change Rj into β ) 2 and make the transition to βR′j ) 22 or change R′j into β ) 1 and make the transition to Rjβ ) 11. • Numerically calculate the nonadiabatic coefficient for the attempted nonadiabatic transition (either U‚dβ,R′j(Qj+1) or U‚dRjβ(Qj+1)). • Define the probability for making a nonadiabatic transition:

V. Simulation Techniques The results reported in this paper are based on simulations similar to those described in ref 100. We provide a brief outline of those techniques in this section for the sake of completeness (a more detailed discussion is provided in ref 100). Simulations were performed with a single AHB 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 bond length of the methyl chloride molecules was constrained at 1.781 Å using the SHAKE118 and RATTLE119 algorithms. 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 Hamiltonian for different instantaneous bath configurations. To this end, we first constructed the protonic Hamiltonian matrix in terms of a basis that consisted of 12 harmonic oscillator basis functions with the first set of 6 basis functions centered at the ionic well of the gas-phase intramolecular potential surface and the other set of 6 basis functions centered at the covalent well. This matrix was diagonalized by using the LAPACK DSYGV routine. Adiabatic dynamics was simulated by solving the corresponding classical equations of motion via the velocity Verlet algorithm with a time step of 1 fs. The initial conditions in all the simulations reported below correspond to equilibrium on either the ground state or excited-state potential energy surface. Equilibrium configurations were generated via a constant temperature molecular dynamics simulation, using a Nose´Hoover thermostat.120 The results reported in the following section were obtained by averaging over 16 000-20 000

Π(RjR′j f βR′j) )

|Uj+1‚dβR′j(Qj+1)|∆t 1 + |Uj+1‚dβR′j(Qj+1)|∆t

(11)

or

Π(RjR′j f Rjβ) )

|Uj+1‚dRjβ(Qj+1)|∆t 1 + |Uj+1‚dRjβ(Qj+1)|∆t

(12)

and use it to sample the transitions. • In the case when the transition is rejected, Rj+1R′j+1) RjR′j, modify the cumulative factor as follows:

Fj+1 ) Fj × WRjR′j ×

1 1-Π

(13)

In the case where the transition is accepted, Rj+1R′j+1 ) βR′j or Rj+1R′j+1 ) Rjβ, perform a momentum jump, Pj+1 f Pj+1 + ∆P (see ref 100 for an explicit expression for ∆P) and modify the cumulative factor to

Fj+1 ) 2Fj × WRjR′j ×

1 × (Uj+1‚dβR′j)∆t Π

(14)

Fj+1 ) 2Fj × WRjR′j ×

1 × (Uj+1‚dRjβ)∆t Π

(15)

or

respectively. • Calculate the excited-state population at the jth time step by averaging over a large number of trajectories:



P2(j) ) sum of Fj from trajectories where (RjR′j) ) (22) (16) overall number of trajectories VI. Results and Discussion As discussed in Section III, VER in the AB model involves adiabatic solvation and nonadiabatic population relaxation, which are distinctly different processes. A nonadiabatic transition can take place at any point during solvation on the excitedstate surface. Nevertheless, the fact that red-shifting the transition frequency increases the probability for a nonadiabatic transition suggests that population relaxation will follow the excited-state solvation in a sequential manner. We will therefore first consider the solvation process that follows impulsive photoexcitation from the ground vibrational state to the first excited vibrational state. To this end, we simulated the adiabatic nonequilibrium solvation dynamics without allowing for nonadiabatic transitions to take place. The ensemble-averaged vibrational energy of the adiabatic excited-state during solvation is shown in Figure 4. It is seen to exhibit a damped oscillatory behavior. The ∼80 fs period of the oscillations is suggestive of coupling with the lowfrequency symmetric stretch of the complex. Indeed, the oscillations disappear upon constraining the A-B distance to its equilibrium value (see upper panel of Figure 4), and can therefore be attributed to intramolecular energy exchange between the hydrogen stretch and the A-B stretch. The damping can be attributed to energy loss to the solvent, which is seen to occur on a time scale of ∼0.5 ps. It should be noted that the damping is somewhat faster in the case of the vibrating complex, which can be attributed to a VER pathway that involves energy transfer to the A-B stretch. Further insight into the excited-state solvation process can be obtained by considering the solvation of the ionic and covalent forms separately (see lower panel of Figure 4). Interestingly, the solvation of those two forms follow very different trends. More specifically, while the solvation of the covalent form exhibits a damped oscillatory behavior, that of the ionic form involves a more or less monotonic increase of the vibrational energy. Those opposite trends can be attributed to the corresponding strength of the polar solute-solvent interactions. In the case of the covalent form, one starts out with a solute that has a relatively small dipole that increases during solvation. Thus, solvation is accompanied by stabilization and energy release to the solvent. Furthermore, the fact that the dipole is rather small leads to weak coupling to the solvent, which explains the more pronounced energy exchange with the intramolecular A-B stretch. In contrast, the solute in the ionic form starts out with a relatively large dipole, which decreases during solvation, thereby leading to destabilization and the absorption of energy from the solvent. Finally, the fact that the dipole is rather large in this case leads to strong coupling to the solvent, which explains the lack of energy exchange with the intramolecular A-B stretch. One way of probing the solvation is by measuring the dynamical Stokes shift, that is, by monitoring the emission spectrum at different points along the solvation process. As shown in Figure 5, the (inhomogeneously broadened) emission line shape, which coincides with the (inhomogeneously broadened) absorption line shape right after photoexcitation, relaxes into the equilibrium emission spectrum within ∼0.5 ps. The lineshapes remain very broad throughout the solvation process. However, a damped oscillation originating from the abovementioned intramolecular energy exchange can be detected in

Figure 4. Upper panel: The (normalized) ensemble-averaged vibrational energy of the adiabatic excited-state as a function of time during solvation. Shown are the results in the case of a flexible A-B distance (black line) and in the case where the A-B distance is constrained to its equilibrium distance (blue line). Lower panel: The contributions of the ionic (green line) and covalent (red line) subpopulations to solvation on the excited-state potential surface.

Figure 5. The inhomogeneously broadened time-resolved emission line shape during solvation on the excited-state potential energy surface (dynamical Stokes shift).

Figure 6. The average emission frequency during solvation on the excited-state potential energy surface.

the time dependence of the average vibrational frequency during the solvation process (see Figure 6). We next turn to consider the nonadiabatic transitions between the excited and ground vibrational states. As expected, the strength of the nonadiabatic coupling between the ground and excited vibrational states is strongly dependent on the transition


Hanna and Geva

Figure 7. A scatter plot of the correlation between the transition frequency and the strength of nonadiabatic coupling as measured by |U‚d12|.

frequency (see Figure 7). More specifically, large nonadiabatic coupling coefficients are only observed at transition frequencies below ∼200 cm-1, which correspond to solvent polarizations between 0.0115 and 0.0165eC/Å. We therefore only account for nonadiabatic transitions within this range of polarizations. We also restrict the maximum number of nonadiabatic transitions per trajectory to eight, and filter the value of Fj along a trajectory so that when it exceeds some threshold value, Fj is set equal to the threshold value.100 This prevents the excessive accumulation of the Monte Carlo weight factors at long times, which magnifies the sign problem associated with averaging an oscillatory integrand. We have carried out tests to confirm that relaxing those restrictions does not alter the results significantly. In principle, nonadiabatic transitions between the vibrational states can take place at any point along the excited-state solvation process. However, the fact that the transition frequency decreases during the process of solvation suggests that the nonadiabatic transitions will mainly occur after the solvation process is completed and thermal equilibrium on the excitedstate potential surface is established. The red curve in the upper panel of Figure 8 shows the relaxation of the excited-state population when the initial ensemble of bath configurations corresponds to equilibrium on the excited-state potential. The relaxation in this case is seen to occur on a time scale of ∼0.8 ps, and is distinctly nonexponential, which is suggestive of a non-Markovian behavior. For comparison, the black curve in the upper panel of Figure 8 shows the relaxation of the excitedstate population when the initial ensemble of bath configurations corresponds to equilibrium on the ground-state potential. The relaxation is still nonexponential in this case, and is seen to occur on a somewhat slower time scale of ∼1.0 ps, which can be attributed to the fact that the nonadiabatic relaxation process is preceded by solvation on the excited-state potential. This also implies that the nonadiabatic transition rate is insensitive to whether the complex is in its covalent form or ionic form when initially excited. The strength of the nonadiabatic coupling is also seen to be strongly dependent on the A - B distance within the range of distances between 2.6 and 2.8 Å sampled by the complex (see Figure 9). Indeed, a significant acceleration in the nonadiabatic relaxation rate is observed when the A-B distance is constrained to its equilibrium value, 2.7 Å, where the strength of nonadiabatic coupling is seen to peak (see lower panel of Figure 8). In contrast, the nonadiabatic relaxation rate is seen to slow down by more than a factor of 2 when the A-B distance is constrained to a value of 2.6 Å for which the nonadiabatic coupling strength is significantly weaker (see lower panel of Figure 8). This observation is also consistent with the significant increase in

Figure 8. Upper panel: The nonadiabatic relaxation of the excitedstate population, as obtained via the MQCL method when the initial ensemble of bath configurations corresponds to equilibrium on the excited-state potential surface (black line), and to equilibrium on the ground-state potential surface (red line). Lower panel: The dependence of the rate of nonadiabatic relaxation on the A-B distance. The initial ensemble of bath configurations corresponds to equilibrium on the ground-state potential surface. Shown are the relaxation rates obtained in the case of a flexible A-B distance (black line), and in the cases where the A-B distance is constrained to 2.6 Å (red line) and to 2.7 Å (green line).

Figure 9. A scatter plot of the correlation between the nonadiabatic coupling coefficient and the A-B distance.

the free energy gap between the ground and excited states in going from an A-B distance of 2.7 Å to 2.6 Å (see Figure 1). Insight into the mechanism underlying the nonadiabatic population relaxation can be obtained by considering the number of particles that significantly contribute to the nonadiabatic coupling and their positions relative to the proton. In Figure 10 we show the average top 20 largest contributions to the nonadiabatic coupling coefficients whose value is above a threshold of 20 ps-1, as obtained from a 75 ps adiabatic equilibrium simulation on the excited-state potential surface. As can be seen, the nonadiabatic transitions are induced by simultaneous interactions with a relatively large number of particles. It should be noted that a similar analysis in the case of iodine in liquid xenon by Li and Thompson64 led those authors to the distinctly different conclusion that the nonadiabatic transition is induced by a single solvent atom, which is also consistent with the instantaneous-pair theory of Larsen and


Figure 10. The average top 20 largest contributions to the nonadiabatic coupling, above a threshold of 20 ps-1, as obtained from a 75 ps adiabatic equilibrium simulation on the excited-state potential surface. The contributions were normalized so that the largest one is equal to one.

Figure 11. The position distribution of bath particles with the largest contributions to the nonadiabatic coupling. The distribution is given as a function of R, the distance between the bath particle and the mean position of the proton, and θ, the angle between the vector that points from the mean proton position to the bath particle and the axis of the complex. Each point corresponds to a single time step during a 75 ps equilibrium adiabatic simulation on the excited-state surface, where at least one nonadiabatic coefficient was found to be above a threshold of 20 ps-1.

Stratt.122 The different conclusions can be attributed to the fact that the analysis in ref 64 was performed on a nonpolar solution, which is characterized by weak short-range solvent-solute interactions, whereas the analysis presented in this paper was performed on a system that involves a dipolar solute dissolved in a polar liquid, which is characterized by stronger long-range interactions. It should also be noted that the AHB complex is rather large, and therefore has more solvent molecules in its first solvation shell. In Figure 11, we show the position distribution of the particles that make the largest contributions to the nonadiabatic coupling coefficient, as obtained from the same 75 ps adiabatic equilibrium simulation on the excited-state potential surface. The two peaks at cos(θ) ) (1 correspond to contributions from the A and B groups on both sides of the proton, while the continuous semicircular distribution corresponds to contributions from solvent molecules. We have also observed similar distributions in the case of particles that make smaller, yet significant contributions to the nonadiabatic transitions. The results show that both the A and B groups and the solvent molecules contribute to the nonadiabatic transitions. Furthermore, the solvent molecules that contribute the most to the nonadiabatic transitions are distributed rather uniformly around the complex. It is interesting to contrast this conclusion with that obtained from a similar analysis in the case of iodine in xenon,64 where the main contributions to the nonadiabatic transitions came from


Figure 12. The ensemble-averaged vibrational energy of the adiabatic ground state during solvation (blue line). Also shown are the corresponding contributions from trajectories that end in the covalent well (red line) and ionic well (green line).

solvent atoms arranged at the ends or near the center of the iodine molecule. The different conclusions can be explained as another manifestation of the above-mentioned fundamental difference between nonpolar and polar solutions. Finally, we consider the solvation on the ground state potential surface that follows the nonadiabatic transition. To this end, we simulated the adiabatic nonequilibrium solvation dynamics on the ground-state potential surface with the initial conditions as obtained immediately after the nonadiabatic transition takes place. For convenience, we assumed that all the solvation trajectories start at the same time, although in practice they obviously start at random times depending on when the nonadiabatic transition occurs. The relaxation of the average vibrational energy of the adiabatic ground state during solvation is shown in Figure 12. It exhibits a damped oscillatory behavior, as in the case of solvation on the excited-state potential surface. The oscillations can be attributed to intramolecular energy exchange between the hydrogen stretch and the A-B stretch, while the damping reflects energy exchange with the solvent, which is seen to occur on a time scale of ∼1.0 ps. It should be noted that the energy does not reach a true plateau at the end of this process, since the system is not in a true equilibrium state (the inverse rate constant for interconversion between the ionic and covalent forms is 7.6 ps123 and therefore significantly slower). Further insight into the ground state solvation process can be obtained by comparing the trajectories that end in the covalent well to those that end in the ionic well (see Figure 12). The vibrational energy is seen to increase along the trajectories that end in the covalent well, and to decrease along the trajectories that end in the ionic well. Those trends are consistent with those observed during the excited-state solvation, and can be explained in a similar manner. Another quantity of interest is the ionic/covalent branching ratio at the end of the ground state solvation process. It should be noted that the time scale for interchanging between the ionic and covalent forms is 7.6 ps, which is an order of magnitude slower than the time scale of solvation. Thus, the branching ratio at the end of the solvation process does not have to coincide with its equilibrium value. Indeed, we have found the ionic/ covalent branching ratio to be 0.9 at the end of the solvation process, which is smaller by about a factor of 2 relative to the corresponding equilibrium branching ratio of 1.9.123 Thus, the system can be prepared in a nonequilibrium state that favors the covalent form over the ionic form via the above-mentioned procedure. The above-mentioned nonequilibrium branching ratio is closely related to the phenomenon of vibrationally induced

4056 J. Phys. Chem. B, Vol. 112, No. 13, 2008 proton transfer as described in refs 80 and 81, where it was observed that vibrational excitation gives rise to a somewhat faster interconversion rate between the tautomers. It should be noted that the rate enhancement is diminished by the excitedstate solvation process. Thus, one expects a larger rate enhancement in the case of the subpopulation that results from photoexcitation in the red edge of the hydrogen stretch absorption band, since no solvation needs to take place prior to the nonadiabatic transition.80 VII. Concluding Remarks The AB model represents an important class of strongly coupled systems where the validity of the assumptions underlying the LT regime become questionable. For such systems, the coupling between the relaxing mode and bath of accepting modes may not be weak enough to justify treating it as a small perturbation, there is no clear separation of time scales between VER and the bath fluctuations that induce it, and the transition frequency is ambiguous because of large fluctuations of the vibrational energy levels. As a result, energy and phase relaxation cannot be assumed to be decoupled, and the bath Hamiltonian is strongly dependent on the vibrational state of the relaxing mode. In this paper, we considered an alternative description of VER in such systems in terms of solvation and nonadiabatic transitions. This study involves several new features, which set it apart from previous studies that employed a similar description of VER:63-68,82,81 • To the best of our knowledge, the present study represents the first application of the MQCL method to VER. • To the best of our knowledge, there have been only three previous studies of VER in the AB model, namely, refs 48, 80, and 81. The study in ref 48 was based on a classical treatment of VER within the framework of the LT regime, and assumed that the hydrogen is fixed at its equilibrium position in the covalent form and that the vibrational transition between the ground and first excited vibrational states corresponds to a sharp resonance. Within the framework of those approximations, the authors of ref 48 predicted that VER of the hydrogen stretch occurs on the nanosecond time scale. The discrepancy between this prediction and the one reported in the present paper appears to be consistent with the expectation that this system lies beyond the range of validity of the LT regime. The studies in refs 80 and 81 are similar to the one considered here in the sense that they treat VER as a nonadiabatic nonequilibrium process. However, this is done either within the framework of the surface-hopping method80 or via a golden-rule-based approach,81 which are distinctly different from the MQCL method employed here. Furthermore, the analysis in those previous studies was aimed at elucidating the effect of VER on the rate of proton transfer and, as a result, did not address many aspects of the problem that have been addressed in the present paper. • Our results suggest that VER in this model system clearly lies beyond the LT regime. This should be contrasted with the systems studied in refs 63-68 and 82, where the predictions of the LT approach were found to be in good agreement with those obtained via the nonperturbative treatment. • Most previous studies were focused on VER of diatomic molecules, which can only occur via a single intermolecular pathway, namely, energy exchange with the solvent. In contrast, the hydrogen-bonded complex studied here can be viewed as a simple model of a polyatomic molecule, where VER can occur via either an intramolecular pathway (i.e., energy exchange with the A-B stretch) or an intermolecular pathway (i.e., energy

Hanna and Geva exchange with the polar solvent). As we have shown, both solvation and nonadiabatic dynamics in this system are strongly affected by the intramolecular DOF. Furthermore, the fact that the ionic and covalent tautomers are stable on the time scale of VER implies that the VER pathway can be tautomer-dependent. Indeed, we have observed distinctively different solvation pathways for the different tautomers. The work reported in this paper can be extended in several directions. For example, one may consider the rates of VER between other vibrational states, the effect of isotope substitution (from hydrogen to deuterium), and the dependence of the VER rate on density and temperature. Several limitations of the current study should also be pointed out. First, it should be noted that, despite its richness, the AB model provides a somewhat oversimplified description of hydrogen-bonded complexes. For example, while the AB model allows one to consider VER of the hydrogen stretch via energy transfer to the A-B stretch and the rotational and translational solvent DOF, other potentially important VER pathways that involve energy exchange with the proton bending mode and the vibrational modes of the donor, acceptor, and solvent molecules are left out. Consideration of those pathways can be made possible by incorporating a description of the proton motion in 3D124 and including an explicit description of the internal structure of the donor, acceptor, and solvent molecules.125 However, it should be noted that the pathway considered here is rather efficient and leads to rapid VER on the time scale of a few picoseconds. Thus, this particular pathway may remain significant, even if other pathways are accounted for. In this study, the polarizability of the complex is simply accounted for through the use of an empirical charge switching function that depends explicitly on the position of the proton, but not on the A-B distance, as in ref 97. We would like to point out that a more elaborate, valence-bond-based treatment of this polarizability was reported in ref 101, which accounts for both the protonic position and the A-B distance dependence. Moreover, the fact that the solute-solvent interactions are described in terms of nonpolarizable force fields represents another notable limitation of the AB model. Shifting to polarizable force fields126 could significantly modify both the structure and dynamics of the liquid solution and therefore its effect on the solvation and nonadiabatic transition rates would be difficult to predict. Another simplification has to do with the fact that only the hydrogen stretch mode is treated quantum mechanically. Our choice to restrict the quantum-mechanical treatment to this highfrequency mode was motivated by the fact that this degree of freedom is expected to exhibit the most pronounced quantum effects. Quantizing the A-B stretch81 and avoiding the approximation of replacing FW 22(Q,P,0) by its classical limit (see eq 9)69-73 represent natural extensions of the present work, although one expects them to lead to smaller quantum effects. Yet another desirable improvement would be to better account for long-range electrostatic interactions. In this context, it should be noted that a straightforward application of Ewald summations within the framework of a mixed quantum-classical simulation is complicated by the fact that one needs to deal with periodic replicas of the quantum subsystem. Thus, implementing techniques that circumvent this problem, such as those proposed in refs 127-131, would be worthwhile. The removal of these limitations and the extension of this work in the above-mentioned directions is the subject of ongoing research in our group and will be reported in future papers.

VER of Hydrogen-Bonded Complex in a Polar Solvent Acknowledgment. This project was funded by grants from the National Science Foundation and the Petroleum Research Fund. References and Notes (1) Faltermeier, B.; Protz, R.; Maier, M. Chem. Phys. 1981, 62, 377. (2) Brueck, S. R. J.; Osgood, R. M. Chem. Phys. Lett. 1976, 39, 568. (3) Legay-Sommaire, N.; Legay, F. Chem. Phys. 1977, 52, 213. (4) Chesnoy, J.; Richard, D. Chem. Phys. 1982, 67, 347. (5) Chesnoy, J.; Richard, D. Chem. Phys. Lett. 1982, 92, 449. (6) Chateau, M.; Delalande, C.; Frey, R.; Gale, G. M.; Prade`re, F. J. Chem. Phys. 1979, 71, 4799. (7) Delalande, C.; Gale, G. M. J. Chem. Phys. 1979, 71, 4804. (8) Delalande, C.; Gale, G. M. J. Chem. Phys. 1980, 73, 1918. (9) Faltermeier, B.; Protz, R.; Maier, M.; Werner, E. Chem. Phys. Lett. 1980, 74, 425. (10) Chesnoy, J.; Gale, G. M. Ann. Phys. (Paris, Fr.) 1984, 9, 893. (11) Chesnoy, J.; Gale, G. M. AdV. Chem. Phys.1988, 70 (part 2), 297. (12) Harris, C. B.; Smith, D. E.; Russell, D. J. Chem. ReV. 1990, 90, 481. (13) Owrutsky, J. C.; Raftery, D.; Hochstrasser, R. M. Annu. ReV. Phys. Chem. 1994, 45, 519. (14) Elsaesser, T.; Kaiser, W. Annu. ReV. Phys. Chem. 1991, 42, 83. (15) Calaway, W. F.; Ewing, G. E. J. Chem. Phys. 1975, 63, 2842. (16) Laubereau, A.; Kaiser, W. ReV. Mod. Phys. 1978, 50, 607. (17) Roussignol, P.; Delalande, C.; Gale, G. M. Chem. Phys. 1982, 70, 319. (18) Heilweil, E. J.; Doany, F. E.; Moore, R.; Hochstrasser, R. M. J. Chem. Phys. 1982, 76, 5632. (19) Heilweil, E. J.; Casassa, M. P.; Cavanagh, R. R.; Stephenson, J. C. Chem. Phys. Lett. 1985, 117, 185. (20) Heilweil, E. J.; Casassa, M. P.; Cavanagh, R. R.; Stephenson, J. C. J. Chem. Phys. 1986, 85, 5004. (21) Harris, A. L.; Brown, J. K.; Harris, C. B. Annu. ReV. Phys. Chem. 1988, 39, 341. (22) Paige, M. E.; Russell, D. J.; Harris, C. B. J. Chem. Phys. 1986, 85, 3699. (23) Owrutsky, J. C.; Kim, Y. R.; Li, M.; Sarisky, M. J.; Hochstrasser, R. M. Chem. Phys. Lett. 1991, 184, 368. (24) Moustakas, A.; Weitz, E. J. Chem. Phys. 1993, 98, 6947. (25) Kliner, D. A. V.; Alfano, J. C.; Barbara, P. F. J. Chem. Phys. 1993, 98, 5375. (26) Zimdars, D.; Tokmakoff, A.; Chen, S.; Greenfield, S. R.; Fayer, M. D. Phys. ReV. Lett. 1993, 70, 2718. (27) Pugliano, N.; Szarka, A. Z.; Gnanakaran, S.; Hochstrasser, R. M. J. Chem. Phys. 1995, 103, 6498. (28) Paige, M. E.; Harris, C. B. Chem. Phys. 1990, 149, 37. (29) Salloum, A.; Dubost, H. Chem. Phys. 1994, 189, 179. (30) Tokmakoff, A.; Sauter, B.; Fayer, M. D. J. Chem. Phys. 1994, 100, 9035. (31) Tokmakoff, A.; Fayer, M. D. J. Chem. Phys. 1995, 103, 2810. (32) Urdahl, R. S.; Myers, D. J.; Rector, K. D.; Davis, P. H.; Cherayil, B. J.; Fayer, M. D. J. Chem. Phys. 1997, 107, 3747. (33) Owrutsky, J. C.; Li, M.; Locke, B.; Hochstrasser, R. M. J. Phys. Chem. 1995, 99, 4842. (34) Laenen, R.; Rauscher, C.; Laubereau, A. Phys. ReV. Lett. 1998, 80, 2622. (35) Woutersen, S.; Emmerichs, U.; Nienhuys, H.; Bakker, H. J. Phys. ReV. Lett. 1998, 81, 1106. (36) Myers, D. J.; Urdahl, R. S.; Cherayil, B. J.; Fayer, M. D. J. Chem. Phys. 1997, 107, 9741. (37) Myers, D. J.; Chen, S.; Shigeiwa, M.; Cherayil, B. J.; Fayer, M. D. J. Chem. Phys. 1998, 109, 5971. (38) Sagnella, D. E.; Straub, J. E.; Jackson, T. A.; Lim, M.; Anfinrud, P. A. Proc. Natl. Acad. Sci. U.S.A. 1999, 96, 14324. (39) Hamm, P.; Lim, M.; Hochstrasser, R. M. J. Chem. Phys. 1997, 107 (24), 15023. (40) Oxtoby, D. W. AdV. Chem. Phys. 1981, 47 (part 2), 487. (41) Oxtoby, D. W. Annu. ReV. Phys. Chem. 1981, 32, 77. (42) Oxtoby, D. W. J. Phys. Chem. 1983, 87, 3028. (43) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2002, 117, 5827. (44) Sibert, E. L., III; Rey, R. J. Chem. Phys. 2002, 116, 237. (45) Whitnell, R. M.; Wilson, K. R.; Hynes, J. T. J. Phys. Chem. 1990, 94, 8625. (46) Whitnell, R. M.; Wilson, K. R.; Hynes, J. T. J. Chem. Phys. 1992, 96, 5354. (47) Benjamin, I.; Whitnell, R. M. Chem. Phys. Lett. 1993, 204, 45. (48) Bruehl, M.; Hynes, J. T. Chem. Phys. 1993, 175, 205. (49) Gnanakaran, S.; Hochstrasser, R. M. J. Chem. Phys. 1996, 105, 3486. (50) Chorny, I.; Vieceli, J.; Benjamin, I. J. Chem. Phys. 2002, 116, 8904.

J. Phys. Chem. B, Vol. 112, No. 13, 2008 4057 (51) Rey, R.; Hynes, J. T. J. Chem. Phys. 1996, 104, 2356. (52) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 118, 264. (53) Ferrario, M.; Klein, M. L.; McDonald, I. R. Chem. Phys. Lett. 1993, 213, 537. (54) Morita, A.; Kato, S. J. Chem. Phys. 1998, 109, 5511. (55) Rey, R.; Hynes, J. T. J. Chem. Phys. 1998, 108, 142. (56) Gulmen, T. S.; Sibert, E. L., III. J. Phys. Chem. A 2004, 108, 2389. (57) Gulmen, T. S.; Sibert, E. L., III. J. Phys. Chem. A 2005, 109, 5777. (58) Laage, D.; Demirdjian, H.; Hynes, J. T. Chem. Phys. Lett. 2005, 405, 453. (59) Miller, D. W.; Adelman, S. A. Int. ReV. Phys. Chem. 1994, 13, 359. (60) Stratt, R. M.; Maroncelli, M. J. Phys. Chem. 1996, 100, 12981. (61) Deng, Y.; Stratt, R. M. J. Chem. Phys. 2002, 117, 1735. (62) Deng, Y.; Stratt, R. M. J. Chem. Phys. 2002, 117, 10752. (63) Herman, M. F. J. Chem. Phys. 1987, 87, 4779. (64) Li, S.; Thompson, W. H. J. Chem. Phys. 2003, 107, 8696. (65) Terashima, T.; Shiga, M.; Okazaki, S. J. Chem. Phys. 2001, 114, 5663. (66) Sato, M.; Okazaki, S. J. Chem. Phys. 2005, 123, 124508. (67) Sato, M.; Okazaki, S. J. Chem. Phys. 2005, 123, 124509. (68) Sato, M.; Okazaki, S. J. Mol. Liq. 2005, 119, 15. (69) Shi, Q.; Geva, E. J. Phys. Chem. A 2003, 107, 9059. (70) Shi, Q.; Geva, E. J. Phys. Chem. A 2003, 107, 9070. (71) Ka, B. J.; Shi, Q.; Geva, E. J. Phys. Chem. A 2005, 109, 5527. (72) Ka, B. J.; Geva, E. J. Phys. Chem. A 2006, 110, 9555. (73) Ka, B. J.; Geva, E. J. Phys. Chem. A 2006, 110, 13131. (74) Ando, K.; Hynes, J. T. J. Phys. Chem. A 1999, 103, 10398. (75) Zwanzig, R. J. Chem. Phys. 1961, 34, 1931. (76) Landau, L.; Teller, E. Z. Phys. Z. Sowjetunion 1936, 34, 10. (77) Lawrence, C. P.; Skinner, J. L. J. Chem. Phys. 2003, 119, 3840. (78) Bakker, H. J. J. Chem. Phys. 2004, 121, 10088. (79) Fujisaki, H.; Zhang, Y.; Straub, J. E. J. Chem. Phys. 2006, 124, 144910. (80) Hammes-Schiffer, S.; Tully, J. C. J. Phys. Chem. 1995, 99, 5793. (81) Kim, H. J.; Staib, A.; Hynes, J. T. In Femtochemistry and Femtobiology: Ultrafast Reaction Dynamics at Atomic Scale Resolution; Sundsto¨rm, V., Ed.; Imperial College Press: London, 1997; p 510. (82) Herman, M. F. J. Chem. Phys. 1987, 87, 4794. (83) Martens, C. C.; Fang, J. J. Chem. Phys. 1997, 106, 4918. (84) Donoso, A.; Martens, C. C. J. Chem. Phys. 1998, 102, 4291. (85) Donoso, A.; Martens, C. C. J. Chem. Phys. 2000, 112, 3980. (86) Donoso, A.; Kohen, D.; Martens, C. C. J. Chem. Phys. 2000, 112, 7345. (87) Schu¨tte, C. Partial Wigner transforms and the quantum classical Liouville equation. Preprint SC 99-10; Konard-Zuse-Zentrum fu¨r Informationstechnik: Berlin, 1999. (88) Santer, M.; Manthe, U.; Stock, G. J. Chem. Phys. 2001, 114, 2001. (89) Kapral, R.; Ciccotti, G. J. Chem. Phys. 1999, 110, 8919. (90) Nielsen, S.; Kapral, R.; Ciccotti, G. J. Chem. Phys. 2000, 112, 6543. (91) Wan, C.; Schofield, J. J. Chem. Phys. 2000, 113, 7047. (92) Wan, C.; Schofield, J. J. Chem. Phys. 2000, 112, 4447. (93) Kapral, R.; Ciccotti, G. In Bridging Time Scales: Molecular Simulations for the Next Decade; Nielaba, P.; Mareschal, M.; Ciccotti, G., Eds.; Springer-Verlag: Berlin, 2002; pp 445-472. (94) Kapral, R. Annu. ReV. Phys. Chem. 2006, 57, 129. (95) Shi, Q.; Geva, E. J. Chem. Phys. 2004, 121, 3393. (96) Horenko, I.; Salzmann, C.; Schmidt, B.; Schutte, C. J. Chem. Phys. 2002, 117, 11075. (97) Azzouz, H.; Borgis, D. J. Chem. Phys. 1993, 98, 7361. (98) Lippincott, E. R.; Schroeder, R. J. Chem. Phys. 1955, 23, 1099. (99) Hammes-Schiffer, S.; Tully, J. C. J. Chem. Phys. 1994, 101, 4657. (100) Hanna, G.; Kapral, R. J. Chem. Phys. 2005, 122, 244505. (101) Staib, A.; Borgis, D.; Hynes, J. T. J. Chem. Phys. 1995, 102, 2487. (102) Antoniou, D.; Schwartz, S. D. J. Chem. Phys. 1999, 110, 7359. (103) Antoniou, D.; Schwartz, S. D. J. Chem. Phys. 1999, 110, 465. (104) McRae, R. P.; Schenter, G. K.; Garrett, B. C.; Svetlicic, Z.; Truhlar, D. G. J. Chem. Phys. 2001, 115, 8460. (105) Kim, S. Y.; Hammes-Schiffer, S. J. Chem. Phys. 2003, 119 (8), 4389. (106) Yamamoto, T.; Miller, W. H. J. Chem. Phys. 2005, 122, 044106. (107) Hanna, G.; Kapral, R. Acc. Chem. Res. 2006, 39, 21. (108) Marcus, R. A.; Sutin, N. Biochim. Biophys. Acta 1985, 811, 265322. (109) Warshel, A. J. Phys. Chem. 1982, 86, 2218-2224. (110) Kiefer, P. M.; Hynes, J. T. Solid State Ionics 2004, 168, 219224. (111) Bagchi, B.; Oxtoby, D. W.; Fleming, G. R. Chem. Phys. 1984, 86, 257. (112) Maroncelli, M.; Fleming, G. R. J. Chem. Phys. 1988, 89, 875. (113) Barbara, P. F.; Jarzeba, W. AdV. Photochem. 1990, 15, 1. (114) Rossky, P. J.; Simon, J. D. Nature 1994, 370, 263.

4058 J. Phys. Chem. B, Vol. 112, No. 13, 2008 (115) Horng, M. L.; Gardecki, J. A.; Papazyan, A.; Maroncelli, M. J. Phys. Chem. 1995, 99, 17311. (116) Stephens, M. D.; Saven, J. G.; Skinner, J. L. J. Chem. Phys. 1997, 106, 2129. (117) Nielsen, S.; Kapral, R.; Ciccotti, G. J. Chem. Phys. 2001, 115, 5805. (118) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (119) Andersen, H. C. J. Comput. Phys. 1983, 52, 24. (120) Evans, D. J.; Holian, B. L. J. Chem. Phys. 1985, 83, 4069. (121) MacKernan, D.; Kapral, R.; Ciccotti, G. J. Phys. Condens. Matter 2002, 14, 9069. (122) Larsen, R.; Stratt, R. M. J. Chem. Phys. 1999, 110, 1036. (123) Hanna, G.; Kapral, R. To be submitted for publication, 2008.

Hanna and Geva (124) Laria, D.; Ciccotti, G.; Ferrario, M.; Kapral, R. J. Chem. Phys. 1992, 97, 378. (125) Staib, A.; Borgis, D. Chem. Phys. Lett. 1997, 271, 232. (126) Cabral, B. J. C.; Rivail, J. L.; Bigot, B. J. Chem. Phys. 1987, 86, 1467. (127) Staib, A.; Borgis, D. J. Chem. Phys. 1995, 103, 2642. (128) Yang, C.-Y.; Wong, K. F.; Skaf, M.; Rossky, P. J. J. Chem. Phys. 2001, 114, 3598. (129) Turi, L.; Borgis, D. J. Chem. Phys. 2002, 117, 6186. (130) Nicolas, C.; Boutin, A.; Le´vy, B.; Borgis, D. J. Chem. Phys. 2003, 118, 9689. (131) Borgis, D.; Rossky, P. J.; Turi, L. J. Chem. Phys. 2006, 125, 064501.

Vibrational Energy Relaxation of a Hydrogen-Bonded Complex

Recommend Documents