Quantum dynamical characterization of unimolecular resonances

Abstract

We give a selective review of quantum mechanical methods for calculating and characterizing resonances in small molecular systems, with an emphasis on recent progress in Chebyshev and Lanczos iterative methods. Two archetypal molecular systems are discussed: isolated resonances in HCO, which exhibit regular mode and state specificity, and overlapping resonances in strongly bound HO₂, which exhibit irregular and chaotic behavior. Future directions in this field are also discussed.

---

1 Introduction

Resonances are temporarily trapped meta-stable states—the analogue of bound vibrational states in the continuum part of the molecular spectrum. Resonances can be formed by bound-free excitations or by collisions between reactants, hence they are manifested in both unimolecular reactions and bimolecular complex-forming reactions. Unless stabilized by collisional relaxation or spontaneous infrared emission, they will decay into products on a finite time scale that is characterized by their lifetime. In essence, resonances are quantum mechanical phenomenon because they occur at discrete energies (resonance positions), but unlike bound states they have a finite width (resonance width).

Though the phenomenon of resonances has long been recognized and qualitatively understood, the quantitative determination of resonances started to appear only during past two decades. Such quantitative studies are presently limited mostly to triatomic systems. Conceptually, there are two basic approaches to determining resonance energies and widths via computation. The first approach formulates the resonance problem in terms of an eigenvalue/eigenvector calculation, and has many affinities to bound-state calculations. In this approach, the resonance positions and widths are determined as the real and imaginary parts of complex eigenvalues associated with the Hamiltonian under dissipative boundary conditions. The dissipative boundary conditions can be imposed in a number of ways, principal among them being the use of a complex absorbing potential (CAP), e.g.ref. 1, a complex scaled Hamiltonian, e.g.ref. 2, or an incrementally damped matrix recursion, e.g.ref. 3. The second approach is a scattering one which relies on the calculation of the scattering S matrix.^4–6 Resonance states are associated with the complex poles of the S matrix and thus all S matrix related quantities such as the lifetime matrix or scattering probabilities will reflect the resonance structures in their energy-dependent profiles. Analysis of such profiles can enable the determination of resonance energies as well as widths.

A third quantity which is characteristic of a resonance is the product state distributions arising from its decay. The energy of a resonance is essentially determined by the potential energy surface (PES) in the inner region, while its width depends on the coupling between the inner region and the exit channel. The product state distributions reflect scattering from the resonance into product states through the transition state region of the PES, and thus contain additional clues about the intra- and inter-molecular dynamics of the system. To specify product state distributions arising from resonance decay, the resonance eigenfunctions of the dissipative Hamiltonian must be calculated and analyzed for their amplitudes in different product channels. For the case of a bimolecular complex-forming reaction, the scattering waves associated with specific incoming reactant channels must be computed and analyzed in the asymptotic region to extract information about product state distributions. These product state distributions may comprise of contributions from one or more (overlapping) resonances which are accessed at the given scattering energy, together with possible direct scattering components. In order to fully characterize the reaction dynamics, it is necessary to consider resonance energies, widths and product state distributions for as many resonances as possible.

In recent years, quantum calculations based on iterative methods have become increasingly common. These methods have better scaling properties than direct methods because they do not require explicit storage of the Hamiltonian matrix, rather only the multiplication of the Hamiltonian onto a vector. When combined with a sparse representation of the Hamiltonian such as a discrete variable representation (DVR),⁷ both memory and CPU time can be reduced dramatically. In this perspective we will focus on an overview of the application of such iterative methods, highlighting in particular two of the most promising approaches: Chebyshev methods, e.g.ref. 8–12, and Lanczos methods, e.g.ref. 13–19. Other quantum methods^2,5,6,20 based on direct matrix diagonalisation techniques have also continued to develop. Due to limited space, we shall only summarise such approaches briefly. We choose two molecular systems which are dynamically quite different for illustration of the type of information that can be obtained from such detailed computations: the HCO molecule, e.g.ref. 4, 21–23, and the HO₂ molecule, e.g.ref. 3, 6, 15–18 and 24. Not all related references have been listed here, see ref. 25 and 26 for more references for HCO. The HCO molecule represents essentially a regular system, whereas HO₂ represents essentially a chaotic system. Most resonances for HCO are isolated and can be assigned normal mode or local mode labels. In contrast, most resonances in HO₂ are overlapping and defy spectroscopic assignment. Although the quantum fluctuations appear for both systems, the underlying mechanisms are quite different.

Advances in this topic have been in large part due to interplay between experiments and theories. In recent years tremendous progress has been made in detecting resonances in the fully state resolved level for both unimolecular dissociation and bimolecular reactions. Excellent reviews of the experimental studies in this field have appeared, e.g.ref. 27–29, and we limit ourselves here to quantum theoretical investigations. Additionally, we refer the curious reader to several excellent related theoretical reviews on this topic, e.g.ref. 26, 30 and 31. This review is arranged as follows: in Section 2 we describe advances in quantum iterative methods. The two case studies are discussed in Section 3, emphasizing the underlying resonance mechanisms that govern the observables. A summary with outlook to possible future developments is provided in Section 4.

2 Theoretical methods

In the context of the iterative methods which we review in this perspective, it is appropriate identify two general approaches which have been found useful for obtaining stable representations of the Hamiltonian [or equivalently the causal Green operator G⁺(E)] under scattering boundary conditions (see, e.g.ref. 1 and 32). The first involves the incorporation into the Hamiltonian of a (negative imaginary) complex absorbing potential (CAP). The resulting augmented Hamiltonian operator can then be stably represented in terms of a real basis (e.g., a discrete variable grid representation), yielding a complex-symmetric matrix. This approach has many affinities to earlier complex scaling methods. The second general approach involves representation of the real Hamiltonian in terms of a complex basis which satisfies the outgoing-wave scattering boundary conditions, which again leads to a complex-symmetric matrix representation. The two classes of iterative techniques which we review below are, respectively, illustrations of these two general approaches to incorporation of the scattering boundary conditions. The most important contribution of the developments in iterative matrix scattering methods is the ability to implement the essential physics within the context of an efficient low-storage algorithm with scaling properties superior to conventional direct matrix methods.

The “direct” approach to computing quantum resonances involves solving the homogeneous time-independent Schrödinger equation

(E − H′)ψ_E = 0. (1)

Here H′ is a complex-symmetric matrix representation of the Hamiltonian satisfying the outgoing-wave boundary conditions and ψ_E is the resonance wavefunction. This is a bound-state-like eigenproblem. Its solution yields a complex spectrum of which the eigenvalues close to the real axis are typically associated with resonances. For such resonances, the real part of the eigenvalue specifies the resonance energy while the complex part specifies the resonance width which is inversely proportional to its lifetime. The filter-diagonalization approach, pioneered originally by Neuhauser and colleagues (see, e.g.ref. 33), has proven to be a particularly useful method for obtaining converged eigenvalues in targeted spectral windows. The principal is quite simple: one first generates a family of filtered states, filtered about a set of reference energies scanning through the targeted spectral window. Provided the number of filtered states is significantly greater than the actual number of eigenvalues in this window, this set of filtered states should serve as an excellent basis for constructing a subspace representation of the Hamiltonian. Solution of a small generalized eigenproblem then yields the converged eigenvalues within the targeted spectral window (plus some others around the edges or outside the window which will be poorly converged and can be discarded). The primary reason for attempting to target specific spectral windows in this way is related to storage: the width of the targeted spectral windows can be adjusted so that the number of filtered basis functions is not too large. Hence, the dense (generalized) eigenproblem to be solved is always easily manageable. There are many other advantages which have also been recognized and exploited, amongst which are (i) the subspace Hamiltonian matrices (i.e., within the nominated spectral windows) can be constructed directly as the Chebychev or Lanczos recursions proceed without the need to compute the filtered states explicitly in the primary representation—this is a large saving in terms of storage requirements, and (ii) when implemented within Chebychev or Lanczos recursions with the above-mentioned low-storage features, these are “all energy” techniques in the sense that a large number of spectral windows can be managed simultaneously within a single recursion. Both the Lanczos and Chebychev iterative methods which we review below utilize this general approach for determining resonance eigenvalues.

If the process of interest is a bimolecular complex-forming reaction, then one can solve the wavepacket-Lippmann–Schwinger equation⁹ and obtain resonance information by characterizing the poles of the S matrix,

Here |η_α〉 is the initial wavepacket corresponding to the α channel, and ψ^α(E) is the scattering wave function at energy E. This is a scattering theory solution, and again iterative matrix methods such as the Lanczos or Chebychev recursions represent powerful low-storage methods for implementing the action of the Green operator on a wavepacket. From eqn. (2), one typically obtains both resonance scattering and direct scattering information, whereas from eqn. (1) only resonance (or bound) state information is obtained.

2.1 Lanczos subspace method

We begin with the complex-symmetric Lanczos methods which involve the CAP-augmented Hamiltonian. The Lanczos diagonalisation method^34,35 has a long history of use for the calculation of eigenvalues in physics, chemistry, engineering and indeed almost every applied field one can imagine. It is a powerful iterative method for computing eigenvalues and eigenvectors of large matrices (real symmetric, Hermitian or complex symmetric). From a seed vector it generates using a three-term matrix-vector recursion a sequence of states which yield a tridiagonal representation of the original matrix. Its attractive scaling properties stem from the fact that typically just three recursion vectors need be kept in core memory at any one time, and the resulting tridiagonal matrix can be easily diagonalized to yield a spectrum which will contain both converged, unconverged and some spurious or ghost eigenvalues. The tridiagonal representation of the original matrix can also be utilized to iteratively solve linear systems, leading to Lanczos-based iterative linear system solvers such as MINRES³⁶ or QMR³⁷ which can be applied to the solution of eqn. (2). Conventionally, if physically-relevant information such as scattering amplitudes or product state distributions need to be obtained from analysis of the eigenvectors (or linear system solutions), then these need to be computed explicitly in the original DVR representation so that the analysis can be carried out. Recently, however, several advances have been made in application of the Lanczos iteration methods which obviate the need to compute explicit solutions in the primary representation, thus greatly enhancing the flexibility and applicability of the Lanczos approach. These developments include Lanczos-based filter diagonalisation, e.g.ref. 15, 17, for the calculation of real or complex eigenvalues; implicit product state distribution analysis of resonance eigenfunctions;^18,38 and quasi-minimum residual (QMR) solution of the subspace time-independent (TI) wavepacket-Lippmann–Schwinger eqn. (2) to obtain reactive scattering probabilities for arbitrary scattering energies.^19,39 We summarize these new theoretical developments below.

Resonance energies and widths: Lanczos subspace filter diagonalisation. Lanczos subspace filter diagonalisation method (see, e.g.ref. 15, 17), was first introduced by Yu and Smith. In this method, three complex vectors need to be stored and the entire calculation is carried out within a single Lanczos subspace. Since the subspace Hamiltonian is tridiagonal, it is relatively easy to generate a set of filtered states for the chosen energy window by solving quasi-minimum residual equations (QMRFD)¹³ or minimum residual equations (MINRESFD)¹⁴ or a simpler homogeneous linear system (LHFD).¹⁷ The differences among different Lanczos subspace FD versions rely on the methods to generate filtered states and here we only give the summary of LHFD as an example:

(i) Choose a normalised, randomly generated initial vector v₁ ≠ 0 and set β₁, v₀ = 0. Then use the 3-term Lanczos algorithm for complex-symmetric matrices

β_k+1v_k+1 = H′v_k − α_kv_k − β_kv_k−1 (3)

to project the non-Hermitian augmented Hamiltonian into a Krylov subspace. The M × M tridiagonal representation of the Hamiltonian, T_M, has diagonal elements α_k = (v_k|Ĥ′|v_k) and subdiagonal elements β_k = (v_k−1|Ĥ′|v_k). Note that a complex-symmetric inner product is used (i.e., bra vectors are not complex conjugated).

(ii) For all j = 1, 2, … , j_max, generate filtered states ϕ(E_j) by solving the homogeneous linear system

Here a backward substitution recursion is employed:

(a) Choose ϕ_M, the M^th element of ϕ(E_j), to be arbitrary (but non-zero; usually set ϕ_M = 1), and calculate

(b) For k = M − 1, M − 2, …, 2, update scalar ϕ_k−1:

β_kϕ_k−1 = E_jϕ_k − α_kϕ_k − β_k+1ϕ_k+1. (6)

(iii) Construct the overlap matrix with elements S_jj′ = (ϕ(E_j)|ϕ(E_j′)) and subspace Hamiltonian matrix with elements W_jj′ = (ϕ(E_j)|T_M|ϕ(E_j′)). Note that W_jj′ can be calculated using a three-term summation:

(iv) Solve the generalised complex-symmetric eigenvalue problem WB = SBε to obtain the complex energies, {ε}.

(v) Span the energy domain by repeating (ii)–(iv) window by window.

To check the convergence of the eigenvalues as well as the quality of the eigenpairs generated by the above iterative methods, one can typically compute the error norm about the eigenenergy E,

where the Lanczos subspace eigenvector ζ(E) can be cheaply regenerated for each complex eigenenergy using eqn. (4). Alternatively, one can estimate the error norm within the subspace by using a simpler analytic expression. The readers are referred to ref. 16 for more details. Clearly, true eigenvalues should have small error norms and can thus be distinguished from any unconverged/spurious eigenvalues.

Product state distributions from resonance decay. If the resonance eigenfunction is computed out into the asymptotic region then one can determine product state distributions by analyzing its amplitude in the different accessible product channels. Thus, the overlap integrals between the asymptotic channel eigenfunctions and the resonance states must be calculated. Through an elegant construction of the algorithm, both the complex eigenvalues and the product state distributions associated with each individual resonant state can be calculated from a single Lanczos iteration.

One starts with the observation that the wave functions can be transformed between the primary representation and the tri-diagonal Lanczos representation through:

ζ(E_R) = V^Tψ_ER (9)

where V = [v₁, v₂, … , v_M] is the column-orthonormal Lanczos vector matrix. Thus, the overlapping integrals can be re-expressed asThe vectors χ⁽ⁿ⁾ = V^Tφ_n are the subspace projections of the internal eigenfunctions, and can be accumulated in the Lanczos iteration. Therefore, from a single Lanczos recursion, one can calculate the energies, widths and the product state distributions for all converged resonances. This method also has the advantage of low storage requirement because the size of χ⁽ⁿ⁾ vectors is relatively small compared with their counterparts in the primary representation.

Lanczos subspace scattering calculations. Recently, a Lanczos subspace time-independent wavepacket method for calculating scattering probabilities has been developed by Smith et al.^19,39,40 The basic idea is to transform TI wavepacket-Lippmann–Schwinger equation from the primary representationinto a tridiagonal Lanczos representation:

(E − T_M)|ϕ(E)〉 = e₁ . (13)

In eqn. (12), which is equivalent to eqn. (2) above, −iÛ is the absorbing potential. In eqn. (13), e₁ is the first column of the M × M identity matrix, and ϕ(E) is the subspace scattering wave function. By performing the QR factorization,³⁷ the linear system can be solved by the quasi-minimal residual algorithm³⁷where Q_M+1 is a unitary matrix (⁺ indicates the Hermitian adjoint). Only the elements of the upper-triangular matrix R_M with bandwidth 3 and the vector t̃_M+1 = Q_M+1e₁ need to be stored. The subspace scattering wavefunction is obtained via:

Since the Lanczos subspace is independent of a constant energy shift, E, we are able to solve the linear system in eqn. (13) and obtain the desired subspace scattering wave functions for arbitrary E from a single Lanczos recursion.⁴¹ Once again, in a manner similar to that outlined for product state distribution analysis in the previous section, all the overlap integrals necessary for subsequent scattering amplitude analysis can be easily accumulated as the Lanczos recursion progresses. This implies that the scattering wave solutions to eqn. (2) need not be computed explicitly but rather just their subspace analogues from eqn. (13).

2.2 Real Chebyshev recursion methods

The Chebyshev iteration method was first introduced into the reaction dynamics calculations by Tal-Ezer and Kosloff⁸ to expand the evolution operator exp(−iĤt/ħ). Subsequently, Kouri and co-workers⁹ derived a new time-independent (TI) wavepacket-Lippmann–Schwinger equation and presented Chebyshev expansion expressions for both the Green operator and the Dirac delta function. Mandelshtam and Taylor³ introduced a real damping scheme into the Chebyshev recursion, which made the real algorithms possible for dissipative systems. The real Chebychev propagation method can be viewed in an alternative way as a modification of the time-dependent Schrödinger equation. In this respect, two related (discrete/continuous time) forms of the modified equations have been proposed by Chen and Guo⁴² and more generally by Gray and Balint-Kurti.¹² Real Chebyshev recursion methods have grown rapidly in recent years, and have been successfully applied to different fields such as bound or resonance states calculations, e.g.ref. 10, reactive scattering, e.g.ref. 43, and surface scattering, e.g.ref. 44.

The basic Chebyshev iteration is also a three-term recursion:

Where Ĥ_norm = (Ĥ − H̄)/ΔH with H̄ = 0.5(H_max + H_min), and ΔH = 0.5(H_max − H_min). Φ(0) is an initial wavepacket, and e^−r̂ is a damping operator. As will become apparent below, relevant information such as resonance energies, widths, product state distributions or other scattering probabilities can be extracted from either autocorrelation functions c_n = (ξ^ŷ₀,ξ^ŷ_n) or correlation functions , which can be accumulated during the Chebyshev recursion.

Energies and widths: low storage filter diagonalisation. An auto-correlation function based filter diagonalisation scheme was first proposed by Neuhauser et al.³³ In this scheme eigenvalues in the energy windows can be obtained from a single sequence of auto-correlation functions without the need to explicitly construct the filtered states. As indicated above, this eliminates a major memory bottleneck. Mandelshtam et al. combined this scheme with the modified real Chebyshev propagation approach, and implemented a box-like low storage filter diagonalisation (LSFD) formalism.¹⁰ In this method, a partial Fourier transform of the Chebyshev vector is indeed the filtered state. The evolution operator is employed to set up a small-size generalised eigenproblem, and all the matrix elements can be expressed in terms of autocorrelation functions. The most important advantage associated with this approach is that one can employ a real algorithm with a single, extended Chebyshev vector recursion. This leads to substantially greater efficiency in comparison with step-by-step propagation of a complex wavepacket. Other different LSFD versions based on the Chebyshev recursion exist, such as Neuhauser et al.'s versions,⁴⁵ and Chen and Guo's versions.¹¹ For comparisons of these different LSFD methods with Lanczos subspace FD methods, the reader is referred to ref. 46 and 47.

Product state distributions. Resonance wave function can be calculated through a partial time-energy Fourier transformation of the Chebyshev vectors:Here cosϕ = (E − H̄)/ΔH. Time-dependent and time-independent derivations about ψ(E) do exist, but the final expressions of ψ(E) are the same⁴⁸ (see the following part for more details). Like the time-independent Lanczos subspace method, we can arrive at a similar expression for product state distributions in terms of overlapping integrals after obtaining the resonance wave functions.¹⁸ Thus, analogous to the Lanczos algorithm described above, a single sequence of Chebyshev iterations allows us to calculate both energies and product state distributions for different resonances. The difference between the Chebyshev and Lanczos methods is that the coefficient in ψ(E) is analytical in the former case, thus there is no need to solve a linear system as we do for the Lanczos method.

Chebyshev scattering calculations. Eqn. (17) is also the central equation for scattering calculations, which can be understood from both time-dependent and time-independent points of view.^9,18 In the standard time-dependent framework, one can expand the time evolution operator (for the standard time-dependent Schrödinger equation) in terms of Chebyshev polynomials, . Here J_k is the Bessel function. The energy-domain wave function can be obtained by a partial Fourier-transformation of the TD wave function Φ(t), which is indeed eqn. (17). In a novel time-dependent framework,¹² a modified time-dependent Schrödinger equation is utilised to propagate the real part of the wavepacketBy choosing , the real part of the wavepacket Re(Φ′(t)) = Re(Φ′(kτ)) = ξ^ŷ_k can be propagated as a damped Chebyshev iteration. Thus, by performing Fourier transformation of the real part of Φ′(t) and changing the variable from the angle domain to the energy domain, one can again obtain the same energy-dependent wave function as in eqn. (17). In the time-independent framework, one can use Chebyshev polynomials to expand the Green operator , where the expansion coefficient is . Since energy domain wave function is related to the Green operator via, once again we are leading to eqn. (17).

Over the past few years, the applications of real Chebyshev methods into reactive scattering have grown rapidly. For example, Neuhauser et al. combine their FD scheme with damped real Chebyshev algorithm, and extend the real Chebyshev method into molecule-surface scattering.⁴⁴ Gray et al. derived the S matrix formulation from real wavepacket propagation¹² and applied the approach to more challenging systems such as H₂ + O(¹D), (see e.g.ref. 49). Smith et al.^18,48 have recently compared different real Chebyshev versions and the Lanczos method for resonance calculations and scattering calculations. One of their conclusions is that the Chebyshev method can converge the broad outline of the scattering probabilities quickly, whereas the resolution of the individual resonances takes much longer. This indicates that real Chebyshev method is suitable for the large-scale nonzero total angular momentum calculations in reactive scattering for complex forming reactions, in which the fine structures of resonances are not so important to obtain experimental observables. Such J > 0 calculations have been reported recently for half state-resolved probabilities using exact real wavepacket methods, e.g.ref. 43, which require the repeated calculations for many J values. One can predict that exact J > 0 calculations including those at a fully state-resolved level using the real Chebyshev method will grow rapidly in the near future, though they are still very challenging with today's computational resources.

2.3 Related quantum methods

There are several other quantum mechanical methods which have been utilized to deal with resonances. While space does not permit us to discuss them in any detail, we mention in particular the modified version of the log-derivative Kohn's variational principle,⁶ and the truncation/recoupling method,²⁰ which have been very effectively used for resonance calculations in recent years.

3 Case studies

Case 1 HCO—a regular system

HCO is a relatively weakly bound molecule (dissociation energy 0.834 eV and dissociation threshold approximately 0.133 eV), with a low density of resonance states and weak intramolecular coupling. In this case the discrete level structure is important in the dissociation dynamics. There have been many experimental and theoretical studies on HCO, e.g.ref. 4, 21–23, and several review articles related to this molecule have appeared.^26,29 The potential energy surfaces most studies used for this system are RLBH PES⁵⁰ and WKS PES.²¹ The PES has a shallow well and a small barrier separates the inner region from the exit channel. Various quantum methods have been used for this molecule and the results may be summarized as follows:

(i) Most of the resonances for HCO are isolated, and can be unambiguously assigned to three normal mode quantum numbers (ν₁, ν₂, ν₃). Here v₁, v₂, and v₃ represent the H–C stretch mode, the C–O stretch mode, and the bending mode, respectively. Only in some cases does strong coupling between two energetically close-lying zero-order states disturb this clear picture. The assignments can be done by analyzing the nodal structures of the resonance wave functions and several examples are given in Fig. 1 for illustration.⁵¹ The non-vanishing tails of the wave functions toward large values of the dissociation coordinate R are only visible in (e) and (f). The assignment is clear in (a)–(c) and (e). The wave function in (d) seems to follow a particular periodic classical orbit and the state in (f) seems to be unassignable.

Fig. 1 Typical resonance wave functions for HCO. The boxes show the assignments, the resonance positions in eV, and the linewidths in cm⁻¹. Reproduced with permission of the American Institute of Physics from ref. 51.

(ii) The resonance widths (rates) show particularly the striking dependence on the normal mode quantum numbers (ν₁, ν₂, ν₃), i.e., the widths follow the pattern Γ_CH > Γ_bend > Γ_CO. Fig. 2 shows resonance widths for three progressions as a function of the CO stretching quantum number ν₂.²² Resonances with pure excitation in the CO stretching mode (0, ν₂, 0) have the smallest widths (the longest lifetime). This is due to the inefficient energy transfer from r to R, such that the system needs a long time before enough energy is accumulated in the dissociation coordinate to permit dissociation. The bending mode is more strongly coupled to R and addition of one or two quanta of the bending motion can significantly increase the widths. Excitation in R allows a rather rapid bond rupture, and therefore the resonances (ν₁, 0, 0) have the largest widths (the shortest lifetimes). There is quite satisfactory agreement between the quantum results²² and the experimental ones.²³

Fig. 2 Resonance widths as a function of CO stretching quantum number for three progressions. Open squares are calculated results,²² while the full dots are measured ones.²³ Reproduced with permission of the American Institute of Physics from ref. 22.

(iii) The nascent product state distributions of CO also exhibit remarkable mode and state specificity. The CO rotational state distributions are non statistical and highly structured and depend sensitively on the modes excited and on the number of quanta in each mode. Each of the progressions yields distinct distributions. Fig. 3 shows several examples of the final state distributions of CO following the decay of some (0, ν₂, 0) resonances.²¹ These distributions all look rather similar, having a bimodal shape, which reflects the similar shape of the corresponding (0, ν₂, 0) resonance wave functions. In addition, the mode specificity is also shown in the nascent CO vibrational state distributions. For example, the vibrational state distributions from (0, ν₂, 0) resonances are consistently inverted.²¹

Fig. 3 Final state distributions of CO following the decay of (0, ν₂, 0) resonances with ν₂ = 4–8. The vibrational state of CO is the state with the largest probability. Reproduced with permission of the American Institute of Physics from ref. 21.

In summary, the unimolecular reaction of HCO is a clear, illustrative example of mode- and state-specific decomposition. HCO exhibits the isolated resonances that can be assigned normal or local mode labels. Thus the dissociation of HCO is essentially regular. For this system, experimental and theoretical results are in general agreement. This is one example of the great success of quantum mechanical calculations carried out on accurate PES, even at a state resolved level.

Case 2 HO₂—an irregular system

Reactions in which a strongly-bound intermediate exists often support a large number of resonance states, and in this case the reaction dynamics can be completely controlled by resonance formation. Quantum dynamical modeling of such reactions is very challenging since large basis sets are needed to be able to represent the degree of vibrational excitation involved. One option available is to make recourse to theories such as Ramsperger–Rice–Kassel–Marcus (RRKM) theory, e.g.ref. 52, or phase space theory (PST), e.g.ref. 53, which invoke certain statistical approximations (rapid ergodic mixing, no recrossing of the transition state, etc.). Deep-well triatomic molecules such as HO₂ provide a good opportunity to compare the statistical theories with fully quantum simulations. Here we choose HO₂ as an example and present some computational results pertaining to its dissociation rates and product state distributions.

The HO₂ PES (e.g., DMBE IV⁵⁴) has a considerably deeper potential well (D_e = 2.37 eV) than HCO and exhibits no pronounced potential barrier (i.e., no saddle point) during bond rupture. For this system, a meaningful assignment of the resonance states is impossible, and in fact the bound states close to the threshold are already irregular.⁵⁵Fig. 4 depicts the resonance energies and rates calculated from both Lanczos and Chebyshev methods.¹⁸ Note that the rates are widths divided by ħ, and therefore are true unimolecular decay rates only for narrow isolated resonances. From Fig. 4 we see that resonance energies acquired via the two methods are in fairly sound agreement (i.e., 3 or 4 digits of relative accuracy). The resonance rates predicted by both methods agree within 20% of one another. This is because the absorbing boundary conditions are different for the two different methods, and the resonance widths are highly sensitive on the details of absorbing potential or damping function. Strictly speaking, a stabilisation procedure⁵⁶ must be employed to determine more accurate resonance widths. Of course, such stabilisation calculations are still very challenging for large molecules. Analysis of the resonance widths shows that only at very low energies are these resonances are isolated. With increasing energy the resonances begin to overlap. The quantum rates (widths) show a large fluctuation, with the general tendency that the rates increase with energy. Just above threshold, the rates vary over three orders of magnitude, and then the extent of fluctuation decreases with increasing energy. This fluctuating behaviour has been predicted by other quantum calculations on the HO₂ system, e.g.ref. 3, 6.

Fig. 4 Plot of the logarithmic resonance rates, log₁₀(k), versus resonance energy in the energy range 0.09 to 0.47 eV. Circles represent the results from Lanczos subspace filter diagonalisation method and squares denote the results from Chebyshev low storage filter diagonalisation method. The unit of decay rate is s⁻¹. Reproduced with permission of the American Institute of Physics from ref. 18.

Fig. 5(a–d) shows the O₂ rotational distributions for the ground vibrational state at different resonance energies.¹⁸ Firstly, the number of occupied rotational channels increases steadily with energy, which is simply the result of energy conservation. Secondly, the distributions show a very complicated oscillatory behaviour, with the number of oscillations generally increasing with energy. Thirdly, the fluctuations in the distributions seem to be random and unpredictable from resonance to resonance. The rotational state distributions of the fragments reflect the angular dependence of wavefunction at the translational state and the anisotropy of the PES in the exit channel. The rotational state distributions for the HO₂ dissociation are complicated mainly due to the translational–rotational coupling being weak but not negligible. The distributions cannot, therefore, be explained by a simple model such as the Franck–Condon mapping picture or the rotational reflection principle.⁵⁷ The vibrational distributions also fluctuate from resonance to resonance,¹⁸ as can be expected from the rotational state distributions. The general trend is that the ratio increases with energy, approaching unity for the highest resonance energies. For all but several resonances, no population inversion was achieved.

Fig. 5 The O₂ (ν = 0) rotational distributions, at E = (a) 0.154; (b) 0.252; (c) 0.264; and (d) 0.362 eV. Symbols are the same as Fig. 4. All distributions are normalised. Reproduced with permission of the American Institute of Physics from ref. 18.

Whether or not the statistical theories correctly predict the average rates is a question at the heart of unimolecular dissociation. For this system, the decay rates, on average, can be well described by the statistical RRKM theories⁶ (recently, we compared quantum rates with PST for H₂S dissociation and found reasonable agreement³⁸). However, the fluctuations in the rates and product state distributions cannot be produced by the statistical theories. Those fluctuations are purely quantum effects due to the interference of the overlapping resonances. Experimentally, the fluctuations have been observed for NO₂ system in individual NO product distributions as well as in the product vibrational distributions.⁵⁸ Unfortunately, so far there are still no experimental data available for HO₂ system. Thus direct comparisons with experiment for the HO₂ dissociation are still not possible at this time.

Finally, we emphasize the differences of the fluctuation mechanisms in HCO and HO₂ systems. For HO₂, most resonances are overlapping ones and when resonances overlap, interferences will dominate. Thus, the wild fluctuations in the rates and in the product state distributions are a manifestation of prominent quantum interference effects between overlapping resonances. This is manifestly different from the mechanism of HCO fluctuations, for which mode-specificity is responsible. Fluctuations in HCO basically reflect systematically the different types of wave functions of resonance states with different normal or local mode labels. Due to the strong mode specificity, RRKM or PST (which assume rapid ergodic mixing after initial excitation) is seemingly inappropriate for a molecule such as HCO.⁵⁹

Recently, the real Chebyshev recursion and the complex-symmetric Lanczos recursion have been applied to very challenging scattering calculations for complex-forming reactions. As we are primarily concerned with the unimolecular resonances in this review, we give just one example of state-to-state reaction probabilities for the complex forming reaction H + O₂ → HO + O, which also show the signature of the complex resonance structure of the intermediate HO₂. Fig. 6 depicts the state-to-state reactive probabilities associated with the ground O₂ reactant state (ν₀ = 0, j₀ = 1) and ground OH product state (ν = 0, j = 0).¹⁹ Here the resonance energies are well above the unimolecular resonances discussed above. It is apparent that the reaction probabilities are dominated by resonances, most of them being overlapping ones. Analysis of some resonance widths indicates the lifetime of the resonances is of the order of picoseconds. These resonances are clear manifestations of the complex forming and decaying process, and differ from unimolecular resonances only by the way the energized complex is prepared. For more details, the reader is referred to ref. 19.

Fig. 6 The calculated state-to-state reactive probabilities for H + O₂ (ν₀ = 0, j₀ = 1) → O + OH (ν = 0, j = 0) reaction from Lanczos subspace TI wavepacket method. Reproduced with permission of the American Institute of Physics from ref. 19.

4 Concluding remarks

In recent years, resonances in unimolecular and bimolecular reactions have received much attention and tremendous progress has been made in our understanding of resonance-mediated processes at a quantum state-to-state level. Quantum iterative methods such as the Lanczos or Chebyshev recursions have become powerful tools in calculating and characterizing resonances at a detailed level. Such progress is due in part to the increasingly close interplay between experimental and theoretical investigations. So far the most well studied systems are molecules having low threshold energies and weak internal coupling. Resonances in such systems show striking mode and state specificity, which can only be understood by the detailed knowledge of the full PES and the dynamical behavior which governs the unimolecular dissociation in the exit channel. Quantum dynamical methods have proven a triumph in this respect, and statistical theories based on assumptions of rapid ergodic coupling between modes fail.

For complex-forming reactions with a deep well, the overlapping resonances interfere and the resulting fluctuations dominate in rates and in the product state distributions. Statistical theories can predict average rates reasonably well, but have had limited success in capturing the detail of product state distributions. The agreement between theory and experiment has not reached the same level as for weakly bound systems: exact quantum methods are very time-consuming (especially when stabilizing the resonance widths) and the fine details of resonance structure prove quite sensitive to the parameters of the computation (basis set, iterative recursion, method of imposing absorbing boundary conditions, etc.). In light of this, one can make the observation that it will be very difficult to reproduce experimental data at the state-to-state level for such reactions. Going further, one could provocatively query the value of pursuing the theory at all, given that the numerics are clearly so very difficult. In reality, we are not so pessimistic! The fact that one can perform the quantum calculations and observe internally consistent properties and trends (even though the absolute accuracy may still be uncertain) is already an important advance. Even if theory and experiment are not gelling together quantitatively, there are many important issues which can be explored with both approaches.

For the future, there are a number of areas that require further investigation. The role of total angular momenta J > 0 is a major target, especially for bimolecular reactions. The role of IVR in the exit channel in determining the disposal of excess energy into the products is another active area of investigation. For the HCO system several exact or approximate J > 0 calculations^60,61 have been performed and are used to rationalize discrepancies or agreement between experimental and theoretical results. For HO₂, QCT calculations⁶² have inferred that rotational effects will play a significant role. However, while some initial-state-resolved J > 0 scattering calculations have been performed (e.g.ref. 43), no exact J > 0 resonance calculations have been available until very recently.⁶³ In addition, iterative matrix methods have yet to be extended to deal with differential cross sections, a necessary step if the full richness of data from crossed beam experiments is to be modeled. Clearly, it will be important to develop more efficient computational algorithms if questions such as these are to be addressed for such challenging molecular systems.

Using unimolecular resonances to control chemical reactions is another intriguing direction for future investigation. By using intense coherent femtosecond pulses it is possible to excite resonances coherently. For most molecules the laser wave-length needed to reach the unimolecular dissociation limit is in the near-visible range, so it is more practical and economical to manipulate these resonances than to control the dynamics electronically excited states via UV radiation. The use of van der Waals resonances as precursors to a chemical reaction has been proposed several years ago,⁶⁴ and we believe there is much to be done in this area.

Acknowledgements

We are grateful to the Australian Research Council for supporting this work (Discovery Project No. DP0211019). We thank Dr. Anthony Rasmussen and Dr. Hua-Gen Yu for helpful discussions.

References

1. D. Neuhauser and M. Baer, J. Chem. Phys., 1989, 90, 4351.
2. N. Moiseyev, Phys. Rep., 1998, 302, 211.
3. V. A. Mandelshtam, T. P. Grozdanov and H. S. Taylor, J. Chem. Phys., 1995, 103, 10074.
4. K. T. Lee and J. M. Bowman, J. Chem. Phys., 1986, 85, 6225.
5. H. W. Jang and J. C. Light, J. Chem. Phys., 1995, 102, 3262.
6. A. J. Dobbyn, M. Stumpf, H.-M. Keller and R. Schinke, J. Chem. Phys., 1996, 104, 8357.
7. J. C. Light, I. P. Hamilton and J. V. Lill, J. Chem. Phys., 1985, 82, 1400.
8. H. Tal-Ezer and R. Kosloff, J. Chem. Phys., 1984, 81, 3967.
9. D. J. Kouri, M. Arnold and D. K. Hoffman, Chem. Phys. Lett., 1993, 203, 166.
10. V. A. Mandelshtam and H. S. Taylor, J. Chem. Phys., 1997, 107, 6756.
11. R. Chen and H. Guo, J. Chem. Phys., 1999, 111, 464.
12. S. K. Gray and G. G. Balint-Kurti, J. Chem. Phys., 1998, 108, 950.
13. H. G. Yu and S. C. Smith, Chem. Phys. Lett., 1998, 283, 69.
14. H. G. Yu and S. C. Smith, Ber. Bunsen-Ges. Phys. Chem., 1997, 101, 400.
15. H. G. Yu and S. C. Smith, J. Chem. Phys., 1997, 107, 9985.
16. S. C. Smith, H. G. Yu, C. Leforestier and J. C. Rayez, Phys. Chem. Chem. Phys., 1999, 1, 1311.
17. H. Zhang and S. C. Smith, Phys. Chem. Chem. Phys., 2001, 3, 2282.
18. H. Zhang and S. C. Smith, J. Chem. Phys., 2001, 115, 5751.
19. H. Zhang and S. C. Smith, J. Chem. Phys., 2002, 116, 2354.
20. S. Skokov, J. X. Qi, J. M. Bowman, C. Y. Yang, S. K. Gray, K. A. Peterson and V. A. Mandelshtam, J. Chem. Phys., 1998, 109, 10273.
21. H. J. Werner, C. Bauer, P. Rosmus, H. M. Keller, M. Stumpf and R. Schinke, J. Chem. Phys., 1995, 102, 3593.
22. H.-M. Keller, H. Floethmann, A. J. Dobbyn, R. Schinke, H.-J. Werner, C. Bauer and P. Rosmus, J. Chem. Phys., 1996, 105, 4983.
23. J. D. Tobiason, J. R. Dunlop and E. A. Rohlfing, J. Chem. Phys., 1995, 103, 1448.
24. B. Kendrick and R. T. Pack, Chem. Phys. Lett., 1995, 235, 291.
25. J. M. Bowman, J. Phys. Chem. A, 1998, 102, 3006.
26. R. Schinke, C. Beck, S. Y. Grebenshchikov and H. M. Keller, Ber. Bunsen-Ges. Phys. Chem., 1998, 102, 593.
27. G. C. Schatz, Science, 2000, 288, 1599–1600.
28. K. Liu, Annu. Rev. Phys. Chem., 2001, 52, 139.
29. S. A. Reid and H. Reisler, Annu. Rev. Phys. Chem., 1996, 47, 495.
30. R. G. Gilbert and S. C. Smith, Theory of Unimolecular and Recombination Reactions, Blackwell Scientific Publishers, Oxford, 1990.
31. T. Baer and W. L. Hase, Unimolecular Reaction Dynamics: Theory and Experiments, Oxford University Press, New York, 1996.
32. W. H. Miller and B. M. D. D. Jansen op de Haar, J. Chem. Phys., 1987, 86, 6213.
33. M. R. Wall and D. Neuhauser, J. Chem. Phys., 1995, 102, 8011.
34. C. Lanczos, J. Res. Natl. Bur. Stand., 1950, 45, 255.
35. J. K. Cullum and R. A. Willoughby, Lanczos Algorithms for Large Symmetric Eigenvalue Computations, Birkhauser, Boston, 1985.
36. C. C. Paige and M. A. Saunders, SIAM J. Numer. Anal., 1975, 12, 617.
37. R. W. Freund, SIAM J. Sci. Stat. Comput., 1992, 13, 425.
38. H. Zhang and S. C. Smith, J. Phys. Chem. A, 2002, 106, 6129.
39. H. Zhang and S. C. Smith, J. Phys. Chem. A, 2002, 106, 6137.
40. A. J. Rasmussen, S. J. Jeffrey and S. C. Smith, Chem. Phys. Lett., 2001, 336, 149.
41. R. W. Freund, in Numerical Linear Algebra, ed. L. Reichel, A. Ruttan and R. S. Varga, De Gruyter, Berlin, 1993, pp. 101.
42. R. Chen and H. Guo, J. Chem. Phys., 1996, 105, 3569.
43. E. M. Goldfield and A. J. H. M. Meijer, J. Chem. Phys., 2000, 113, 11055.
44. D. A. McCormack, G.-J. Kroes and D. Neuhauser, J. Chem. Phys., 1998, 109, 5177.
45. J. W. Pang, T. Dieckmann, J. Feigon and D. Neuhauser, J. Chem. Phys., 1998, 108, 8360–8368.
46. H. Zhang and S. C. Smith, Chem. Phys. Lett., 2001, 347, 211.
47. D. Xie, R. Chen and H. Guo, J. Chem. Phys., 2000, 112, 5263.
48. H. Zhang and S. C. Smith, J. Chem. Phys., 2002, 117, 5174.
49. M. Hankel, G. G. Balint-Kurti and S. K. Gray, J. Phys. Chem. A, 2001, 105, 2330.
50. H. Romanowski, K.-T. Lee, J. M. Bowman and L. B. Harding, J. Chem. Phys., 1986, 84, 4888.
51. H.-M. Keller and R. Schinke, J. Chem. Phys., 1999, 110, 9887.
52. R. A. Marcus, J. Chem. Phys., 1952, 20, 359.
53. P. Pechukas and J. C. Light, J. Chem. Phys., 1965, 42, 3285.
54. M. R. Pastrana, L. A. M. Quintales, J. Brandao and A. J. C. Varandas, J. Phys. Chem., 1990, 94, 8073.
55. A. J. Dobbyn, M. Stumpf, H.-M. Keller and R. Schinke, J. Chem. Phys., 1995, 103, 9947.
56. A. U. Hazi and H. S. Taylor, Phys. Rev. A, 1970, 1, 1109.
57. R. Schinke, Photodissociation Dynamics, Cambridge University Press, Cambridge, UK, 1993.
58. S. A. Reid and H. Reisler, J. Phys. Chem., 1996, 100, 474.
59. M. Stumpf, A. J. Dobbyn, D. H. Mordaunt, H.-M. Keller, H. Fluethmann and R. Schinke, Faraday Discuss., 1995, 102, 193.
60. C.-Y. Yang and S. K. Gray, J. Chem. Phys., 1997, 107, 7773.
61. J. Qi and J. M. Bowman, J. Chem. Phys., 1996, 105, 9884.
62. J. M. C. Marques and A. J. C. Varandas, J. Phys. Chem. A, 1997, 101, 5168.
63. H. Zhang and S. C. Smith, J. Chem. Phys. 2003, submitted.
64. G. Radhakrishnan, S. Buelow and C. Wittig, J. Chem. Phys., 1986, 84, 727–738.

---