Differential Cross Sections and Product Rotational Polarization in A +

Aug 6, 2009 - The state-to-state differential cross sections for some atom + diatom reactions have been calculated using a new wave packet code, MAD-W...
4 downloads 11 Views 5MB Size
14488

J. Phys. Chem. A 2009, 113, 14488–14501

Differential Cross Sections and Product Rotational Polarization in A + BC Reactions Using Wave Packet Methods: H+ + D2 and Li + HF Examples† A. Zanchet,‡ O. Roncero,*,‡ T. Gonza´lez-Lezana,‡ A. Rodrı´guez-Lo´pez,§ A. Aguado,⊥ C. Sanz-Sanz,‡,# and S. Go´mez-Carrasco‡,¶ Instituto de Fı´sica Fundamental, CSIC, Unidad Asociada UAM-CSIC, Serrano 123, 28006 Madrid, Spain, Centro de Supercomputacio´n de Galicia, AV. de Vigo s/n (Campus Sur), 15706 Santiago de Compostela, Spain, and Departamento de Quı´mica Fı´sica, Facultad de Ciencias C-XIV, Unidad Asociada UAM-CSIC, UniVersidad Auto´noma de Madrid, 28049, Madrid, Spain, School of Chemistry, UniVersity of Birmingham, Edbaston, Birmingham B15 2TT, United Kingdom, and Theoretical Chemistry Department, Institute of Physical Chemistry, UniVersity of Heidelberg, Im Neuenheimer Feld 229, D-69120 Heidelberg, Germany ReceiVed: April 28, 2009; ReVised Manuscript ReceiVed: June 24, 2009

The state-to-state differential cross sections for some atom + diatom reactions have been calculated using a new wave packet code, MAD-WAVE3, which is described in some detail and uses either reactant or product Jacobi coordinates along the propagation. In order to show the accuracy and efficiency of the coordinate transformation required when using reactant Jacobi coordinates, as recently proposed [J. Chem. Phys. 2006, 125, 054102], the method is first applied to the H + D2 reaction as a benchmark, for which exact timeindependent calculations are also performed. It is found that the use of reactant coordinates yields accurate results, with a computational effort slightly lower than that when using product coordinates. The H+ + D2 reaction, with the same masses but a much deeper insertion well, is also studied and exhibits a completely different mechanism, a complex-forming one which can be treated by statistical methods. Due to the longer range of the potential, product Jacobi coordinates are more efficient in this case. Differential cross sections for individual final rotational states of the products are obtained based on exact dynamical calculations for some selected total angular momenta, combined with the random phase approximation to save the high computational time required to calculate all partial waves with very long propagations. The results obtained are in excellent agreement with available exact time-independent calculations. Finally, the method is applied to the Li + HF system for which reactant coordinates are very well suited, and quantum differential cross sections are not available. The results are compared with recent quasiclassical simulations and experimental results [J. Chem. Phys. 2005, 122, 244304]. Furthermore, the polarization of the product angular momenta is also analyzed as a function of the scattering angle. I. Introduction The differential cross section (DCS) provides the most detailed information on the reaction dynamics,1-3 even more if its dependence with the rotational polarizations of reactants and products are considered.4-6 Nowadays, there are many theoretical treatments allowing its calculation.7-14 The quasiclassical trajectory (QCT) method is probably the most commonly used7,8 because it is able to describe most of the features of the reactions at a relatively low computational cost. However, the possible importance of quantum effects might preclude classical approaches to describe the dynamics of the process. Until very recently, most quantum simulations of DCSs for reactive collisions were obtained with time-independent (TI) methods based on hyperspherical coordinates.10,15-17 These methods are able to calculate the whole S matrix at once, but the size of the matrices involved grows quadratically with the number of basis set functions, n, required to represent reactants and products. †

Part of the “Vincenzo Aquilanti Festschrift”. * To whom correspondence should be addressed. E-mail: ORoncero@ imaff.cfmac.csic.es. ‡ Instituto de Fı´sica Fundamental, Unidad Asociada UAM-CSIC. § Centro de Supercomputacio´n de Galicia. ⊥ Universidad Auto´noma de Madrid. # University of Birmingham. ¶ University of Heidelberg.

The computational effort, however, increases as n3. Thus, they are limited to systems with few degrees of freedom involved in the reaction, being applied essentially to triatomic and few tetra-atomic systems.18 Wave packet (WP) techniques provide a single column of the S matrix but only need to store vectors and compute matrix/ vector multiplications.19,20 In addition, they provide information over a wide range of energies in a single calculation but require many iterations, especially for slow dynamics occurring in the presence of long-lived resonances. The computational cost involved is high, especially when taking into account that DCSs are typically measured at one or just a few collisional energies. This is why most WP calculations performed until now have been restricted to integral cross sections (ICSs). A recent review21 provides a extensive list of the systems studied with these techniques. Althorpe was the first to calculate an “exact” DCS with WP methods.12 He used the reactant-product decoupling method22,23 to study the H + H2 prototype reaction and isotopic variants.24-26 This method is very efficient because the WP is split into several parts, each of them confined to a single arrangement channel, using the best adapted coordinates. The connection among the different portions of the WP is done through imaginary potentials. This procedure, however, presents problems when

10.1021/jp9038946  2009 American Chemical Society Published on Web 08/06/2009

Cross Sections and Polarization in A + BC Reactions the reaction shows resonances by the presence of wells, which would be strongly distorted by the partition of the configuration space. With the advent of massively parallel computers, it has been possible to perform exact WP calculations using a single set of coordinates for the full propagation.13,14 In most cases, product Jacobi coordinates are used, transforming only the initial WP, avoiding any transformation along the propagation. In these coordinates, the initial WP cannot be factorized. For example, a pure helicity state in reactant coordinates transforms into a superposition of many helicity components in product coordinates when using a body-fixed (BF) frame. The use of reactant Jacobi coordinates to calculate state-tostate reaction probabilities and scattering matrix amplitudes to obtain DCSs has also been proposed.27 In this procedure, the initial WP is represented by a single helicity component. The evaluation of the state-to-state transition amplitudes requires, however, the transformation from reactant to product Jacobi coordinates at each iteration. This is made very efficiently because only a single R′ value (R∞′ ) is required to evaluate the flux28 and because the coordinate transformation is divided into several steps, each of them requiring four loops at most. The computational cost is then only a fraction of the action of the Hamiltonian on the WP. It has been argued, however, that such a transformation could not be accurate enough to produce DCSs.14 The aim of this work is to show that such a transformation method can provide excellent results and to apply it to the study of Li + HF. In ref 27, the comparative efficiency of using reactant or product Jacobi coordinates was discussed in terms of the masses involved; for H + H′L f HH′ + L collisions (where H/L denotes heavy/light atoms), reactant Jacobi coordinates were considered to be the best adapted ones, while product Jacobi coordinates were very well suited for H + LH′ f HL + H′. Here, we shall first check the reliability of the method by applying it to the intermediate case of H + D2 f HD + D as a prototype of direct reaction, where reactant coordinates are slightly more efficient. The insertion reaction H+ + D2 f HD + D+, however, involves the same masses, but the potential presents a longer range because of its ionic character, and product Jacobi coordinates become more efficient. Calculations for some selected total angular momenta are shown, and the state-selected DCSs are obtained within the random phase approximation. The calculated DCSs are compared with recent exact TI calculations29 and experimental results.30 A more detailed study is applied to the Li + HF f LiF(V′, j′) + H reaction. This system constitutes a benchmark for reactive scattering calculations since it is one of the lightest reaction involving three different atoms, allowing quite reliable potential energy surface (PES) calculations, for the ground31-34 and excited states.35,36 Experimentally, there are collisional37-41 as well as transition-state spectroscopy42 studies. Simulations of spectroscopic processes give detailed information about the electronic curve crossings giving rise to reaction barriers; infrared excitation of the complex between reactants43,44 can be used to probe the ground electronic state, electronic spectroscopy probes the upper electronic states,35,36,45-47 and electronic predissociation studies give information on the nonadiabatic couplings.47-50 Theoretically, reactive collisions are the most studied processes in this system, including QCT,41,51-53 quantum TI,32,54-58 as well as time-dependent WP59-64 calculations. Lately,thesestudieshavebeenextendedtoultracoldtemperatures.65,66 The large electric dipoles of the reactants and products make this system ideal to study its stereodynamics experimentally,38,67-69

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14489 motivating many theoretical simulations.6,70-74 Despite all these studies, the DCSs have only been calculated using the QCT approach.41 The existence of many resonances and zero-point effects in Li + HF collisions suggests that quantum mechanical simulations of DCSs are the most adequate to mimic experimental results. One of the main goals of this work is to calculate quantum DCSs for this reaction and study the stereodynamics of this reaction by analyzing the vector correlations on this system. The paper is structured as follows. Section II is devoted to describe the method. In the first part of section III, the method is applied to H + D2 and H+ + D2 reactive collisions to check the accuracy compared with TI methods, discussing the relative efficiency of using the reactants or products case, for such mass combination. The second part of section III, describes the results obtained for Li + HF f LiF + H, well suited for reactant Jacobi coordinates, making emphasis on the final polarization of LiF products and the kk′j correlations. Finally, section IV is devoted to extract some conclusions. II. Wave Packet Method The DCS for AB(V, j) + C f A + BC(V′, j′) reactive collisions, at an energy E, is commonly expressed in a BF frame as9,75,76

∂σV,jfV',j'(E) 1 |T (E, Θ)| 2 ) ∂Θ (2j + 1) Ω,Ω′ V,j,ΩfV',j',Ω′



TV,j,ΩfV',j',Ω′(E, Θ) )

1 2kVj

J J (E)dΩΩ′ (Θ) ∑ (2J + 1)SV,j,ΩfV',j',Ω′ J

(1) where kVj ) (2µ(E - EVj)/p2)1/2 is the wave vector in the entrance J (Θ) are the reduced Wigner channel of energy EVj and dΩΩ′ rotation matrices77 depending on the center of mass (CM) scattering angle Θ. J is the total angular momentum, with projections Ω and Ω′ on the reactant and product BF z-axes, respectively. Hereafter, bold letters are reserved for vectors, quoted variables and quantum numbers refer to products, while nonquoted refer to reactants. The definitions given below for the case of reactant coordinates are equivalent to those for products. The three atoms are in the x-z BF plane, with the z-axis parallel to the vector R joining the CM of the diatomic molecule AB to the third atom. It may be distinguished between internal degrees of freedom R, r, and γ, with r being the AB internuclear vector, and cos γ ) R · r/Rr, and the three Euler angles φ, θ, and χ defining the rotation from the space-fixed (SF) to the BF frames. The m and µ are the reduced masses associated with r and R, respectively. DCS in eq 1 provides direct information about the kk′ correlation. However, it corresponds to an isotropic distribution of the initial state, and the sum over final Ω′ washes out the effects of the initial and final polarization of diatomic rotation on the reactivity. In a BF frame, Ω (Ω′) is the projection of the diatomic angular momentum j of the reactant (j′ of the product) along the direction R (R′), which coincides asymptotically with Vj,V′ j′ k (k′). Thus, the matrix defined as FΩ1Ω2,Ω′1Ω′2(E, Θ) ) TV,j,Ω1fV′,j′,Ω′1(E, Θ)TV,j,Ω2fV′,j′,Ω′2(E, Θ) provides direct information about the effect of the polarization of j on the reactivity and the final polarization of products after the reaction. Since the Vj,V′ j′ number of projections is high, FΩ1Ω2,Ω′1Ω′2(E, Θ) is typically 4,5,78 The analysis can be simplified expanded in state multipoles.

14490

J. Phys. Chem. A, Vol. 113, No. 52, 2009

Zanchet et al.

by defining the orientation (O) and alignment (A ) of the angular momenta of reactants (R) and products (P) as

R,P OVj,V'j' (Θ, E) )



ΩR,P R,P AVj,V'j' (E, Θ) )



ΩR,P R,P

Ω (Θ, E) ) UVj,V'j'

(

ΩR,P

√jR,P(jR,P + 1)

(

)

R,P

Ω UVj,V'j' (Θ, E)

)

3(ΩR,P)2 ΩR,P - 1 UVj,V'j' (Θ, E) R,P R,P j (j + 1)

Vj,V'j' (E, Θ) ∑ |TV,j,ΩfV',j',Ω′(E, Θ)|2/ ∑ FΩΩ,Ω′Ω′

ΩP,R

ΩΩ′

(2) where (jR, ΩR) ≡ (j, Ω) and (jP, ΩP) ≡ (j′, Ω′). For high j, the semiclassical limit can be used, in which cos θj ) Ω/[j(j + 1)]1/2, θj being the angle between j and the z-axis (and similarly for products). Thus, when the orientation OR approaches +1 (or -1), the most probable event corresponds to j pointing in the same (opposite) direction as k. If the quantization axis is taken along k and/or an initial isotropic distribution is assumed, the orientation and all odd terms of the multipole expansion vanish. In photodissociation, the asymmetry is introduced by circularly polarized light, so that integral and differential78 cross sections can be different from 0. The alignment A R, however, is generally nonzero and is the primary quantity determined from polarization of laser-induced fluorescence. It takes the limiting values of -1 and 2, depending on whether the reaction is more probable when j is perpendicular or parallel to k, respectively. Such limits provide a simple procedure to analyze the correlation among the different vector magnitudes. J (Θ) f 0 when Θ f 0 or π, In addition, note that since dΩΩ′ R,P there is only Ω ) 0 contributing at such angular values. Thus, A ) -1 at Θ ) 0, π. As an example, in ref 62, this procedure was used with the total ICS of the Li + HF(V ) 0,1 and j ) 0-3) reactive collisions. Important steric effects were found for V ) 0, showing a clear preference for collinear collisions in which the HF reactant gets vibrationally excited more easily, as a required condition for the reaction to occur. For V ) 1, however, the alignment was close to 0, thus showing no clear preference since the vibrational excitation of HF was already enough to activate the reaction. An alternative study was done on this reaction by Alvarino et al.71,72 using the stereodirected representation of Aquilanti et al.79 In such studies the S matrix for a particular J value was transformed from the V, j representation to a DVR-like one, giving information about the attack and recoil angles. Here, we shall focus our attention on the kk′j′ correlations by studying the polarization of product angular momenta as a function of the scattering angle. The main quantities in eq 1 are the S matrix elements, J (E), which provide all of the information about the SV,j,ΩfV′,j′,Ω′ collision process. In practice, it is more efficient to calculate the S matrix in a parity-adapted basis set and then transform it to the S matrix whose elements are used in eq 1 through the expressions

where f ) [(1 + δ0Ω)(1 + δ0Ω′)]1/2/2, and ( are used to denote the parity under spatial inversion,  ) (1. The calculation of SJ using WP methods in product Jacobi coordinates has been described previously.13,14 Here, it will be described using reactant Jacobi coordinates, as proposed recently,27 making emphasis on the differences. This method has been implemented in the code MADWAVE3. This code uses either reactant or product Jacobi coordinates to obtain state-to-state scattering matrix elements, and it also treats photoinitiated processes. In addition, it considers several coupled electronic states described in a diabatic representation. The code is parallelized using the MPI library, with respect to the helicity components and also the angular grid points. A. Propagation. One of the most commonly used timedependent integrators is the Chebyshev real WP propagator, or slight variants.80-85 In this method, the evolution operator is expanded in modified Chebyshev polynomials,81 so that the WP at time t is defined as

ΨR(t) )

J J SV,j,ΩfV',j',-Ω′ (E) ) SV,j,-ΩfV',j',Ω′ (E) ) J+ J- (-1)JSV,j,ΩfV',j',Ω′ ) f(SV,j,ΩfV',j',Ω′

(4)

k

where R denote the quantum numbers required to specify the state of the reactants. ΨR(k) are obtained using the modified Chebyshev recurrence86

ΨR(k ) 0) ) ΨR(t ) 0) ˆ sΨR(k ) 0) ΨR(k ) 1) ) e-φH -φ ˆ sΨR(k) - e-φΨR(k - 1)} ΨR(k + 1) ) e {2H

(5)

ˆ s, t) in eq 4, are given by87 The coefficients fk(H

ˆ s, t) ) (2 - δk0)e-iE0t/p(-i)kJk(t∆/p) fk(H

(6)

ˆ s ) (H ˆ where Jk(x) are Bessel functions of the first kind. H E0)/∆ is the scaled Hamiltonian, with E0 ) (Emax + Emin)/2 and ∆ ) (Emax - Emin)/2. Emax and Emin are the minimum and maximum energy values of the real Hamiltonian of the system, ˆ , represented in the grid and basis set used. H The real absorbing function, φ, in eq 5, depends on the dissociative coordinates and serves to avoid the boundary problems arising from the use of finite grids, such as reflection and transmission. Typically, φ ) -Ax(x - xI)n for x > xI and 1 elsewhere (with x ) r or R). The absorption parameters xI and Ax are optimized to minimize the reflection and transmission effects that occur when the WP reaches the edges of the grid, leaving nearly unaffected the portion of the WP in the internal region.88 To extract the magnitudes for a fixed total energy, the TI eigenfunctions of the total Hamiltonian, Ψ+R (E), with total energy E, are obtained in terms of the ΨR(k) modified Chebyshev components as89

J J SV,j,ΩfV',j',Ω′ (E) ) SV,j,-ΩfV',j',-Ω′ (E) ) J+ J+ SV,j,ΩfV',j',Ω′ ) f(SV,j,ΩfV',j',Ω′

∑ fk(Hˆs, t)ΨR(k)

Ψ+ R (E) )

(3) with





1 ˆ , E)ΨR(k) c (H 2πpaR(E) k)0 k s

(7)

Cross Sections and Polarization in A + BC Reactions

ˆ s, E) ) ck(H

(2 - δk0)pe-ik arccos(E-E0)/∆

√∆2 - (E - E0)2

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14491

(8)

〈r, R, γ|ΦJ Ω (k)〉 J ) WMΩ(φ, θ, χ) rR Ωg0



(9) where the 〈r,R,γ|ΦJ Ω (k)〉 coefficients depending on the internal degrees of freedom, r, R, and γ, are propagated numerically when eq 9 is introduced into eq 5. The radial variables r and R are described in finite grids of equidistant points, while γ is described by Gauss-Legendre quadrature points. For the Euler angles, a parity-adapted basis set is used, of the form

J WMΩ (φ, θ, χ) )



2J + 1 J* (DMΩ (φ, θ, χ) + 16π2(1 + δ0,Ω) J* (-1)J+ΩDMΩ (φ, θ, χ))

(10)

J* are Wigner rotation matrices.77 When the system where DMΩ presents open-shell fragments and several electronic states participate, it is convenient to include the electronic wave function in the angular basis set to build more general symmetryadapted functions.90 Since we are applying the method for reactive scattering on the ground electronic surface, we shall omit here the electronic state for simplicity. When the Hamiltonian operator is applied to the functions of eq 9, a set of coupled differential equations is obtained.91 The radial kinetic terms are solved using the fast Fourier transform method. The action of the angular momentum operators is performed by multiplying the WP by a matrix91 as

J j2〈Rk, rl, γn |ΦΩ (k)〉 )

∑A

J J J 〈WMΩ′ |l 2 |WMΩ 〉〈Rk, rl, γn |ΦΩ (k)〉

)



l J An'Ω′,nΩ 〈Rk, rl, γn |ΦΩ (k)〉

n'

(11) where j and l are the angular momenta associated with r and R, respectively. In the above expression the A matrices are defined as

∑ √w′ Y

n jΩ(γ′n, 0)j(j

j

l An'Ω′,nΩ )

∑ √w′ Y

+ 1)√wnYjΩ(γn, 0)

√wnYjΩ(γn, 0)

2 J J n jΩ′(γ′n, 0)〈WMΩ′YjΩ′ |l |WMΩYjΩ〉

j

(12) with

The only couplings between adjacent Ω components are due to the nondiagonal terms of l 2. This calculation is very well suited for massively parallel computing;91,92 a Ω value is assigned to each processor, and the exchange of information is only among adjacent processors because of the band structure of the l 2 matrix elements. The present implementation of the method parallelizes with respect to Ω and also on γ with the MPI library, becoming very efficient even for a large number of Ω components. The angular momentum terms are the only ones coupling different γn points, and unfortunately, the corresponding matrices do not have so many zeros. Therefore, the exchange of information is among all processors involved, and the efficiency gets lower. However, for low J, it is convenient to reduce computing times by increasing the number of processors. In the case that R approaches 0, l 2/2µR2 becomes very large. It is necessary to limit it to some finite value, so that Emax does not become too large. Several procedures have been proposed. J J YjΩ| l 2|WMΩ′ Yj′Ω′〉 matrix Lin and Guo14 diagonalized the 〈WMΩ to limit the eigenvalues to Ecut and then transform back to the BF frame. The resulting matrix was completely full, which involves information exchange among all processors having different Ω values. Thus, such a procedure would become inefficient for massive parallelization. The second method13 consists of limiting the diagonal terms of the J J YjΩ|l 2|WMΩ′ Yj′Ω′〉 matrix (in the BF representation) to Ecut, 〈WMΩ setting to 0 the couplings corresponding to the modified diagonal terms with the rest of the Ω’s. This procedure assures the band structure of the matrix, without affecting the parallelization, but introduces some instabilities for high J’s. In this work, we follow a different strategy. A mean value is defined as the average of all eigenvalues in the J J YjΩ|l 2|WMΩ′ Yj′Ω′〉 matrix as a function of R. When it 〈WMΩ becomes equal to a given quantity Ecut, the value of Rcut is set. Thus, for R < Rcut the value of R is substituted by Rcut for a particular J, j pair. In this way, the band structure of the matrix is maintained and the instabilities disappear, provided that Ecut is sufficiently high. C. Initial Wave Packet. In reactant Jacobi coordinates, the initial state of the fragments, R ) J, , V, j, and Ω0, at a sufficiently long distance is taken as

j J n',n〈Rk, rl, γn |ΦΩ (k)〉

n'

j An',n )

δΩ,Ω′(1√1 + δΩ,0 + δΩ′,0√J(J + 1) - ΩΩ′√j(j + 1) - ΩΩ′

(13)

B. Wave Packet Representation. The WP is expressed in Jacobi coordinates, reactants or products, in a BF frame as defined above. The same representation is valid for WP, ΨR(t), TI wave functions, ΨR(E), and Chebyshev components, ΨR(k), and we shall describe it only for the last one, ΨR(k). In a partial wave expansion, these coefficients are calculated for a particular angular momentum J and a given parity under total inversion, , separately, each one being represented as

〈r, R|ΨJM R (k)〉

J J 〈WMΩ YjΩ |l 2 |WMΩ′ Yj'Ω′〉 ) δΩ,Ω′[J(J + 1) + j(j + 1) - 2Ω2] -

ΨR(t ) 0) )

∑W

J J,,V,j,Ω0 (t MΩ(φ, θ, χ)χVj(r)YjΩ(γ)〈R|gΩ

) 0)〉



(14) where χVj(r) are the vibrational eigenfunctions of the AB(V, j) reactant and YjΩ(γ) are normalized associated Legendre func2 2 0(t ) 0) ) δ tions. Initially, gJ,,V,j,Ω Ω Ω,Ω0 exp[-(R - R0) /2Γ 2 1/4 iK0R]/(πΓ ) , that is, an incoming complex Gaussian. The asymptotic distance, where the initial WP is placed, corresponds to that distance where the potential between the two reactants is 0. This implies the use of larger grids for J > 0 since the centrifugal term varies as R-2. In order to avoid this problem, another alternative is followed, which consists of the following: (1) Backward propagatation in time, up to a time at which the average centrifugal energy is negligible. The propagation is performed using a split operator propagator on a single radial

14492

J. Phys. Chem. A, Vol. 113, No. 52, 2009

Zanchet et al.

grid (usually large) for a particular J, , V, j state, including all possible Ω channels. (2) After this propagation, Coriolis couplings produce a WP with several Ω components. All of the components of the WP with Ω * Ω0 are set to 0 to start with a proper Ω0 value. The WP is then renormalized, and the incident flux is evaluated as

aR(E) )

1 2i



µ 2πp2kVj

∫ dR eik R〈R|gΩJVjΩ (-tmax)〉 Vj

0

0

(15)

(3) The WP is propagated forward in time for the same number of steps, so that at the end, the WP is essentially placed at the same position as the initial Gaussian. (4) To start the propagation for the reactive collision, the WP at time zero is reconstructed according to eq 14, with the new 0(t ) 0) obtained in this procedure. To make the initial gJ,,V,j,Ω Ω 0(t ) 0) is added WP real, the complex conjugate of each gJ,,V,j,Ω Ω and renormalized. Thus, ΦJ Ω (k) Chebyshev components remain real along the propagation.81,83-85 When using reactant Jacobi coordinates, after this procedure, the propagation starts. When using product Jacobi coordinates, however, the initial WP is transformed from reactant to product coordinates. As a result, the initial WP expressed in product Jacobi coordinates is a superposition of different Ω′ values, and each component is given by

JR 〈r', R', γ′|ΦΩ′ (k ) 0)〉 ) Ω′

(-1)

J [dΩΩ′ (β)

rR r'R'

∑ Ω



1 + δΩ′0 × 1 + δΩ,0

J + (-1)J+Ωd-ΩΩ′ (β)] ×

〈r, R, γ|ΦJR Ω (k ) 0)〉

|

(16)

Jj J Jj' SV,j,lfV',j',l ′(E)TΩ′,l ∑ ∑ (-1)j+l+j'+l ′TΩ,l ′ l

l′

(17)



k

Ω′

(20) This last quantity is obtained as a transformation from the Chebyshev coefficients to the time-independent wave function, eq 7, followed by a BF-to-SF frame transformation. The coefficients CV′j′Ω′(k) in the above equation are simply the overlap between the WP and the final states of the products at R′ ) R∞′ in product Jacobi coordinates, that is

CV'j'Ω′(k) )

∫ sin γ′dr'dγ′χV'j'(r')Yj',Ω′(γ′, 0)〈r', R′∞, γ′|ΦΩ′J (k)〉

(21)

where the indexes referring to the initial state have been suppressed for clarity. These coefficients are directly obtained by numerical integration when using product coordinates in the propagation. When using reactant coordinates, however, a reactants-to-products transformation must be done at each iteration k which might be very demanding computationally. An efficient method for this transformation has been proposed recently,27 requiring much less effort than the application of HΨ(k). In such a method, the transformation is not made in a single step, as in eq 16, but it is divided into the following three steps: (1) The asymptotic states of the products, χV′j′(r′)Yj′Ω′(γ′), are represented in the mixed coordinates (R′∞, R, γ) (according to possibility A of ref 27). This is done once before starting the propagation, and the transformed functions are kept in a two-dimensional grid as small as possible to save memory storage and computation time. (2) At each iteration k, the Chebyshev component is transformed to an intermediate system of coordinates, that is J 〈r, R, γ|ΦJ Ω (k)〉 f 〈R′∞, R, γ|ΦΩ (k)〉

(22)

The transformation matrix, given in eq 16 of ref 27, required for this transformation is constructed before the propagation and stored. (3) The transformation from reactant to product BF frames is performed according to J 〈R′∞, R, γ|ΦΩ′ (k)〉 )

where Jj TΩ,l

Jj Jj' ˆ s, E)CV'j'Ω′(k)]TΩ′,l [ ∑ ck(H ∑ ∑ TΩ,l ′

(r,R,γ)≡(r',R',γ′)

where cos β ) R · R′/RR′. D. S Matrix Elements. Parity-adapted S matrix elements are obtained by a transformation from the SF frame as J SV,j,ΩfV',j',Ω′ )

J QV,j,l fV',j',l ′(E) )

∑ TΩ′,Ω〈R′∞, R, γ|ΦJΩ(k)〉

Ωg0

(

l j J [1 + (-1)j+l] ) (-1)j-l-Ω√2l + 1 0 Ω -Ω 2

)

(18)

are the SF-to-BF transformation matrix elements. The S matrix in the SF reference system is obtained from the flux of the TI eigenfunctions at the asymptote93 as J SV,j,l fV',j',l ′(E) ) -



i

2kV'j' 1 e-il ′π/2 QJ (E) (19) πµ′ aVjΩ(E) h(2)(k R′ ) V,j,l fV',j',l ′ l V'j' ∞

where hl(2)(x) are spherical Bessel functions of the third kind and

(23) where the T matrix is given by eq 23 of ref 27. (4) Finally, the overlaps CV′j′Ω′(k) are calculated by numerical integration. III. Results and Discussion The WP simulations have been performed using the following PESs: the BKMP2 by Boothroyd et al.94,95 for H + D2 f HD + D, the ARTSP by Aguado et al.96 for H+ + D2 f HD + D+, and the APW by Aguado et al.97 for Li + HF f LiF + H. The parameters used in the WP calculations are listed in Table 1 for the three systems, indicating if reactant or product Jacobi coordinates were used. For the case of the H + D2 reaction, TI calculations with hyperspherical coordinates have been performed for comparison using the ABC code.11 In this last case, the following parameters have been used: the maximum

Cross Sections and Polarization in A + BC Reactions

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14493

TABLE 1: Parameters Used in the WP Calculations in Reactant and/or Product Jacobi Coordinates for the Three Reactions Studied

reactants coord.

products coord.

Li + HF f LiF + H reactants coord.

0.1 8 128 6 0.5 Å-2 0.1 8 128 6 0.5 Å-2 36 3.2 0.6 0.2 5.6 2.5 6

0.3 14 256 11 0.002 Å-4 0.01 18 512 14 0.002 Å-4 70 8.5 0.25 0.12 15 2.5 8

0.25 17.5 400 13 0.017 Å-2 0.5 17.5 512 14.5 0.025 Å-2 50 13 0.2 0.07 10 3 7

H + D2 f HD + D H+ + D2 f HD + D+

rmin/Å rmax/Å Nr rI/Å AR Rmin/Å Rmax/Å NR RI/Å AR Nγ R0/Å E0/eV ∆E/eV R∞′ Vcut/eV l Ecut /eV

rotational quantum number is jmax ) 50, the helicity truncation parameter is Ωmax ) 10, the maximum hyperradius is Fmax ) 12.0 Bohr, the number of log-derivative propagation sectors is 110, and the maximum internal energy in any channel is Emax ) 2.0 eV. A. H + D2 Collision in Reactant Coordinates. This reaction has been studied intensively both experimental and theoretically (see, for example, ref 98), and here, it is used as a benchmark. Calculations for all angular momenta, from J ) 0 up to 25, were done for the H + D2(V ) 0, j ) 0) f HD(V′, j′) + D reaction using the ABC code and the MAD-WAVE3 code with reactant Jacobi coordinates, giving rise to the DCSs shown in Figure 1. Similar calculations were carried out recently by Chu et al.99 and are used here as a benchmark to check the accuracy of the method. The “exact” reactant wave packet (reac-WP) calculations presented here include all possible Ω projections and required 6000 iterations to be converged. The state-to-state DCSs obtained are compared with TI results, showing an excellent agreement. Thus, the use of reactants coordinate, undergoing the coordinate transformation at each iteration, yields accurate results. The DCSs are perfectly converged when the number of helicity components is limited to Ωmax ) 11. The same calculations are carried out but using product Jacobi coordinates (prod-WP), and the parameters used are very similar to those reported previously by Hankel et al.13 The number of angular grid points required (∼50) is larger than that in the reacWP case. A higher number is needed because the permutation symmetry of the D2 reactants is not properly taken into account when using product coordinates. In product coordinates, the initial WPs correspond to a superposition of many Ω′, while in reactant coordinates, there is only one Ω. In fact, this number increases with the position of the initial WP, determined by R0. The Coriolis term couples different Ω′, and the initial distribution in the helicity components is expected to become even broader along the propagation. However, the DCSs do not seem to be too sensitive to such a situation. In fact, the DCSs obtained ′ similar in prod-WP show a degree of convergence with Ωmax to that obtained using the reac-WP method. A comparative analysis of the efficiency of both WP approaches reveals that the method based on reactant coordinates yields slightly cheaper calculations, especially due to the considerable reduction in the number of angular grid points. The number of helicity components required to converge the

Figure 1. State-to-state DCSs obtained with TI and reac-WP calculations for H + D2(V ) 0, j ) 0) f HD(V′ ) 0, j′) + D, for Etrans ) 1.0083 eV.

DCSs is slightly smaller when using product coordinates, but the reduction is so small that it becomes nearly insignificant computationally. The alignment of rotational angular momenta of HD products, shown in Figure 2, is -1 at Θ ) 0 and π, as commented above. For other scattering angles, however, it changes with the final j′ in a nontrivial way, varying notably with the initial translational energy as well. Therefore, the analysis of A can be used as a sensitive magnitude for the description of such a direct reaction. B. H+ + D2 Collision in Product Coordinates. Part of the advantages of using reactant Jacobi coordinates are lost for longrange potential interactions; when R∞ becomes large, dense grids have to be used to represent the product wave functions. This is the case of the ionic H+ + D2(V ) 0, j ) 0) f HD(V′, j′) + D+, which follows an insertion mechanism. The dynamics has been described in detail previously.29 It is found that nearly all helicity components have to be included in the calculations to reproduce correctly the exact results obtained with the TI method.100 Until now, the WP calculations performed on this reaction beyond the centrifugal sudden (CS) approach101,102 have only included very few helicity components (up to Ωmax ) 5), which cannot guarantee converged results. In this work, we use product Jacobi coordinates including all Ω projections. A very large number of iterations is needed to converge the results, ∼105, and the calculation of all partial waves becomes very demanding. As in the previous case, the use of reactant coordinates is slightly more efficient than that of products but only for low J because of the reduction in the number of angular grid points when considering the permutation symmetry. In fact, reac-WP and prod-WP calculations performed up to J ) 30 (not shown

14494

J. Phys. Chem. A, Vol. 113, No. 52, 2009

Zanchet et al.

Figure 2. Alignment parameter obtained using eq 2 from reac-WP calculations for H + D2(V ) 0, j ) 0) f HD(V′ ) 0, j′) + D, for several Etrans.

here) are in excellent agreement. For higher J, however, the situation changes completely. Since the reaction is essentially statistical, the higher j′ are obtained with the highest probability, and therefore, also high Ω′ are very probable. In reactant coordinates, at long R∞′ where the transformation to products is done, the number of angular grid points is reduced in the energetically accessible region because the skewing angle is rather small. Moreover, the accessible angular points approach the linear configuration. Since ΦΩ ∝ sinΩ γ, for γ ) 0 and π, this situation implies that high Ω are very poorly described unless huge angular grids are considered. Those Ω . 0 components are, on the other hand, of great importance given the enormous effect of the Coriolis coupling in such an insertion reaction with such low masses. A reac-WP calculation would require increasingly large angular grids as J increases, with the corresponding computational difficulty. For this reason, product coordinates are more efficient for this kind of reactions. Since this insertion reaction is known to be rather well described by statistical methods, the random phase approximation is expected to work very well.103 In this approximation, the phase of the S matrix is considered to vary very rapidly with J due to the presence of a high number of resonances. As a consequence, the interference term appearing in eq 1 may be neglected, obtaining a good description of the exact DCS. For this reason, only state-to-state reaction probabilities are needed, not the amplitudes. These probabilities are calculated for some J values, 0, 5, 10, ..., 45, and 50 in this case, including all possible Ω′ components. For those J’s which were not directly calculated, the probabilities are interpolated using a variation of the J-shifting approach. The total reaction probabilities obtained for some J’s are shown in Figure 3, showing a highly

Figure 3. Total reaction probabilities obtained using the prod-WP method, including all helicity components Ω′, for the H+ + D2(V ) 0, j ) 0) f HD + D+ reaction, for different total angular momentum J, as a function of the translational energy. SQM results, from ref 29, have been included for comparison.

oscillating shape due to the presence of many resonances. The present prod-WP results are compared with probabilities obtained previously29 by means of a statistical quantum method103-105 (SQM), which, despite not describing the resonant structure, provides a good average reproduction of exact TI results. For high J, the reaction probabilities at the threshold show a maximum, which becomes more significant as J increases. This structure corresponds to the high J maximum of the opacities studied in ref 29 with exact TI and statistical methods. The total reaction cross sections obtained with the prod-WP and the SQM methods are compared in the top panel of Figure 4. Due to the higher zero-point energy of HD products (∼0.05 eV), the cross section shows a sudden increase at this threshold and then a slow decrease. In general, the SQM reaction probabilities are higher than those obtained with the prod-WP method, except for intermediate energies and low total angular momenta. For the two exact TI points (EQM in the figure), a similar situation holds. However, for the lower energy of 0.1 eV, the TI result is higher than the prod-WP result. The very oscillating shapes shown in Figure 3 associated with very longlived resonances are difficult to converge completely, especially at low energies, since it implies very long propagations, perhaps longer than those performed here. The numerical problems associated with this can explain the underestimation of the total

Cross Sections and Polarization in A + BC Reactions

Figure 4. Reaction probabilities as a function of total angular momentum J for the H+ + D2(V ) 0, j ) 0) f HD(V′, j′) + D+ reaction at a translational energy of 524 meV summed over all V′, j′ channels (bottom panel) and ICS for V′ ) 0, j′ (middle panel). Prod-WP results are compared with the exact TI results of ref 29. ICSs for two close energies, 520 and 530 meV, calculated with the prod-WP method are also included for comparison. In the top panel, the total integral cross sections for the H+ + D2(V ) j ) 0) as a function of translational energy obtained with the prod-WP, SQM, and EQM29 methods are compared.

cross section obtained at low energies for the prod-WP. At 0.524 eV, however, the prod-WP results are in excellent agreement with respect to the exact TI results and are a significant improvement with respect to the adiabatic centrifugal sudden approach results reported in ref 29, which yielded a cross section of 11.6 Å2. The total reaction probability as a function of angular momentum for a translational energy of 0.524 eV, in the bottom panel of Figure 4, shows the oscillating structure associated with the resonances appearing in the reaction for all J’s. Since those resonances correspond to very excited quasi-bound states, very dense grids are required to describe precisely their nodal structure. It is difficult to achieve the description of the nodes of all of the resonant structures with the same accuracy when using different coordinates. This is why the opacity function obtained with exact WP propagation does not match completely the exact TI results, not only for the interpolated J’s but also for those J’s at which “exact” numerical calculations were carried out. The rotationally resolved ICS for HD(V′ ) 0, j′) products, in the middle panel of Figure 4, shows a very flat structure, and the agreement between the two methods is rather satisfactory. The TI results at 0.524 eV exhibit some shallow oscillations which are not present in the prod-WP results at the same energy.

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14495

Figure 5. State-resolved DCS for the H+ + D2(V ) 0, j ) 0) f HD(V′, j′) + D+ reaction at a translational energy of 0.524 eV obtained with the prod-WP method compared with the exact TI results from ref 29. The points are the experimental results of ref 30.

However, prod-WP results at close energies display structures similar to those shown by TI results, which could be interpreted as evidence that such structures are due, again, to some particular resonance. The state-resolved DCS obtained with the TI method presents a very oscillating structure, while the prod-WP results do not, as shown in Figure 5. The random phase approximation invoked in the present WP approach precludes considering all of the interference terms in eq 1, introduced in the TI method. Also, the slight asymmetry of some of the angular distributions obtained with the TI approach cannot be reproduced completely by the WP results, which, according with the above-mentioned approximation, are strictly forward-backward symmetric. Despite these differences, the overall agreement is very good. In fact, both theoretical approaches are found to describe quite well the available experimental results.30 The DCS seems to exhibit a monotonic behavior as a function of energy, as revealed in Figure 6, with two equal maxima at Θ ) 0 and 180°, decreasing as the energy increases. The oscillations obtained as a function of energy are due to resonances as well as to the peaks appearing as a function of Θ in Figure 5, which are lower than the TI results because in the RPA, the interference terms are neglected, but they are noticeable especially for Θ close to 0 or π. Apart from these oscillations, especially important at forward and backward scattering directions, the DCSs show a quite simple structure, in nearly perfect agreement with statistical methods.29 As new j′ channels become open, the cross section shows a sharp

14496

J. Phys. Chem. A, Vol. 113, No. 52, 2009

Zanchet et al.

Figure 6. Three-dimensional plot of the state-resolved DCS (in Å2/sr) as a function of the CM scattering angle and translational energy for the H+ + D2(V ) 0, j ) 0) f HD(V′, j′) + D+.

increase up to a maximum whose value depends linearly on j′ since the number of accessible sublevels increases to 2j′ + 1. For this nearly statistical reaction, the alignment of rotational angular momenta of HD products, shown in Figure 7, becomes nearly always close to 0. This behavior only changes at Θ ) 0 and π where A ) -1 for the reasons commented above. The value of A ) 0 clearly indicates that j′ and k′ are randomly aligned. Since the reaction proceeds through long-lived resonances, the system rotates many times, erasing any dependence on the initial conditions. In addition, the resulting mechanism does not present any particular j′k′ correlation because following statistical arguments, all final helicities are formed with the same probability since they correspond to the same energy for a particular V′, j′ channel. This explains the A ≈ 0 obtained. The small oscillations can be associated with the neglect of interference terms invoked in the RPA and are expected to

increase slightly for an exact treatment, as obtained here for H + D2 reaction, in Figure 2. C. Li + HF Collision in Reactant Coordinates. The two reactions discussed above have the same mass combination, whose skewing angle has an intermediate value, and no particular advantages seem to be gained by using either reactant or product Jacobi coordinates.27 Li + HF f LiF + H, as a prototype of H + H′L f HH′ + L reactions, is more efficiently described in reactant coordinates.27 No quantum simulation of the DCS has been reported until now, and one goal of this work is to calculate it. For this purpose, reac-WP calculations have been done here for J ) 0, 1, 2, ..., 45, including Ωmax ) 7, 15, and 31. The vibrationally resolved reaction probabilities (Figure 8) obtained for the three cases are indistinguishable, thus ensuring a good convergence with the number of Ω’s. The total reaction cross section is shown in the top panel of Figure 8 up

Cross Sections and Polarization in A + BC Reactions

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14497

Figure 7. Alignment of the HD rotational products as a function of the scattering angle for the H+ + D2(V ) 0, j ) 0) f HD(V′, j′) + D+ at a translational energy of 0.524 eV.

to 0.25 eV, and it is compared with previous QCT results obtained using the same PES. The experimental results of Ho¨bel et al.,40 also included in Figure 8, in arbitrary units, have been multiplied by an arbitrary scaling factor. The total cross section exhibits oscillations as a function of energy, which are associated with the corresponding ones of the reaction probabilities, in turn explained by bound-state structures at the transition state.62 In addition to these structures, the reaction at low energy seems to be mediated by a large number of resonance peaks. This feature, also found in other PESs,62 would explain the differences found between QCT and quantum results. A better comparison with the experimental results should include a distribution of initial rotational states of the HF reactants. The QCT and reac-WP total DCSs are compared in Figure 9 for some selected collisional energies and present significant differences, which are also attributed to quantum effects, such as the resonances discussed above. The angular distributions display a complicated structure with an overall preference toward forward scattering. A further analysis at a state-to-state level reveals a different behavior of the different final LiF(V′, j′) rovibrational states, as shown in Figure 10. Thus, for V′ ) 0, the DCS corresponds essentially to forward scattering; for V′ ) 1, it is backward, while for V′ ) 2, it is the sideways direction. These differences indicate that there are several mechanisms; the existence of a prominent peak in either the forward or backward direction for V′ ) 0 and 1 is typical of a direct stripping or rebound mechanism, respectively, while the sideways, more isotropic distribution is more typical of an indirect mechanism.106 The indirect pathways for V′ ) 2 could be explained by the resonances, leading the system to rotate several

Figure 8. Reaction probabilities (bottom panels) for the Li + HF(V ) j ) 0, J) f LiF(V′) + H reaction as a function of translational energy. (Top panel) Li + HF(V ) j ) 0) total reaction cross section obtained in this work, compared with previous QCT41 and experimental40 results.

times during the formation of a complex, which induces a much isotropic angular distribution. The direct processes are interpreted in terms of the so-called direct interaction with product repulsion (DIPR) mechanism.107 The ground electronic state is the result of a curve crossing between a covalent and an ionic electronic state correlating with reactants and products, respectively. The saddle point in the ground electronic state is the result of the crossing. The sudden change from covalent to ionic in this region can be interpreted as a “jump” of an electron from the Li to the H atom, following a harpoon-like process, at a bent geometry. The transient HFis in a repulsive potential curve, and the H atoms fly apart rapidly, “pushing” the F atom, leading to rotationally excited Li+F- products. The departing H atom is very light, and it can be assumed that l ′ , j′. Since l > j, it can be argued that l ≈ J ≈ j′. The direct stripping and rebound mechanisms are associated with high and low total angular momenta, which would imply high and low j′ values, respectively. The two distributions obtained for V′ ) 0 and 1 present a double forward/ backward peak at approximately the same j′ ) 25 value in Figure 10. The dominance of the forward/backward peaks depends on J. Thus, when including only J < 20, the velocity distribution for V′ ) 2 does not appreciably change. However, for V′ ) 0 and 1, significant modifications can be observed; with only low J’s, the velocity distribution becomes completely sideways, with

14498

J. Phys. Chem. A, Vol. 113, No. 52, 2009

Zanchet et al.

Figure 10. Rotationally resolved differential cross section (in Å2/sr) for the Li + HF (V ) j ) 0) f LiF(V′, j′) + H reaction for a translational energy of 0.241 eV, for the three final vibrational states, V′ ) 0, 1, and 2.

Figure 9. Total differential cross section for the Li + HF(V ) j ) 0) (summing on all final LiF states) for several translational energies and compared with the QCT results of ref 41.

much lower probabilities. Note that the system can, in principle, rotate more for high J’s, thus leading to broader angular distributions. Sideways distributions are therefore attributed to the resonances appearing at low J’s, leading to indirect mechanisms. The angular distribution of recoil velocities of products, ′ ), is more directly obtained in time-of-flight or similar D(Θ,Ekin experiments. Such a quantity corresponds to an average over the initial relative velocity distribution (of energy Ekin) between reactants, that is

D(Θ, E′kin) )

∑ ∫ dEkin g(Ekin)δ(EVj + Ekin V'j'

EV'j' - E′kin)

∂σV,jfV',j'(Ekin) ∂Θ

(24)

where g(Ekin) is the distribution of relative translational energy, 0 2 which is taken as a Gaussian function, exp[-(Ekin - Ekin ) /∆2]. The need for DCSs at different energies does not imply an additional computational effort since there is information over a wide range of energies in WP calculations. The results obtained here for three selected initial kinetic energies, E0kin, are 0 ) 0.241 eV, the shown in Figure 11 for ∆ ) 5 meV. For Ekin average does not change the discussion made above for the stateresolved cross section shown in Figure 11. Since the energy of LiF(V′ ) 0, j′ ) 0) fragments is 0.08 eV below the HF(V ) 0,

j ) 0), the maximum ring for V′ ) 0 is at ∼0.32 eV. For V′ ) 0, the backward peak is three times lower than the forward one. For V′ ) 1, 0.03 eV above the initial threshold, the rings are shorter, presenting the preference toward the backward direction. Finally, for V′ ) 2 at 0.14 eV, the distribution is clearly sideways and wider. For the other two lower energies, V′ ) 2 is closed, and D(Θ,E′kin) for V′ ) 0 is also mainly forward but shows important sideways and broad wings. For V′ ) 1, it is backward but, in general, much broader. In fact, for these lower energies, the reaction is mediated by resonances,62 as can be seen in the reaction probabilities shown in Figure 8, which explains why the intensity is more isotropically distributed. The rotational alignment of LiF products, shown in Figure 12, is always very close to -1 in the forward/backward directions. This happens in a much broader interval than in the previous two cases and therefore is not simply due to the J (Θ), as geometric aspect introduced by the behavior of dΩΩ′ discussed above. This clearly indicates that j′ is almost perpendicular to k′ for the backward and forward directions. This can be explained by the direct DIPR mechanism since when H leaves, it pushes the F atom, making LiF rotate in a plane containing k′, so that j′ is perpendicular to it. Since this happens for all energies and all j′ values, and considering that j′ ≈ l, it may be concluded that j′ is perpendicular to the kk′ plane. Similar features were already found in electric depletion experiments,108 in agreement with phase-space model simulations.4 For sideways angles, however, the alignment becomes positive again for nearly all j′ and energies considered. For V′ ) 0, the sideways probability is much lower than that for the forward or backward directions but still significant to allow its numerical determination. The corresponding alignment gets closer to 2 for the three energies considered, indicating that j′

Cross Sections and Polarization in A + BC Reactions

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14499

Figure 11. Velocity map distribution of the differential cross section (in Å2/sr) for the Li + HF(V ) 0, j ) 0) f LiF(V′) + H at translational energies of (a) 0.241, (b) 0.136, and (c) 0.097 eV.

The fact that A becomes nearly 2 should then be attributed to resonances associated with Ω′ . 0. IV. Conclusions

Figure 12. Rotational alignment parameters for the products as a function of scattering angle for the Li + HF(V ) j ) 0) f LiF(V′ ) 0, 1, 2, j′) + H reaction at translational energies of 0.241. 0.136, and 0.097 eV.

becomes nearly parallel to k′. If such change is attributed to the existence of indirect processes or resonances, it is interesting that A is not 0, as it is for the H+ + D2 case discussed above. The process does not seem to become statistical probably because the corresponding resonances are much wider, with shorter lifetimes, than those attributed to the transition state.62

In this work, we have presented a wave packet method to calculate state-to-state differential cross sections for A + BC f AB + C reactive collisions. The present method has several things in common with a couple of previous wave packet methods13,14 in which product Jacobi coordinates were used. Apart from some small technical differences, the major innovation of the present method is that reactant Jacobi coordinates are used as described in ref 27 to calculate the S matrix elements. The method has been implemented with a FORTRAN code, MAD-WAVE3, which uses either reactant or product Jacobi coordinates to extract the state-to-state S matrix elements and is parallelized using the MPI library in the helicity components Ω and in the angular grid points. The differential cross sections obtained for the neutral H + D2 reaction obtained with reactant Jacobi coordinates are in excellent agreement with hyperspherical close coupling calculations performed with the ABC code. This demonstrates that the required coordinate transformation procedure is both efficient and accurate. Surprisingly, it has been found that the differential cross section converges for the same number of helicity components, Ωmax ) 11, when using reactant or product coordinates. However, the calculations using reactant coordinates are well converged with a lower number of angular grid points because of the permutation of the two identical nuclei is directly expressed in this set of coordinates. As a consequence, the use of reactant coordinates is more efficient, by nearly a factor of 2.

14500

J. Phys. Chem. A, Vol. 113, No. 52, 2009

For the ionic H+ + D2 reaction, the situation is different. For low J’s, reactant and product coordinates are of comparable efficiency. However, for high J, and because of the important role of Coriolis couplings at longer distances, the number of angular points required to properly describe the products for high Ω values is very high, the skewing angle is rather small, and the product channel described using reactants coordinates is limited to a small angular interval close to γ ) 0. This problem is considered to be less important for systems involving heavier atoms since the Coriolis coupling reduces. The differential cross section has then been calculated for the H+ + D2 reaction using product Jacobi coordinates. Since the calculation of all of the partial waves is very demanding computationally, the random phase approximation is used, allowing one to calculate only a few selected J values and to interpolate the S2 matrix elements for intermediate J’s. The results obtained are in very good agreement with hyperspherical close coupling calculations performed previously on this system for E ) 0.524 eV.29 The only difference is the diminution of the intensity of the existing oscillations in the TI calculations because of the application of the random phase approach. The differential cross section is then calculated on a broader energy range for several HD(V′ ) 0, j′), showing a rather similar structure, close to that obtained with statistical quantum methods.104,105 Finally, the reac-WP method has been applied to the Li + HF(V ) 0, j ) 0) reactive collisions as a prototype of LH + H′ f L + HH′ reactions, which is very well adapted for the use of reactant Jacobi coordinates.27 Thus, the state-to-state differential cross sections have been simulated for the first time using an “exact” quantum treatment. The results are compared with previous QCT calculations,41 and the rotational alignment of the LiF(V′, j′) products indicates that j′ is perpendicular to the k-k′ plane for the forward/backward peaks, while it seems to be nearly isotropically distributed for sideways distributions. The differential cross sections are well converged by using only up to Ω ) 15, and this shows the ability of the present method to describe such reactions. Acknowledgment. We thank Prof. J. Aoiz for interesting discussions and for providing us with the QCT results obtained for the Li + HF reaction. We want to acknowledge Carlos Mourino for his supervision of the parallelization programming with MPI and CESGA for the computing facilities. This work has been supported by the Ministerio de Ciencia e Innovacio´n under Projects CTQ2007-62898 and FIS2007-62006. S.G.C acknowledges the Alexander von Humboldt Foundation for the finatial support. The calculations were performed on the IFF and CESGA computers. References and Notes (1) Casavecchia, P.; Balucani, N.; Volpi, G. G. Annu. ReV. Phys. Chem. 1999, 50, 347. (2) Liu, K. Annu. ReV. Phys. Chem. 2001, 52, 139. (3) Yang, X. Int. ReV. Phys. Chem. 2005, 24, 37. (4) Case, D. R.; Herschbach, D. R. Mol. Phys. 1975, 30, 1537. (5) de Miranda, M. P.; Clary, D. C. J. Chem. Phys. 1997, 106, 4509. (6) Aoiz, F. J.; Martı´nez, M. T.; Sa´ez-Ra´banos, V. J. Chem. Phys. 2001, 114, 8880. (7) Aoiz, F. J.; Herrero, V.; Sa´ez-Ra´banos, V. J. Chem. Phys. 1992, 97, 7423–7436. (8) Balucani, N.; Casavecchia, P.; Banares, L.; Aoiz, F. J.; Gonza´lezLezana, T.; Honvault, P.; Launay, J. M. J. Phys. Chem. A 2006, 110, 817. (9) Schatz, G. C.; Kuppermann, A. J. Chem. Phys. 1976, 65, 4642. (10) Launay, J. M.; Le Dourneuf, M. Chem. Phys. Lett. 1989, 163, 178. (11) Skouteris, D.; Castillo, J. F.; Manolopoulos, D. E. Comput. Phys. Commun. 2000, 133, 128.

Zanchet et al. (12) Althorpe, S. C. J. Chem. Phys. 2001, 114, 1601. (13) Hankel, M.; Smith, S. C.; Allan, R. J.; Gray, S. K.; Balint-Kurti, G. G. J. Chem. Phys. 2006, 125, 164303. (14) Lin, S. Y.; Guo, H. Phys. ReV. A 2006, 74, 022703. (15) Kuppermann, A.; Hipes, P. G. J. Chem. Phys. 1986, 84, 5962. (16) Pack, R. T; Parker, G. A. J. Chem. Phys. 1987, 87, 3888. (17) Castillo, J. F.; Manolopoulos, D. E.; Stark, K.; Werner, H.-J. J. Chem. Phys. 1996, 104, 6531. (18) Pogrebnya, S. K.; Palma, J.; Clary, D. C.; Echave, J. Phys. Chem. Chem. Phys. 2000, 2, 693. (19) Balakrishnan, N.; Kalyanaraman, C.; Sathyamurthy, N. Phys. Rep. 1997, 280, 79. (20) Balint-Kurti, G. G. AdV. Chem. Phys. 2004, 128, 249. (21) Althorpe, S. C.; Clary, D. C. Annu. ReV. Phys. Chem. 2003, 54, 493. (22) Peng, T.; Zhang, J. Z. H. J. Chem. Phys. 1996, 105, 6072. (23) Althorpe, S. C.; Kouri, D. J.; Hoffman, D. K. J. Chem. Phys. 1997, 107, 7816. (24) Althorpe, S. C.; Ferna´ndez-Alonso, F.; Bean, B. D.; Ayers, J. D.; Pomerantz, A. E.; Zare, R. N.; Wrede, E. Nature 2002, 416, 67–70. (25) Juanes-Marcos, J. C.; Althorpe, S. C.; Wrede, E. Science 2005, 309, 1227. (26) Koszinowski, K.; Goldberg, N. T.; Zhang, J.; Zare, R. N.; Bouakline, F.; Althorpe, S. C. J. Chem. Phys. 2007, 127, 124315. (27) Go´mez-Carrasco, S.; Roncero, O. J. Chem. Phys. 2006, 125, 054102. (28) Balint-Kurti, G. G.; Dixon, R. N.; Marston, C. C. Int. ReV. Phys. Chem. 1992, 11, 317. (29) Carmona-Novillo, E.; Gonza´lez-Lezana, T.; Roncero, O.; Honvault, P.; Launay, J.-M.; Bulut, N; Aoiz, F. J.; Banares, L.; Trottier, A.; Wrede, E. J. Chem. Phys. 2008, 128, 014304. (30) Song, H.; Dai, D.; Wu, G.; Wang, C. C.; Harich, S. A.; Hayes, M.; Wang, X.; Gerlich, D.; Yang, X.; Skodje, R. T. J. Chem. Phys. 2005, 123, 074314. (31) Aguado, A.; Sua´rez, C.; Paniagua, M. Chem. Phys. 1995, 201, 107. (32) Parker, G. A.; Lagana`, A.; Crocchianti, S; Pack, R. T J. Chem. Phys. 1995, 102, 1238. (33) Aguado, A.; Paniagua, M.; Lara, M.; Roncero, O. J. Chem. Phys. 1997, 107, 10085. (34) Burcl, R.; Piecuch, P.; Spirko, V.; Bludsky´, O. Int. J. Quantum Chem. 2000, 80, 916. (35) Jasper, A. W.; Hack, M. D.; Truhlar, D. G.; Piecuch, P. J. Chem. Phys. 2002, 116, 8353. (36) Aguado, A.; Paniagua, M.; Sanz, C.; Roncero, O. J. Chem. Phys. 2003, 119, 10088. (37) Becker, C. H.; Casavecchia, P.; Tiedemann, P. W.; Valentini, J. J.; Lee, Y. T. J. Chem. Phys. 1980, 73, 2833. (38) Loesch, H. J.; Stienkemeier, F. J. Chem. Phys. 1993, 98, 9570. (39) Ho¨bel, O.; Mene´ndez, M.; Loesch, H. J. Phys. Chem. Chem. Phys. 2001, 3, 3633. (40) Hobel, O.; Bobbenkamp, R.; Paladini, A.; Russo, A.; Loesch, H. J. Phys. Chem. Chem. Phys. 2004, 6, 2198. (41) Bobbenkamp, R.; Paladini, A.; Russo, A.; Loesch, H. J.; Mene´ndez, M.; Verdasco, E.; Aoiz, F. J; Werner, H.-J. J. Chem. Phys. 2005, 122, 244304. (42) Hudson, A. J.; Oh, H. B.; Polanyi, J. C.; Piecuch, P. J. Chem. Phys. 2000, 113, 9897. (43) Paniagua, M.; Aguado, A.; Lara, M.; Roncero, O. J. Chem. Phys. 1998, 109, 2971. (44) Topaler, M.; Piecuch, P.; Truhlar, D. G. J. Chem. Phys. 1999, 110, 5634. (45) Chang, X. Y.; Ehlich, R.; Hudson, A. J.; Piecuch, P.; Polanyi, J. C. Faraday Discuss. 1997, 108, 411. (46) Topaler, M. S.; Truhlar, D. G.; Chang, X. Y.; Piecuch, P.; Polanyi, J. C. J. Chem. Phys. 1998, 108, 5378. (47) Aguado, A.; Lara, M.; Paniagua, M.; Roncero, O. J. Chem. Phys. 2001, 114, 3440. (48) Zeiri, Y.; Katz, G.; Kosloff, R.; Topaler, M. S.; Truhlar, D. G.; Polanyi, J. C. Chem. Phys. Lett. 1999, 300, 523. (49) Jasper, A. W.; Hack, M. D.; Chakraborty, A.; Truhlar, D. G.; Piecuch, P. J. Chem. Phys. 2001, 115, 7945. (50) Aguado, A.; Paniagua, M.; Sanz, C.; Roncero, O. In Recent AdVances in the Theory of Chemical and Physical Systems; Progress in theoretical chemistry and Physics Series; Springer: Dordrecht, The Netherlands, 2006; Vol. 15, p 385. (51) Alvarino, J. M.; Casavecchia, P.; Gervasi, O.; Lagana`, A. J. Chem. Phys. 1982, 77, 6341. (52) Aoiz, F. J.; Martı´nez, M. T.; Mene´ndez, M.; Sa´ez-Ra´banos, V.; Verdasco, E. Chem. Phys. Lett. 1999, 299, 25. (53) Aoiz, F. J.; Verdasco, E.; Sa´ez-Rabanos, V.; Loesch, H. J.; Mene´ndez., M.; Stienkemeier, F. Phys. Chem. Chem. Phys. 2000, 2, 541. (54) Miller, D. H.; Wyatt, R. J. Chem. Phys. 1987, 86, 5557.

Cross Sections and Polarization in A + BC Reactions (55) Baer, M.; Garcı´a, E.; Lagana, A.; Gervasi, O. Chem. Phys. Lett. 1989, 158, 362. (56) Baer, M.; Last, I.; Loesch, H. J. J. Chem. Phys. 1994, 101, 9648. (57) Lagana, A.; Bolloni, A.; Crocchianti, S. Phys. Chem. Chem. Phys. 2000, 2, 535. (58) Skouteris, D.; Crocchianti, S.; Lagana, A. Chem. Phys. 2008, 349, 170. (59) Gogtas, F.; BalintKurti, G. G.; Offer, A. J. Chem. Phys. 1996, 104, 7927. (60) Aguado, A.; Paniagua, M.; Lara, M.; Roncero, O. J. Chem. Phys. 1997, 106, 1013. (61) Zhu, W.; Wang, D.; Zhang, J. Z. H. Theor. Chem. Acc. 1997, 96, 31. (62) Lara, M.; Aguado, A.; Roncero, O.; Paniagua, M. J. Chem. Phys. 1998, 109, 9391. (63) Lara, M.; Aguado, A.; Paniagua, M.; Roncero, O. J. Chem. Phys. 2000, 113, 1781. (64) Xie, D. Q.; Li, S.; Guo, H. J. Chem. Phys. 2002, 116, 6391. (65) Weck, P.; Balakrishnan, N. J. Chem. Phys. 2005, 122, 154309. (66) Weck, P.; Balakrishnan, N. J. Chem. Phys. 2005, 122, 234310. (67) Loesch, H. J.; Stenzel, E.; Wu¨stenbecker, B. J. Chem. Phys. 1991, 95, 3841. (68) Loesch, H. J.; Stienkemeier, F. J. Chem. Phys. 1993, 99, 9598. (69) Loesch, H. J. Annu. ReV. Phys. Chem. 1995, 46, 555. (70) Noorbatcha, I.; Sathyamurthy, N. Chem. Phys. Lett. 1982, 93, 432. (71) Alvarino, J. M.; Aquilanti, V.; Cavalli, S.; Crocchianti, S.; Lagana`, A.; Martı´nez, T. J. Chem. Phys. 1997, 107, 3339. (72) Alvarino, J. M.; Aquilanti, V.; Cavalli, S.; Crocchianti, S.; Lagana`, A.; Martı´nez, T. J. Phys. Chem.A 1998, 102, 9638. (73) de Miranda, M. P.; Aoiz, F. J.; Banares, L.; Sa´ez-Ra´banos, V. J. Chem. Phys. 1999, 111, 5368. (74) Skouteris, D.; Crocchianti, S.; Lagana, A. Chem. Phys. Lett. 2007, 440, 1. (75) Zhang, J. Z. H.; Miller, W. H. J. Chem. Phys. 1989, 91, 1528. (76) Child, M. S. Molecular Collision Theory; Dover Publications, Inc.: Mineola, NY, 1996. (77) Zare, R. N. Angular Momentum; John Wiley and Sons, Inc.: New York, 1988. (78) Siebbeles, L. D. A.; Glass-Maujean, M.; Vasyutinskii, O. S.; Beswick, J. A.; Roncero, O. J. Chem. Phys. 1994, 100, 3610. (79) Aquilanti, V.; Cavalli, S.; Grossi, G.; Anderson, R. W. J. Phys. Chem. 1991, 95, 8182. (80) Huang, Y.; Kouri, D. J.; Hoffman, D. K. J. Chem. Phys. 1994, 101, 10493. (81) Mandelshtam, V. A.; Taylor, H. S. J. Chem. Phys. 1995, 103, 2903.

J. Phys. Chem. A, Vol. 113, No. 52, 2009 14501 (82) Huang, Y.; Iyengar, S. S.; Kouri, D. J.; Hoffman, D. K. J. Chem. Phys. 1996, 105, 927. (83) Kroes, G.-J.; Neuhauser, D. J. Chem. Phys. 1996, 105, 8690. (84) Chen, R.; Guo, H. J. Chem. Phys. 1996, 105, 3569. (85) Gray, S. K.; Balint-Kurti, G. G. J. Chem. Phys. 1998, 108, 950. (86) Mandelshtam, V. A.; Taylor, H. S. J. Chem. Phys. 1995, 102, 7390. (87) Tal-Ezer, H.; Kosloff, R. J. Chem. Phys. 1984, 81, 3967. (88) Riss, U. V.; Meyer, H.-D. J. Chem. Phys. 1996, 105, 1409. (89) Balint-Kurti, G. G.; Gonza´lez, A. I.; Goldfield, E. M.; Gray, S. K. Faraday Discuss. 1998, 110, 169. (90) Dubernet, M.-L.; Hutson, J. J. Chem. Phys. 1994, 101, 1939. (91) Roncero, O.; Caloto, D.; Janda, K. C.; Halberstadt, N. J. Chem. Phys. 1997, 107, 1406. (92) Goldfield, E. M.; Gray, S. K. Comput. Phys. Commun. 1996, 84, 1. (93) Balint-Kurti, G. G.; Dixon, R. N.; Marston, C. C. J. Chem. Soc., Faraday Trans. 1990, 86, 1741. (94) Boothroyd, A. I.; Keogh, W. J.; Martin, P. G.; Peterson, M. R. J. Chem. Phys. 1991, 95, 4343. (95) Boothroyd, A. I.; Keogh, W. J.; Martin, P. G.; Peterson, M. R. J. Chem. Phys. 1996, 104, 7139. (96) Aguado, A.; Roncero, O.; Tablero, C.; Sanz, C.; Paniagua, M. J. Chem. Phys. 2000, 112, 1240. (97) Aguado, A.; Paniagua, M.; Werner, H.-J. 2004, Unpublished. (98) Dai, D.; Wang, C. C.; Harich, S. A.; Wang, X.; Yang, X.; Chao, S. D.; Skodje, R. T. Science 2003, 300, 1730. (99) Chu, T.-S.; Han, K.-L.; Hankel, M.; Balint-Kurti, G. G. J. Chem. Phys. 2007, 126, 214303. (100) Honvault, P.; Launay, J.-M. Chem. Phys. Lett. 1998, 287, 270. (101) Gonza´lez-Lezana, T.; Aguado, A.; Paniagua, M.; Roncero, O. J. Chem. Phys. 2005, 123, 194309. (102) Chu, T.-S.; Han, K.-L. J. Phys. Chem. A 2005, 109, 2050. (103) Rackhan, E. J.; Gonza´lez-Lezana, T.; Manolopoulos, D. E. J. Chem. Phys. 2003, 119, 12895. (104) Rackham, E. J.; Huarte-Larranaga, F.; Manolopoulos, D. E. Chem. Phys. Lett. 2001, 343, 356. (105) Gonza´lez-Lezana, T. Int. ReV. Phys. Chem. 2007, 26, 29. (106) Herschbach, D. R. AdV. Chem. Phys. 1966, 10, 319. (107) Mestdagh, J.-M.; Soep, B.; Gaveu, M.-A.; Visticot, J.-P. Int. ReV. Phys. Chem. 2003, 22, 285. (108) Hsu, D. S. Y.; Weinstein, N. D.; Herschbach, D. R. Mol. Phys. 1975, 29, 257.

JP9038946