Wave packet dynamics of internal rotational ... - ACS Publications

in water, waterxrythritol, and water:glycerol solutions (36pages). Ordering information is given on any current masthead page. Wave Packet Dynamics of...
0 downloads 0 Views 2MB Size
4200

J . Phys. Chem. 1992, 96,4200-4209

near- Tg spectral modes is seemingly excluded by our measurements. Acknowledgment. The Partial SuPPrt of this work by the National Science Foundation under Grant DMR90-43885 is gratefully acknowledged. Registry NO. H20, 7732-18-5; glycerol, 56-81-5; erythritol, 149-32-6.

Supplementary Material Available: Tables of diffusion coefficient D, spectral variance V, and T / q (q being the solvent viscosity) for polystyrene sphere probes in water, water:erythritol, water:glycerol, and pure glycerol; figures showing D and Vagainst TI,, (q being the solvent viscosity) for polystyrene sphere probes in water, watererythritol, and water:glycerol solutions (36 pages). Ordering information is given on any current masthead page.

Wave Packet Dynamics of Internal Rotational Predissociation in NeHF and NeDF William J. Hovingh and Robert Parson* Department of Chemistry and Biochemistry, University of Colorado, and Joint Institute for Loboratory Astrophysics, University of Colorado and National Institute of Standards and Technology, Boulder, Colorado 80309-0440 (Received: November 18, 1991; In Final Form: February 17, 1992)

We investigate the process of internal rotational predissociation in the NeHF and NeDF van der Waals complexes using time-dependent quantum mechanical calculations. The computations are performed in three degrees of freedom and for arbitrary values of the total angular momentum; they use a combination of the Fourier wave packet propagation method in the dissociative coordinate with a body-fixed basis set expansion in the intramolecular angular coordinates. Calculated values of the predissociation lifetimes of several states of these complexes are found to be in excellent agreement with the results of time-independent quantum calculations. Product distributions are easily extracted from the evolving wave packets. The fragmentation mechanisms of several metastable states are further explored by examination of the time-evolving probability distribution.

I. Introduction Time-dependent quantum mechanical methods are now available as powerful tools in the study of the dynamics of van der Waals molecules. In the laboratory, advances in high-resolution far-infrared and near-infrared vibrational spectroscopies have allowed the study of low-frequency modes of many van der Waals complexes, notably in dimers consisting of a rare gas atom At the same time, ab initio and a hydrogen halide chemists have been able to provide increasingly reliable potential energy surfaces for such systems, fueling theoretical investigations of the unusual dynamics of these weakly bound, floppy complexes. Most dynamical studies to date have involved time-independent quantum mechanical calculations or classical trajectory simulations. However, rapid increase in available computing speed and concurrent advances in computational techniques6 have now opened up the possibility of fully quantum mechanical, time-dependent studies of the dynamics of small van der Waals systems. A feature of particular interest in these weakly bound complexes is the ability to see evidence of low-energy metastable states. These states predissociate, and transitions to them appear as lifetimeIn broadened lines in spectra of the van der Waals comple~.~,~.' order to employ these observations in understanding the behavior of the molecule, it should be helpful to have in mind clear dynamical pictures of the dissociation process-exactly the kind of information that can best be obtained from a time-dependent calculation. A dissociation mechanism of particular interest is internal rotational predissociation, in which rotational energy of one fragment in the complex is converted into relative translational motion of the fragments. This process has been observed experimentally in ArHD? NeHF? HeHF, and HeDF' and has been (1) Nesbitt, D. J. Chem. Rev. 1988, 88, 843. (2) Nesbitt, D. J.; Lovejoy, C. M.; Lindeman, T. G.; ONeil, S. V.; Clary, D. C. 1. Chem. Phys. 1988, 91, 722. (3) Lovejoy, C. M.; Nesbitt, D. J. J . Chem. Phys. 1991, 94, 208. (4) Lovejoy, C. M.; Nesbitt, D. J. J . Chem. Phys. 1990, 93, 5287. (5) Saykally, R. J. Acc. Chem. Res. 1989, 22, 295. (6) Kosloff, R. J. Phys. Chem. 1988, 92,2087. (7) Ashton, C. J.; Child, M. S.; Hutson, J. M. J . Chem. Phys. 1983, 78, 4025. (8) McKellar, A. R. W. J . Chem. Phys. 1974,61,4636. McKellar, A. R. W. Faraday Discuss. Chem. Sot. 1982, 73, 89.

0022-365419212096-4200$03.00/0

studied with time-independent quantum mechanical methods in ArHCl,' ArHD,9 and HeHF,'O among others. This process depends strongly on the anisotropy of the intermolecular potential and can therefore act as a sensitive probe of that feature of the interaction. In this paper, we report a time-dependent quantum mechanical study of internal rotational predissociation in NeHF and NeDF. We have adopted an approach, developed by Gray and Wozny,l' which combines wave packet propagation in the intermolecular separation coordinate with a basis set expansion in the remaining degrees of freedom; the method is adapted to treat a complex made up of an atom and a rigid rotor. The calculations reported agree well with experiment and fully with other calculations; they are also fully converged, without requiring excessive computational effort. The rest of the paper is organized as follows. In section 11, we discuss the method used and present details of the calculations. In section 111, we briefly review what is already known, from experiments and previous calculations, about the specific molecules to be treated. In section IV, we present the results of our calculations and compare the computed predissociation lifetimes of selected states with the results of experiment and time-independent calculations. We also show the time evolution of the probability density in several cases, and in section V we examine that time evolution and discuss the mechanistic pictures of the fragmentation dynamics that emerge from that examination. 11. Method

-

In thissection, we outline the method used to study the quantum dynamics of the fragmentation process XBC X + BC, where X is an atom bound by a weak, weakly anisotropic potential to BC, a diatom. In particular, we consider the case in which the diatom can be treated as a rigid rotor and the dissociation of the triatomic complex is due to the transfer of energy from rotational motion of the diatom to X-BC stretching motion. Our method is in the spirit of Gray and Wozny's application of the closecoupled wave packet formalism to the vibrational predissociation (at total J = 0) of NeCl, and related complexes,"J* with this (9) Hutson, J. M.; LeRoy, R. J. J. Chem. Phys. 1983, 78, 4040. (10) Gianturco, F. A.; Palma, A.; Villareal, P.; Delgado-Barrio, G. Chem. Phys. Lett. 1984, I 1 I , 399. (11) Gray, S. K.; Wozny, C. E. J. Chem. Phys. 1989, 91, 7671.

0 1992 American Chemical Society

The Journal of Physical Chemistry, Vo1. 96, No. 11, 1992 4201

Internal Rotational Predissociation in NeHF and NeDF difference: we freeze the diatom stretching motion, while including all angular degrees of freedom. Thus, the problem is solved in three dimensions and for arbitrary values of total angular momentum. In what follows, we describe the representation and propagation of the time-dependent wave function and discuss the selection of an initial state; we then show how observable quantities are extracted from the calculation and how the results are tested for convergence. A. Constructing the Wave Packet Equations. Using the closecoupled wave packet (CCWP) formalism developed by Kouri and ~ e w o r k e r sand ~ ~Jackson and Metiu,14we will obtain solutions of the time-dependent Schrijdinger equation i-a+ at

with respect to a space-fixed frame. 5 , k is a spherical harmonic that describes the orientation of the diatom in the body-fmed frame of the complex. Since the body-fixed frame is defined by the condition that the orbital angular momentum has no projection along the body-fixed z axis, it follows that J, = j , and therefore that the D and Y functions are indexed by the same k quantum number. Since the third Euler angle in the rotation matrix element coincides with the body-fmedazimuthal angle (p, it has no physical significance and is omitted in defining the basis functions.18 We now substitute eqs 2-5 into eq 1 to obtain a set of equations for the time dependence of the channel packets C,,k(R,t):

&+

(Note that we employ atomic units, such that h = 1.) In Jacobi coordinates, the Hamiltonian for the interaction of an atom with a rigid rotor is

&=

--!-(

+ B? + 2 2 + V(R,B)

Rl-$R)

(2)

2pR2

2P

where p is the atom-diatom reduced mass; B is the rotational constant of the diatom and j its angular momentum; C = IJ - jl, the orbital angular momentum of the complex; and the interaction potential Vdepends on the length of R (the vector from the diatom center of mass to the atom) and the angle, 6,that R makes with the diatom axis. We will expand V(R,B) in a Legendre series V(R,B) = CVA(R) P ~ ( C O0S) A

(3)

When the anisotropy of the potential is large compared to the spacing between end-over-end rotational levels of the complex, but both are small compared to the rotational level spacing of the diatom itself, we expect that two nearly good quantum numbers will bej, the magnitude of the angular momentum of the diatom, and k,its projection on the bodyfixed z axis (that is, its projection along R). Thus, to obtain the solutions of eq 1 with this Hamiltonian, we expand the wave function in a set of basis functions in body-fixed coordinates

Here, (@,a)are the polar coordinates of the vector R in space, and (B,q)are the polar coordinates of the diatom axis in the body frame. The superscripts J, M,and p'denote the three rigorously good quantum numbers for this problem: J and M are respectively the total angular momentum and its projection on the space-fixed z axis; p' denotes the symmetry of the wave function under reflection through the plane of rotation of the complex, which can be even or odd, corresponding to the labels (e&. The factor R' is included in the expansion to simplify the form of the radial kinetic energy term in the close-coupled equations. We choose the basis functions q , k to be eigenfunctions of the body-fmed angular momentum operators j and j , as well as of J, J,, and reflection symmetry:Is

Here, d M , k is a rotation matrix element with the convention of Brink and Satchler;I6J7it describes the orientation of the complex (12) Gray, S.K.; Wozny, C. E. J . Chem. Phys. 1991, 94, 2817. (13) Kouri, D. J.; Mowrey, R. C. J . Chem. Phys. 1983,52, 3 5 . Sun, Y.; Kouri, D. J. J . Chem. Phys. 1988,89,2958. Sun, Y.; Judson, R. S.;Kouri, D. J. J. Chem. Phys. 1989, 90, 241. (14) Jackson, B.; Metiu, H. J . Chem. Phys. 1987, 86, 1026. (15) Schatz, G. C.; Kupperman, A. J . Chem. Phys. 1976,65, 4642, Appendix A.

These are the close-coupled wave packet equations for an atomrigid-rotor system in body-fmed coordinates. They have the form

a-

i-C(R,t) at

1

a2

+ Q(R)

= -- -Z(R,t) 2 r aR2

e(R,t)

(7)

where the Cj,k(R,t) are written as a column vector. The matrix Q incorporates the matrix elements of all terms in the Hamiltonian which involve angles and angular momenta. Specifically, this matrix includes the rotational kinetic energy of the diatom and a Coriolis term that is diagonal in 0.k); there are two Coriolis terms that are diagonal in j but which kinetically couple neighboring k levels; finally, the matrix elements of the potential are diagonal in k but couple different channels with different j . (In a space-fixed formulation, where the basis functions are eigenfunctions of j and C, all kinetic terms are diagonal, but the potential matrix elements are more complicated and have offdiagonal terms in both j and C.) B. Propagating the Wave Packet Equations. We solve eq 7 using the split operator method of Fleck et al.I9 in which the propagation at each time step is approximated by 1

a2

( 'if ) exp(-) exp(-yQ(R))e(R,t) 2 r aR2

= exp --Q(R)

iAt

(8)

In practice, the Cik(R,t) are defined on an evenly spaced grid in R. The matrix exp(-iA@) is computed at each grid point; the exponentiated partial derivative term is evaluated using Fourier transforms

and the leftmost exponential term is incorporated into the prop agation at the next time step. Thus, the computation at each step (16) Zare, R. N. Angular Momentum; Interscience: New York, 1988. (17) Brink, D. M.; Satchler, G. R. Angular Momentum, 2nd ed.;Clarendon Press: Oxford, 1968. (1 8) Brocks, G.; van der Avoird, A.; Sutcliffe, B. T.; Tennyson, J. Mol. Phys. 1983, 50, 1025. (19) Fleck, J. A.; Morris, J. R.; Feit, M. D. Appl. Phys. 1976, IO, 129.

4202 The Journal of Physical Chemistry. Vol. 96, No, 11. 1992 consists of a matrix multiplication at each grid point and two fast Fourier transform operations in each of the 0,k) channels. The method is accurate to third order in At at each time step and is unconditionally stable. However, care must be taken to ensure that the grid spacing is small enough that spatial variations of the wave function are adequately represented and that the time step is small enough to ensure that unacceptable phase errors do not accumulate over the course of the calculation. In addition, absorbing boundary conditions20,21are applied to damp outgoing amplitude before it reaches the edge of the grid at large R; otherwise, because the fast Fourier transform assumes a periodic form for C(R), amplitude will come out that end and spuriously enter the other. (No damping is needed a t small R, because the hard wall of the potential prevents amplitude from reaching that edge of the grid.) We implement the damping by multiplying, at each time step, by an attenuation function A(R) = 1, R < RD A(R) = exp(-a(R

- RD)2At),

R

> RD

(10)

C. Constructing the Initial Wave Function. Before a wave function can be propagated we must, of course, choose an initial state. In studying the fragmentation of XBC, we suppose that the complex is initially in some quasi-bound state of the Hamiltonian. However, because we work in a finite range in R, we can at best approximate any quasi-bound state with some L2 function. To construct such an approximation, a simple stabilization procedure22is employed. The channel packets Cj,k(R,t) are expanded in an .C2 basis set; in this case we choose a distributed Gaussian basis,23so that the full wave function is expressed as iJ,k

(Evenly spaced Gaussians of equal width give reasonable convergence.) We then repeatedly diagonalize the Hamiltonian in the combined Gaussian/body-fixed basis, varying either the Gaussian spacing or the number of basis Gaussians; eigenvalues that are stable under these variations correspond to the quasibound states we seek to approximate. Choosing the wave function corresponding to one of these eigenvalues, we have constructed an initial state which can now be propagated. We should make it clear at this point that the eigenfunction from the stabilization calculation, which we use as an initial condition for solution of the timedependent Schrijdinger equation, is an approximation both to the wave function of the predissociating state prepared by the experiment and to a scattering resonance wave function with outgoing asymptotic boundary conditions. It is not exactly the experimentally prepared wave function; that depends sensitively on the interaction between the molecule and the laser field, and its computation is beyond the scope of our methods. On the other hand, the stabilized wave function cannot be an exact resonance wave function, since the latter has infinite extent in the scattering coordinate, while our calculations are performed on a finite range in R. Nonetheless, the stabilization calculation will produce a close approximation to the resonance wave function in the interaction region, and the time dependence of this approximation in the limit of long times should closely approximate that of the true resonance wave function. This expectation is borne out by our ability to extract resonance lifetimes from the decay of the autocorrelation function. Furthermore, the stabilized wave function adequately represents the “important” component of the wave function of the complex in a CW laser experiment, in this sense: The radiation field excites the complex to a superposition of states with energies tightly grouped around the resonance energy in question;24each of these is finite in its extent in R (since the wave function before laser (20) Kosloff, R.; Kosloff, D. J . Comput. Phys. 1986, 63, 363. (21) Kosloff, D.; Kosloff, R. J. Comput. Phys. 1983, 52, 35. (22) Hazi, A. U.; Taylor, H . S.Phys. Reo. A 1970, 1, 1109. (23) Light, J. C.; Hamilton, I. P. J . Chem. Phys. 1986, 84, 306. (24) Shapiro. M.; Bersohn, R. Annu. Reu. Phys. Chem. 1982, 33, 409.

Hovingh and Parson excitation, in the vibrational ground state of the complex, is also bounded in R). In the limit of an isolated resonance, the contributions to this superposition which are not the resonance wave function will be continuum wave functions that rapidly decay into the asymptotic region, leaving behind the wave function of a (relatively) slowly decaying metastable state-essentially the same wave function that is left after the nonresonance contributions to the stabilized wave function have disappeared in our calculations. Consequently, we have confidence that, at least in the absence of overlapping rcsonances, the stabilization procedure used here to generate initial states is perfectly adequate. In the simulations reported here, we have found that all significant results are completely independent of the details of the initial wave function. D. Extracting Observable Quantities. Ultimately, we wish to determine the values of observable quantities from the time-dependent wave function. In particular, we can find the lifetime and resonance energy of the quasi-bound state whose dissociation we model. The distribution of the fragments among the open product channels is also of interest. One way to measure the lifetime is to monitor the probability that survives within a large spherical shell in R

(12)

where R*lies well beyond the range of the potential. If the state represented by rl, dissociates with a single lifetime T , the inverse lifetime I’= 1/ T should be the slope of a plot of In (Psur)against t . A second measure of lifetime depends on our decision to choose as the initial condition a wave function that closely approximates a quasi-bound state of the complex. We assume that the initial wave function is nearly an eigenfunction of the Hamiltonian with a complex e i g e n v a l ~ e . ~That ~ is, rl,(O) = where

The autocorrelation function will then have the simple time dependence (+(O)(+(t))

= e-iEJe-rJ/2

(14)

and the slope of a plot of In I(rl,(O)lrl,(t) )I2 against t will again yield the inverse lifetime, rr. (Departures from the simple Breit-Wigner form will appear as nonexponential contributions to the correlation function.) In addition, the resonance energy E, can be estimated by Fourier transform analysis of the autocorrelation function. The resonance energy so estimated should be approximately equal to the stabilized eigenvalue found in the procedure to choose the initial state, although it will differ somewhat, since the true resonance wave function includes spatially extended components that the stabilized wave function cannot. Gray and WoznyI2 have used the method of Baht-Kurti et to determine the distribution, among available product states, of the products of a unimolecular dissociation process studied with time-dependent methods. In the present work, we have adopted a simpler procedure, which we find provides useful results that agree well with the method of ref 26. We monitor the wave function at some point on the grid in R,far removed from the interaction region (but also removed from the damping region at the outside edge of the grid). Transforming the channel amplitudes at those points into a space-fured basis, we calculate the probability density, at that point, in each available product channel. Normalizing the total probability density to unity, we then plot the time dependence of the relative flux density and observe that, after an initial treatment, the probability in each channel settles down and oscillates about a nearly constant value. (25) Taylor, J. R. Scattering Theory; Interscience: New York, 1972; Chapter 13. (26) Balint-Kurti, G.C.; Dixon, R. N.; Marston, C. G.J. Chem. SOC., Faraday Trans. 1990, 86, 1741.

The Journal of Physical Chemistry, Vol. 96, No. 11, 1992 4203

Internal Rotational Predissociation in NeHF and NeDF TABLE I: Testing Convergence with Respect to Time Step Size“

A%,m, MHz

overlap from survival from correlation amplitude A f , au probabilityb functionC erroP 2000 1000 500 250 125

:::::: ::::::

1031.96 883.40 885.62 886.17 886.31

1031’73

0.033341 0.002 283 0.000574 0.000 298

overlap phase error 0.63544 0.143 95 0.035 33 0.008 79

Comparison of the results of five propagations with identical initial states and parameters, varying only the size of the time step. The state computed here is the Z bend resonance of NeHF on the semiempirical surface of ref 3. bEquation 12. CEquation14. “Equation 16. Finally, the great strength of time-dependent methods is our ability to examine the time development of the wave function itself. In section IV, we present plots of a reduced probability density function

Here, besides integrating out the orientation angles (Y and P, we have averaged over the body-fixed azimuthal angle cp. This averaging leaves a function of R and 8, reflecting the dependence of the potential on these coordinates and its independence of cp; the (relatively) simple surface and contour plots of p will help us to elucidate the dynamics of the predissociation process. E. Testing Convergence of the Results. In testing the results of these calculations for precision, we considered the convergence, with respect to the spacing of the grid in R and the size of the time step (At), of two kinds of numerical measures: the overall amplitude and phase of the wave function at the end of the calculation, and the value of the observable quantities discussed in section 1I.D. We have found that it is possible to obtain excellent convergence without excessive computational effort. For this reason, we have not used any of the more powerful propagation schemes that are available. To assess the convergence of the wave function itself, we seek measures of the amplitude and phase differences between two wave functions. Following Leforestier et aL2’ we define the overlap errors in amplitude, t12,and phase, aI2,between two wave functions and qZto be

61,2

= mod

(($11+2))

(16)

(The form of the definition o f t differs slightly from that in ref 27 to account for the exponential decrease in the norm of our dissociative wave functions.) We then propagate the same initial wave function, over the same total time span, with a series of different time step or grid spacing parameters and observe the behavior of the amplitude and phase errors between successive final wave functions in the series. Similarly, we can monitor the behavior of the inverse lifetime (found from the slopes of semilog plots of the survival probability and squared autocorrelation function) as the same parameters are varied. In Table I, we present the results of such a test, in which the time step in the problem is reduced by successive powers of 2 from 1000 to 25 au; the wave function is propagated for lo5 au (1 au = 2.42 X p). Note that the lifetimes calculated are converged to within 0.1% at a time step of 500 au. Note also that the overlap amplitude error is of the same order of magnitude as the error in lifetime between successive propagations. We have found this to hold true independent of the state studied. Similar convergence tests were performed by varying the grid spacing and number of body-fixed channels; all of the calculations reported in this paper were performed with time steps of 500 au and a grid spacing of ~

~~

(27) Leforestier, C.; Bisseling, R.; Cerjan, C.; Feit, M. D.; Friesner, R.; Guldberg, A.; Hammerich, A.; Jolicard, G.;Karrlein, W.; Meyer, H.-D.; Lipkin, N.; Roncero, 0.;Kosloff, R. J . Compur. Phys. 1991, 94, 59.

0.25 au and included all channels with j I5. We very conservatively estimate that the error in our reported numerical lifetimes (due to finite grid spacing, time step, and basis set size) is less than 0.2%. Another feature of the convergence test to note is that the ratio of overlap phase errors between successive propagations is approximately 4, suggesting that the measure of error falls off as (At)2. The propagation is third order in A? at each time step, but the total number of time steps in a given calculation increases as l/At; consequently, the error at the end of the propagation is second order in At. Finally, we point out that the bulk of the computing time used in the calculations reported here was consumed in the stabilization calculations used to generate the initial states; the propagations themselves were relatively inexpensive. 111. Bound and Quasi-Bound States of NeHF and NeDF

The Ne-HF van der Waals complex provides an ideal system on which to apply the method described here. This molecule satisfies the criterion that the potential anisotropy is larger than the rotational energy level spacing of the complex but smaller than the spacing of HF rotational levels; thus, as discussed by HutsoqB the choice of a body-fixed basis set is appropriate. A high-quality a b initio potential energy surface exists,29 and high-resolution infrared spectroscopic studies of the complexl and its deuterated analogue) have made it possible to construct a refined (semiempirical) potential surface and have suggested the presence of novel rotational predissociation dynamics. In this section, we briefly describe the observed bound and quasi-bound states of these complexes; more detailed discussion can be found in ref 2. High-resolution infrared spectra of NeHF were first reported by Nesbitt et al. in 1988.2 Nearly simultaneously with the experimental work, Nesbitt and Clary performed time-independent quantum mechanical calculations on the complex, using an ab initio potential energy surface of ONeil et and predicted the spectrum with remarkable success. From this combined work, a picture of unusual vibrational dynamics emerged. On the ab initio potential energy surface, the strongest binding is -65 cm-I, found at the linear NeHF geometry; a secondary well of -39 cm-’ is in the linear NeFH configuration. Between these minima, there is a saddle point at -27 cm-I. Thus, the anisotropy in the potential is smaller than twice the HF rotational constant, 28 = 40 cm-’ (the spacing between the first two rotational levels of free HF). Consequently, one expects that the energy levels of the complex (at a given total angular momentum) will look like vibrational levels in the van der Waals stretch superimposed on the rotational energy levels of a nearly free rotor. In Hutson’s nomenclature2*(adapted from Bratoz and Martin30), this is “case 2“ coupling: the diatom angular momentum is strongly coupled to the intermolecular axis, so that a body-fixed basis is appropriate, while the rotation of the diatom is only slightly hindered by the potential. The infrared spectrum of this complex consists of transitions from the ground state of the molecule to states in which the HF has one quantum of vibrational excitation. Since the frequency of this vibration is much higher than any other frequency in the problem, this degree of freedom is decoupled from all the others, and the excitation process essentially puts the complex onto a new vibrationally adiabatic potential energy surface. The accuracy of this picture is confirmed by the small frequency shift from the monomer transition as well as by the absence of lifetime broadening that can be attributed to vibrational predissociation; this is true for all of the rare gas + HF complexes. Experimentally, two principal excitation bands are observed. The first, which lies very near the HF monomer u = 1 0 origin, is the fundamental band in which none of the low-frequency degrees of freedom are excited. This band is 10 times less intense than a combination band that lies 43 cm-l to the blue, which is

-

(28) Hutson, J. M. In Advances in Molecular Vibrations and Collision Dynamics, in press. (29) ONeil, S. V.;Nesbitt, D. J.; Rosmus, P.; Werner, H.-J.; Clary, D. C. J . Chem. Phys. 1989, 91, 71 1. (30) Bratoz, S.; Martin, M. L. J . Chem. Phys. 1965, 42. 1051.

4204 The Journal of Physical Chemistry, Vol. 96, No. 11, 1992

identifiable as a set of transitions to states in which the HF has acquired one quantum of rotational excitation. (This unusual ratio of intensities provides spectroscopic evidence that the HF undergoes nearly free internal rotation in the complex.) The calculations predict that these internal rotation excited states should be unbound by about 22 cm-'. In addition, they predict the existence of two bound states in which the van der Waals stretch is excited but HF rotation is not and an additional internal rotation excited state at about 15 cm-I above the R zero of energy. The stretch excited states are expected to have a negligible transition dipole moment to the ground state; neither they nor the state a t +15 cm-' is observed in NeHF. It is in the HF rotation excited states, none of which is bound in NeHF, that the novel dynamics appear. These are, in fact, quasi-bound states which suffer some kind of predissociation dynamics. The first three were predicted by the calculations of ref 29 to be characterized by the approximate quantum number j = 1 and to correlate with the triply degenerate] = 1 state of free HF. The lowest-lying of the three states is dominated by the j = 1, k = 0 channel (in the notation of eqs 4-6), which is coupled by the anisotropy of the potential to the dissociative j = k = 0 channel. It is predicted to predissociate so rapidly that the corresponding transitions will be lifetime broadened beyond the limits of detection; in fact, these transitions are not observed. This state is denoted 2 in the discussion of Nesbitt et al., since the wave function for the HF rotor in the body frame is dominated by the azimuthally symmetric k = 0 channel. The next two quasi-bound states are dominated by t h e j = lkl = 1 channels but are distinguished by the reflection-symmetry quantum number; one is even and the other (at slightly lower energy), odd. Just as the first quasi-bound state is denoted 2, these two states are denoted IIc and IIf,respectively. Since all states (on this vibrationally excited surface) that lie below the IIr state have even reflection symmetry, and since that is a rigorously good quantum number, this state has no available exit channel on this surface and can only dissociate through the negligible coupling between the HF and van der Waals stretches; hence, it can be considered to be bound. Transitions to this state appear as the instrumentally narrow Q branch lines in the intense combination band of the spectrum. In contrast, transitions to the IIe state appear as lifetimebroadened lines in the P and R branches of the same band, with line widths that depend on the total angular momentum of the final state and increase as J(J I), approximately. This J dependence of the line widths was understood by Nesbitt et al. to arise from the following predissociation mechanism: The dominant channel in this state, j = lkl = 1, is not directly coupled by the Hamiltonian to the only available exit channel, j = k = 0. It is, rather, coupled to the j = 1, k = 0 channel by a Coriolis term that depends on J as [ J ( J l)ll/*. As noted earlier, i,t is this second channel that is coupled by the potential to the exit channel. Consequently, the predissociation rate depends on both the Coriolis coupling term, with its dependence on J , and the potential coupling of t h e j = 1, k = 0 channel to the exit channel. This state of affairs can be illustrated by making plots of the radially adiabatic effective potential (the so-called Born-Op penheimer angulal-radial separation, or BOARS,potential). This function of R is computed by finding the eigenvalues of the effective potential matrix, Q in eq 6,at a given value of R. When the first few BOARS curves are plotted against R, in Figure la, we see that the lowest is an isolated curve that forms an attractive well at small R and tends to zero from below as R grows large. (The zero of energy is defined as R with no angular momentum in the diatom.) The lowest eigenstate of radial motion in this well corresponds approximately to the ground state on this potential surface. Higher-energy eigenstates in this well correspond to van der Waals stretch excited states; all of these are bound. At about 40 cm-I higher in energy (that is, 2BHF)we find a trio of wells that asymptotically approach the energy of HF with one quantum of angular momentum; even the lowest eigenstates in these wells all have positive energy eigenvalues and correspond to the 2 , If', and IIc bend excited states discussed above. The

Hovingh and Parson

-

j=l

j=O

-40L-.----80 4

8

12

R

16

20

/ a.u.

j=2

j=l

j=O

f2 g -20- l 0 I

q,

-30

+

4

,

,

,

12

18

20

-80

4

8

R / a.u.

+

-

'

404

\

3'

1

c

4.

t

k

4

8

04

r7

J=U

-10 -20

::lL -50

-80 4

8

R /12a.u.

18

20

Figure 1. BOARS curves for (a, top) NcHF at J = 2, (b, middle) NeDF at J = 2, and (c, bottom) NeDF at J = 16.

Internal Rotational Predissociation in NeHF and NeDF

The Journal of Physical Chemistry, Vol. 96, No. 11, I992 4205 TABLE II: Comparisons of Line Widths of the neStates in NeHF, As Measured Experimentally' and Calculated with Time-Independentb and Time-DewndenF Methods on the ab Initio Potentiald

J

exDeriment

time-independent calculation

this work

1 2 3 4

61 192 334 492

40 112 214 334

38.5 112 213 332

'Reference 2. bCorrected values, reported in ref 3. 'This work. dReference 29. -0.20

I

I

I

I

I

I

I

I

0

20

40

60

80

100

120

140

time / psec

TABLE 111: Comparison of Resonance Energies and Line Widths of II' States in NeHF, from Calculations on the Semiempirical

Figure 2. Semilog plots of the survival probability and squared autocorrelation function for the nestate of NeHF at J = 3.

adiabatic surfaces with even reflection symmetry are, in this radially adiabatic picture, coupled by kinetic energy terms in R. Hence, since none of the eigenstates on the surfaces of excited internal rotation is bound (i.e., lies below the energy zero), the two with the same symmetry as the lowest surface, and ne,must in fact be metastable states of the system. The situation is somewhat different in the NeDF van der Waals complex. Because the rotor constant of DF is about half that of HF, all of the states characterized by the approximate quantum numberj = 1 are bound (at least at low J). Nesbitt and Lovejoy have observed transitions to those states, as well as a van der Waals stretch combination bandn3They also observed that the rotational progressions of transitions to the even-symmetry internally excited states all terminate abruptly at modest values of final J. They attributed this to rapid, direct rotational predissociation of the complex through the centrifugal barrier. (That is, the energy of end-over-end rotation of the complex exceeds the binding energy of the complex; since this motion is directly coupled to motion in the dissociative coordinate, the complex dissociates rapidly.) Figure l b shows the first four BOARS surfaces for NeDF at total J = 2. Because the rotation constant of DF is half that of HF, thejDF= 1 surfaces are now nested deeply inside the jDF =0 surface, and their radial ground states are bound. When the total angular momentum is raised to 16 (Figure IC), however, all of the surfaces are raised by the centrifugal pseudopotential, and even the ground state lies above the energy zero and so is, in fact, a shape resonance rather than a truly bound state. Analysis of the J dependence of this direct predissociation behavior allowed Nesbitt and Lovejoy to determine the binding energy of NeDF.3 This information, along with the vibrational frequencies and rotation constants of NeDF, was used to modify the ab initio potential to bring its predictions into better agreement with experiment. The resulting semiempirical potential has minima at -86 cm-l (in the linear Ne-H-F geometry) and -55 cm-' (Ne-F-H) and an intervening saddle point at -39 cm-'.

IV. Results A. Lifetimes of Bend Excited States of NeHF. We have used the CCWP method, as described in section 11, to extract a variety of quantitative and qualitative results for the NeHF and NeDF complexes. The first numerical results that we obtain from our calculations are lifetimes 7 for the dissociation of the quasi-bound bend excited states of NeHF. We present these results as line widths, Avfwh = ' / t m ,for ease of comparison with other sources. Figure 2 shows semilog plots of the survival probability (eq 12) and the square of the autocorrelation function (eq 14) for the p r e d i i a t i o n of the NeHF nestate at J = 3. In both plots, there is initially some transient, nonexponential behavior, which we expect since the initial wave function is not exactly that of a quasi-bound state. This transient behavior persists for a few picoseconds, after which both functions decay exponentially with the same time constant. We find that the shape of the transient depends on the details of the initial state, but the exponential behavior at long times does not.

J 1 2 3 4

5

time-independent calculitionb Er,dcm-l L\Yfwhmr MHz 11.2668 11.9025 12.8541 14.1189 15.6937

time-dependent calcuiation' E,: cm-I AYywhm, MHZ 11.274 11.909 12.859 14.122 15.695

74.2 f 3.3 205.2 f 6.8 394 f 15 6 2 9 f 16 889 f 18

68.7 0.1 201.9 f 0.4 392.2 i 0.8 624.1 f 1.2 883.2 f 1.8

'Reference 3. bLovejoy, C. M., private communication. 'This work. dEnergies referenced to the energy of HF with j = 0 and u = 1. eEnergy of the stabilization eigenfunction used as an initial condition in the lifetime calculation.

=

o0 0

, 20

2 40

~60

80 -

100;

120;

140 :

time / psec Figure 3. Time-dependent probability amplitude at R = 20 au in the three available exit channels in the dissociation of the Ae resonance of NeHF at J = 5. The total flux is normalized to unity at each time step, and the relative cross sections calculated with the method of ref 26 are indicated by the arrows at the right of the plot.

In Table 11, we present our calculated line widths for the ne state, comparing them with experimental line widths and with the time-independent scattering calculations of Lovejoy and N ~ b i t t . ~ Both sets of calculations were performed using the ab initio surface, which has been fit to a Maitland-Smith potential in the first five terms in the Legendre expansion of eq 3. Clearly, the two kinds of calculations agree well, while neither agrees particularly well with the experimental data. This discrepancy is to be expected, since no ab initio potential surface is likely to match perfectly, in all of the details which sensitively affect predissociation lifetimes, the potential surface of the actual complex; moreover, this particular surface is known to misestimate certain features, such as the well depth, of the potential. Nonetheless, the theoretical results match experiment reasonably well in magnitude and quite well in terms of the J dependence in the line widths. In Table 111, we present a similar comparison of timedependent and time-independent calculations using the semiempirical intramolecular potential; we also include calculations of the resonance energies of the quasi-bound states. For our timedependent calculations, the energies presented are the energy eigenvalues of the stabilized wave functions used as initial conditions for the

4206

The Journal of Physical Chemistry, Vol. 96, No. 11, 1992

Hovingh and Parson t = 19.35 psec

t = 19.35 psec

1

oa

0.8

P 0.4

oa

0

i

i io

;1

ii

in

i

in

R

t = 19.83 psec 1

Od

ao

P 0.4

oa

0

a

io

a

I

p

:

R

t = 20.32psec 1

Od

0.0

P OA

oa

0

I4

t

=

10

IO

10

16

20.80psec

1

(u

P a4

oa

I

0 4

6

6

10

I?,

14

SD

R

Figure 4. Snapshots of the reduced probability density (eq 15), p(R,B), in the dissociation of the Z bend excited state of NeHF at J = 0.

Internal Rotational Predissociation in NeHF and NeDF calculation. Because the full contribution of the continuum Cannot be included in the stabilization calculation, we expect that these energies should be somewhat higher than the resonance energies calculated in a time-independent scattering calculation. An examination of these results suggests that the two methods are complementary: the CCWP method obtains resonance lifetimes with high precision, while the time-independent method is more precise in calculating resonance energies. (It is possible, in principle, to obtain accurate resonance energies from CCWP calculations by Fourier analyzing the autocorrelation f~nction.~' However, to do so with precision comparable to that of time-independent methods requires following the dynamics for impractically long times.) This complementarity, and the excellent agreement in line width results, has held true for calculations on a variety of other quasi-bound states in NeHF and NeDF, with lifetimes ranging from 50 ps to 10 ns. B. product State Dlstributiom. Another observable of interest is the distribution of products among the available exit channels, when more than one is available. A rough estimate of these distributions can be obtained simply by monitoring the time-dependent probability density in each space-fixed channel at some fued radius in the asymptotic region, as discussed in section 1I.D. In Figure 3 we have plotted the time dependence of the probability ICj,t(RJ)(2for a predissociating resonance with multiple exit channels-in this case, the NeHF resonance dominated by t h e j = lkl = 2, even symmetry, channel, with total J = 5 . Here, R , = 20 au. The three available exit channels a r e j = 0, C = 5 ; j = 1, C = 4; and j = 1, C = 6. (This "Ae" resonance is not experimentally observed because of its small transition moment with the ground state.) After some initial transient behavior, the three probabilities (which have been normalized at each time step to sum to 1) settle down to small oscillations about fixed values. Clearly, this method only gives numbers accurate to several percent, but such results could provide adequate comparison with experiment in a system that possessed an observable predissociating resonance with multiple exit channels. More precise estimates can be found using the methods of ref 26,if the resonance energy is known exactly. The arrows in Figure 3 represent the relative cross sections for dissociation into each of the three exit channels, as calculated by the method of ref 26;the result is clearly in reasonable agreement with our simpler scheme. As can be Seen in the figure, the products of this particular dissociation process are predominantly in the j = 1, C = 6 state. We find that there is a consistent, strong propensity for this Ae resonance to dissociate into the j = 1, C = J 1 channel. We plan to discuss this result in depth in a future publication. C. Dynamical Pictures. The value of solving these problems in a timedependent manner is, obviously, that the results generated yield time-dependent (that is, dynamical) information. In particular, our method allows us to observe the behavior of the wave function over the course of the dissociation process. We will do so here by looking at the reduced probability density function, p(R,6) defined in eq 15. Figure 4 shows a series, evenly spaced in time, of surface and contour plots of p(R,B) for the J = 0 2 bend excited state of NeHF. The reduced density is dominated at small R by two lobes at 6 = 0 and 6 = II ( p = IP,o(cos 6)12, since the wave function is dominated by t h e j = 1, k = 0 channel). These lobes are cut off at 1/50 of the peak value in order to show the interesting dynamics of the probability leaking out into the asymptotic region. At large R, the outgoing probability appears as a sequence of wavelets whose profiles are independent of 8 (since the only exit channel h a s j = k = 0, so that p = IPoo(cos@)I2).In between these two regions, at R 9 au, is a ridge of probability density near B =: 0 (that is, near the N e H F configuration of the complex). This ridge appears to bob up and down as the outgoing wavelets escape from the potential well. Figure 5 shows the shape of the function ICj=l,k=Ol as a function of R, at a fixed time in the decay of the IIc resonance; this am-

+

-

~~

~~

(31) Feit, M.D.:Fleck, J. A.;Steiger, A. J . Compul. Phys. 1982, 47,412.

The Journal of Physical Chemistry, Vol. 96, No. 11, 1992 4207

v.1

,

4.0

8.0

12 0

16.0

20.0

R Figure 5. Amplitude of the channel function in the j = 1, k = 0 channel in the neresonance of NeHF at J = 2. Superimposed is a dashed curve representing the product of the channel packet function C1,l and the square of the Coriolis matrix element coupling the two channels, J(J +

1)/R2. plitude is observed to decay with the same lifetime as the overall resonance but maintains its shape throughout. Superimposed over this shape is an appropriately normalized curve obtained by multiplying the amplitude in the j = k = 1 channel (which dominates the wave function) by the Coriolis term which couples the two channels. The curves are nearly identical. Since this Coriolis term is much smaller than the V l ( R )term that couples to the exit channel, we can see the Coriolis coupling as limiting the rate of decay of the resonance. Figure 6 shows another series of snapshots of the probability density, this time for the dissociation of the lowest-energy odd reflection symmetry state in NeDF which suffers internal rotational predissociation. The resonance wave function is dominated by t h e j = 2, k = 1 channel, and the only available exit channel is j = k = 1. The escaping lump of probability density appears to weave back and forth in 6 as its escapes to large R. Figure 7 shows a single snapshot of the reduced probability density in the dissociation of the Ae state, (The product state distribution in this state is discussed in section 1V.B.) As in the NeDF IIrstate, wavelets of probability appear to weave back and forth about 0 = 7r/2 as they fly off to large R. In this case, however, the wave function in the asymptotic region is not constrained by symmetry to vanish at 6 = 0 and 6 = 7r; consequently, there is substantial amplitude throughout the entire range of 6. A notable feature of Figure 7 is the oscillation of the wings; that is, the amplitude at the ends of the range in 6 appears to oscillate between one end and the other. V. Discussion As we have seen, this time-dependent method has been effectively applied to study the quantum dynamics of internal rotational predissociation in an atom-plus-rigid-rotor complex. We are able to calculate dissociation lifetimes, ranging from 50 ps to 10 ns, which agree extremely well with the results of proven time-independent techniques; this suggests that we are able to solve the quantum mechanical problem with similar accuracy. As we noted earlier, this time-dependent method is complementary to timeindependent techniques, in the sense that we can calculate lifetimes with greater precision, but our calculated resonance energies are less accurate; the two kinds of calculation require similar levels of computational effort. If our only purpose were to compute resonance energies and lifetimes, we would have gained little by extending the CCWP formalism to this set of problems. Again, as we noted earlier, the real value of this technique is its ability to reveal information about the dynamics of the dissociation process. The remainder of this discussion will be devoted to exegesis of the dynamical pictures offered in section IV. First, consider the dissociation of the L: resonance, for which time snapshots of p(R,6) are shown in Figure 4. Besides the outgoing spherical wavelets of probability escaping into the as-

4208

Hovingh and Parson

The Journal of Physical Chemistry, Vol. 96, No. I I , I992 t = 241.88 ps I-

Od

-

0.6

-

0.4

-

P

oa -

-.

a-.

4

i

i

io

in

i4

io

io

nip

i4

ti

is

io

l4

ii

ii

io

i4

io

in

nip

R

t

= 242.36ps

1-

0.8 -

0.6

-

0.4

-

oa

-

P

O T

-i

i

i

IO

ip

R

~~

t = 242.85 ps 1

011

0.6

P 0.4

03

-i

i

i

io

in R

t = 243.33 ps I-

0.0-

0.6

-

0.4

-

oa

-

P

0

74

i

i

io

h R

Figure 6. Snapshots of the reduced probability density in the dissociation of the lowest-energy, rotationally predissociating, odd reflection symmetry state in NeDF. The resonance is dominated by the j = 2, k = 1 channel; the only available channel is j = k = 1.

The Journal of Physical Chemistry, Vol. 96, No. 11, 1992 4209

Internal Rotational Predissociation in NeHF and NeDF

0.0500-

0.0375-

0

0.0250-

.r(

2

0.0125 -

0.0

ri

2

d

y

38

Figure 7. Snapshot of the reduced probability density in the dissociation of the even reflection symmetry state dominated by j = k = 2 (Ae) in NeHF at J = 5. The wave function at large R is a linear combination of the j = 1 and j = 0 channels, with time-varying relative phase; thus, amplitude appears to oscillate back and forth between the two ends of the diatom.

ymptotic region, the most obvious dynamical element of this picture is the ridge of probability density at R * 8 au which appears to move periodically up and down, once for each time an outgoing wavelet escapes. We emphasize that the appearance of this feature is independent of the exact nature of the starting wave function, and it persists for many dissociation lifetimes. The profile of the ridge in 8 (large near 8 = 0 and small, but nonzero, near 8 = T ; that is, a sum of P,(cos 0 ) and Po(cos 8)) suggests that we are seeing evidence of the mechanism that couples the j = 1, k = 0 channel, which dominates the resonance, to the (isotropic) product channel. In more physical language, it is evidence of a dissociation mechanism in which a wavelet of probability is allowed to escape each time the amplitude for the Ne-H-F linear configuration peaks, just outside the range of the attractive well. In other words, the proton acts like the end of a hockey stick that, upon striking the neon atom “puck”, sends out a small wave of escaping amplitude. In the case of the dissociation of the nebend states of NeHF, the time-dependent formulation allows us to elucidate the dependence of the predissociation line widths on total J. We find that the dissociation can be viewed as a two-step process: first, amplitude must be drained from the dominant (j = lkl = 1) channel into the intermediatej = 1, k = 0 channel by Coriolis coupling, which depends on J as [ J ( J + 1)]1/2;then amplitude can escape into the asymptotic region by the same process as in the 2 dissociation discussed above. Of these two steps, the Coriolis step is “rate limiting”, and the total rate of escape of probability is proportional to the product of the I: dissociation rate (which is weakly dependent on J) and the square of the Coriolis term. This relationship is clearly seen in Figure 8 where the ratio of T to 2 lifetimes is plotted against J(J + 1) for J = 1-5; the direct proportionality of the two lifetimes is evident. Now consider the dissociation of the NeDF resonance dominated by the j = 2, lkl = 1 channel with negative reflection symmetry (n’?, illustrated in Figure 6. (Note that the wave function is constrained by symmetry to vanish at 8 = 0 and 8 = T.) Here, we imagine that the diatom, as it escapes to large R, continues to rotate in a plane oriented at some angle, to the molecular axis; the angle between the diatom axis and the molecular axis will then oscillate between ~ / +2{and 7r/2 - In fact, we observe that the peak of each outgoing wavelet oscillates about 7r/2 having originated near the small-8 end of the diatom. It is as if the interaction of the neOn and the deuteron “kicks” out a wavelet, which then wobbles out to large R as the diatom

r,

r.

6.0

12.0

18.0

24.0

30.0

J(J+l) Figure 8. Ratio of the inverse lifetimes of the Z and IIestates of NeHF (as calculated on the semiempirical potential surface) as a function of J(J + 1). The linear relationship is obvious.

continues to rotate on its way out. We can observe a similar wobbling behavior in the dissociation of the NeHF Ae resonance, illustrated in Figure 7. But we see an additional feature, the oscillation of amplitude between the two linear orientations of the complex (Ne-H-F and Ne-F-H). This can be understood in mathematical terms as arising from the interference of the amplitudes in t h e j = 1 and j = 0 channels. Amplitude spills out in each of these channels but propagates in R with varying relative phase; thus, the amplitude, as a function of 8 at a given value of R, is a linear combination of Poo(cose), Plo(cose), and Pll(cos 8). It should come as no surprise that the square of this quantity is sometimes larger at 8 = 0 and sometimes larger at 8 = T . There is also, however, an appealing physical picture of what is taking place. In the j = k = 0 and j = 1, k = 0 channels, the diatom has substantial amplitude to lie in the plane of the overall rotation of the complex. From the point of view of a diatom lying in this plane, the escaping neon atom’s orbital motion will appear to bring it past each end of the diatom in turn as it flies away; that is, the hydrogen and the fluorine will successively be more likely to be pointing at the neon. The overall picture that emerges from these interpretations is of a time-dependent dissociation process in which probability amplitude “flows” from one channel to another under the influence of two kinds of interaction. First, overall rotation of the complex causes no loss of diatom angular momentum but forces amplitude to flow into channels with different projections of that angular momentum on the molecular axis; that is, Coriolis forces shift the orientation of the diatom’s plane of rotation in the body-fixed frame. Second, the interaction of the diatom with the atom couples channels of different j, leading to the transfer of energy from internal rotation of the diatom into dissociative translation. Altogether, it is clear that a time-dependent formulation of the quantum mechanical problem of unimolecular dissociation gives us a view of the physics of such processes that is less accessible from other methods. It allows a mechanistic interpretation of the dissociation process, while preserving the physical correctness of a quantum formulation. Moreover, such studies can now be performed at a computational cost not much greater than that of time-independent quantum methods. Overall, it appears that the accuracy, flexibility, and efficiency of such methods show great promise for obtaining greater physical insight into dynamical molecular processes. Acknowledgment. It is a pleasure to thank S. K. Gray and C. M. Lovejoy for many stimulating discussions. This research was supported by the donors of the Petroleum Research Fund, administered by the American Chemical Society, and by the National Science Foundation (NSF PHY90- 12244).