Harmonic inversion of time cross-correlation functions: The optimal way to perform quantum or semiclassical dynamics calculations

Computing vibrational energy levels by solving linear equations using a tensor method with an imposed rank

Present day computers do not have enough memory to store the high-dimensional tensors required when using a direct product basis to compute vibrational energy levels of a polyatomic molecule with more than about five atoms. One way to deal with this problem is to represent tensors using a tensor format. In this paper, we use the canonical polyadic (CP) format. Energy levels are computed by building a basis from vectors obtained by solving linear equations. The method can be thought of as a CP realization of a block inverse iteration method with multiple shifts. The CP rank of the tensors is fixed, and the linear equations are solved with an method. There is no need for rank reduction and no need for orthogonalization, and tensors with a rank larger than the fixed rank used to solve the linear equations are never generated. The ideas are tested by computing vibrational energy levels of a 64-D bilinearly coupled model Hamiltonian and of acetonitrile (12-D).

On-the-Fly ab Initio Semiclassical Calculation of Glycine Vibrational Spectrum

We present an on-the-fly ab initio semiclassical study of vibrational energy levels of glycine, calculated by Fourier transform of the wavepacket correlation function. It is based on a multiple coherent states approach integrated with monodromy matrix regularization for chaotic dynamics. All four lowest-energy glycine conformers are investigated by means of single-trajectory semiclassical spectra obtained upon classical evolution of on-the-fly trajectories with harmonic zero-point energy. For the most stable conformer I, direct dynamics trajectories are also run for each vibrational mode with energy equal to the first harmonic excitation. An analysis of trajectories evolved up to 50 000 atomic time units demonstrates that, in this time span, conformers II and III can be considered as isolated species, while conformers I and IV show a pretty facile interconversion. Therefore, previous perturbative studies based on the assumption of isolated conformers are often reliable but might be not completely appropriate in the case of conformer IV and conformer I for which interconversion occurs promptly.

CALCULATION OF TRANSITION AMPLITUDES WITH A SINGLE LANCZOS PROPAGATION

We review in this article a recently proposed energy-global method that is capable of calculating the entire transition amplitude matrix with a single Lanczos propagation. This method requires neither explicit computation nor storage of the eigenfunctions, rendering it extremely memory efficient. Procedures are proposed to handle situations where "spurious" eigenvalues aggregate around true eigenvalues due to round-off errors. This method is amenable to both real-symmetric and complex-symmetric Hamiltonians. Applications to molecular spectra and reactive scattering are presented. Its relationships with other methods are also discussed.

FURTHER GENERALIZATION AND NUMERICAL IMPLEMENTATION OF PSEUDO-TIME SCHRÖDINGER EQUATIONS FOR QUANTUM SCATTERING CALCULATIONS

We review and further develop the recently introduced numerical approach [Phys. Rev. Lett. 86, 5031, (2001)] for scattering calculations based on a so called pseudo-time Schrödinger equation, which is in turn a modification of the damped Chebyshev polynomial expansion scheme [J. Chem. Phys. 103, 2903, (1995)]. The method utilizes a special energy-dependent form for the absorbing potential in the time-independent Schrödinger equation, in which the complex energy spectrum is mapped inside the unit disk Ek → uk, where uk are the eigenvalues of some explicitly known sparse matrix U. Most importantly for the numerical implementation, all the physical eigenvalues uk are the extreme eigenvalues of U (i.e. |uk| ≈ 1 for resonances and |uk| = 1 for the bound states), which allows one to extract these eigenvalues very efficiently by harmonic inversion of a pseudo-time autocorrelation function y(t) = ϕ T Ut ϕ using the filter diagonalization method. The computation of y(t) up to time t = 2T requires only T sparse real matrix-vector multiplications. We describe and compare different schemes, effectively corresponding to different choices of the energy-dependent absorbing potential, and test them numerically by calculating resonances of the HCO molecule. Our numerical tests suggest an optimal scheme that provide accurate estimates for most resonance states using a single autocorrelation function.

ON HARMONIC INVERSION OF CROSS-CORRELATION FUNCTIONS BY THE FILTER DIAGONALIZATION METHOD

Harmonic inversion of Chebyshev correlation and cross-correlation functions by the filter diagonalization method (FDM) is one of the most efficient ways to accurately compute the complex spectra of low dimensional quantum molecular systems. This explains the growing popularity of the FDM in the past several years. Some of its most attractive features are the predictable convergence properties and the lack of adjusting parameters. These issues however are often misunderstood and mystified. We discuss the questions relevant to the optimal choices for the FDM parameters, such as the window size and the number of basis functions. We also demonstrate that the cross-correlation approach (using multiple initial states) is significantly more effective than the conventional autocorrelation approach (single initial state) for the common case of a non-uniform eigenvalue distribution.

Non-normal Lanczos methods for quantum scattering

This article presents a new complex absorbing potential (CAP) block Lanczos method for computing scattering eigenfunctions and reaction probabilities. The method reduces the problem of computing energy eigenfunctions to solving two energy dependent systems of equations. An energy independent block Lanczos factorization casts the system into a block tridiagonal form, which can be solved very efficiently for all energies. We show that CAP-Lanczos methods exhibit instability due to the non-normality of CAP Hamiltonians and may break down for some systems. The instability is not due to loss of orthogonality but to non-normality of the Hamiltonian matrix. While use of a Woods-Saxon exponential CAP-as opposed to a polynomial CAP-reduced non-normality, it did not always ensure convergence. Our results indicate that the Arnoldi algorithm is more robust for non-normal systems and less prone to break down. An Arnoldi version of our method is applied to a nonadiabatic tunneling Hamiltonian with excellent results, while the Lanczos algorithm breaks down for this system.

Long time wave packet dynamics from energy eigenfunctions: nonuniform energy resolution via adaptive bisection fast Fourier transformation

This article presents a new approach to long time wave packet propagation. The methodology relies on energy domain calculations and an on-the-surface straightforward energy to time transformation to provide wave packet time evolution. The adaptive bisection fast Fourier transform method employs selective bisection to create a multiresolution energy grid, dense near resonances. To implement fast Fourier transforms on the nonuniform grid, the uniform grid corresponding to the finest resolution is reconstructed using an iterative interpolation process. By proper choice of the energy grid points, we are able to produce results equivalent to grids of the finest resolution, with far fewer grid points. We have seen savings 20-fold in the number of eigenfunction calculations. Since the method requires computation of energy eigenfunctions, it is best suited for situations where many wave packet propagations are of interest at a fixed small set of points--as in time dependent flux computations. The fast Fourier transform (FFT) algorithm used is an adaptation of the Danielson-Lanczos FFT algorithm to sparse input data. A specific advantage of the adaptive bisection FFT is the possibility of long time wave packet propagations showing slow resonant decay. A method is discussed for obtaining resonance parameters by least squares fitting of energy domain data. The key innovation presented is the means of separating out the smooth background from the sharp resonance structure.

Reduced-Order Modeling, Error Estimation, and the Role of the Start-Vector: The Recursive Residue Generation Method Revisited

The recursive residue generation method (RRGM) [Wyatt, R. E. Adv. Chem. Phys. 1989, 73, 231] is re-derived using the formalism of reduced-order modeling [Bai, Z. Appl. Numerical Math. 2002, 43, 9]. A stopping criteria for the RRGM recursions is proposed, on the basis of an expression for an upper bound to the absolute error [Bai, Z.; Ye, Q. Electron. Trans. Numerical Anal. 1998, 7, 1]. It is further pointed out that, in general, the start-vector has a negligible effect on the convergence of the RRGM.

Dissipative dynamics of a system passing through a conical intersection: Ultrafast pump-probe observables

The dynamics of a system incorporating a conical intersection, in the presence of a dissipative environment, is studied with the purpose of identifying observable ultrafast spectroscopic signatures. A model system consisting of two vibronically coupled electronic states with two nuclear degrees of freedom is constructed. Dissipation is treated by two different methods, Lindblad semigroup formalism and the surrogate Hamiltonian approach. Pump-probe experimental expectation values such as transient emission and transient absorption are calculated and compared to the adiabatic and diabatic population transfer. The ultrafast population transfer reflecting the conical intersection is not mirrored in transient absorption measurements such as the recovery of the bleach. Emission from the excited state can be suppressed on the ultrafast time scale, but the existence of a conical intersection is only one of the possible mechanisms that can provide ultrafast damping of emission.

Classical, semiclassical, and quantum investigations of the four-sphere scattering system

A genuinely three-dimensional system, viz. the hyperbolic four-sphere scattering system, is investigated with classical, semiclassical, and quantum mechanical methods at various center-to-center separations of the spheres. The efficiency and scaling properties of the computations are discussed by comparisons to the two-dimensional three-disk system. While in systems with few degrees of freedom modern quantum calculations are, in general, numerically more efficient than semiclassical methods, this situation can be reversed with increasing dimension of the problem. For the four-sphere system with large separations between the spheres, we demonstrate the superiority of semiclassical versus quantum calculations, i.e., semiclassical resonances can easily be obtained even in energy regions which are unattainable with the currently available quantum techniques. The four-sphere system with touching spheres is a challenging problem for both quantum and semiclassical techniques. Here, semiclassical resonances are obtained via harmonic inversion of a cross-correlated periodic orbit signal.

Semiclassical tunneling splittings from short time dynamics: Herman-Kluk-propagation and harmonic inversion

We investigate a recently proposed method [J. Chem. Phys. 108, 9206 (1998)] to obtain tunneling splittings from short time cross-correlation matrices that were propagated according to the semiclassical propagator of Herman and Kluk. The energy levels were extracted by harmonic inversion of the cross-correlation matrix using the filter diagonalization technique. The aim of this study is twofold: First, the short time behavior of the Herman-Kluk-propagator and the meaning of using cross-correlation matrices rather than autocorrelation functions is addressed. Numerical examples are given for one- and two-dimensional model potentials. Second, the performance of the method is investigated for a system with considerable anharmonicity and coupling. Here the proton transfer in 3,7-dichlorotropolone is considered using an ab initio reaction surface Hamiltonian approach. For this example also the extension to more dimensions is critically discussed.

Locating Pollicott-Ruelle resonances in chaotic dynamical systems: a class of numerical schemes

A class of numerical methods to determine Pollicott-Ruelle resonances in chaotic dynamical systems is proposed. This is achieved by relating some existing procedures that make use of the Padé approximants, and interpolating exponentials to both the memory function techniques used in the theory of relaxation and the filter diagonalization method used in the harmonic inversion of time correlation functions. This relationship leads to a theoretical framework in which all these methods become equivalent and which allows for new and improved numerical schemes.

Pseudotime Schrödinger equation with absorbing potential for quantum scattering calculations

The Schrödinger equation (Hpsi) (r) = [E+u(E)W(r)]psi(r) with an energy-dependent complex absorbing potential -u(E)W(r), associated with a scattering system, can be reduced for a special choice of u(E) to a harmonic inversion problem of a discrete pseudotime correlation function y(t) = phi(T)U(t)phi. An efficient formula for Green's function matrix elements is also derived. Since the exact propagation up to time 2t can be done with only approximately t real matrix-vector products, this gives an unprecedently efficient scheme for accurate calculations of quantum spectra for possibly very large systems.

Semiclassical calculation of chemical reaction dynamics via wavepacket correlation functions

Calculation of chemical reaction dynamics is central to theoretical chemistry. The majority of calculations use either classical mechanics, which is computationally inexpensive but misses quantum effects, such as tunneling and interference, or quantum mechanics, which is computationally expensive and often conceptually opaque. An appealing middle ground is the use of semiclassical mechanics. Indeed, since the early 1970s there has been great interest in using semiclassical methods to calculate reaction probabilities. However, despite the elegance of classical S-matrix theory, numerical results on even the simplest reactive systems remained out of reach. Recently, with advances both in correlation function formulations of reactive scattering as well as in semiclassical methods, it has become possible for the first time to calculate reaction probabilities semiclassically. The correlation function methods are contrasted with recent flux-based methods, which, although providing somewhat more compact expressions for the cumulative reactive probability, are less compatible with semiclassical implementation.

Semiclassical spectra and diagonal matrix elements by harmonic inversion of cross-correlated periodic orbit sums

Semiclassical spectra weighted with products of diagonal matrix elements of operators A(alpha), i.e., g(alphaalpha')(E)= summation operator(n)<nA(alpha)n><nA(alpha')n>/(E-E(n)), are obtained by harmonic inversion of a cross-correlation signal constructed of classical periodic orbits. The method provides highly resolved semiclassical spectra even in situations of nearly degenerate states, and opens the way to reducing the required signal lengths to shorter than the Heisenberg time. This implies a significant reduction of the number of orbits required for periodic orbit quantization by harmonic inversion.