One of the most accurate methods for solving the time-dependent Schrödinger equation uses a combination of the dynamic Fourier method with the split-operator algorithm on a tensor-product grid. To reduce the number of required grid points, we let the grid move together with the wavepacket but find that the naïve algorithm based on an alternate evolution of the wavefunction and grid destroys the time reversibility of the exact evolution. Yet, we show that the time reversibility is recovered if the wavefunction and grid are evolved simultaneously during each kinetic or potential step; this is achieved by using the Ehrenfest theorem together with the splitting method. The proposed algorithm is conditionally stable, symmetric, and time-reversible and conserves the norm of the wavefunction. The preservation of these geometric properties is shown analytically and demonstrated numerically on a three-dimensional harmonic model and collinear model of He–H2 scattering. We also show that the proposed algorithm can be symmetrically composed to obtain time-reversible integrators of an arbitrary even order. We observed 10 000-fold speedup by using the tenth-order instead of the second-order method to obtain a solution with a time discretization error below 10−9. Moreover, using the adaptive grid instead of the fixed grid resulted in a 64-fold reduction in the required number of grid points in the harmonic system and made it possible to simulate the He–H2 scattering for six times longer while maintaining reasonable accuracy. Applicability of the algorithm to high-dimensional quantum dynamics is demonstrated using the strongly anharmonic eight-dimensional Hénon–Heiles model.
One of the most accurate methods for solving the time-dependent Schrödinger equation uses a combination of the dynamic Fourier method with the split-operator algorithm on a tensor-product grid. To reduce the number of required grid points, we let the grid move together with the wavepacket but find that the naïve algorithm based on an alternate evolution of the wavefunction and grid destroys the time reversibility of the exact evolution. Yet, we show that the time reversibility is recovered if the wavefunction and grid are evolved simultaneously during each kinetic or potential step; this is achieved by using the Ehrenfest theorem together with the splitting method. The proposed algorithm is conditionally stable, symmetric, and time-reversible and conserves the norm of the wavefunction. The preservation of these geometric properties is shown analytically and demonstrated numerically on a three-dimensional harmonic model and collinear model of He–H2 scattering. We also show that the proposed algorithm can be symmetrically composed to obtain time-reversible integrators of an arbitrary even order. We observed 10 000-fold speedup by using the tenth-order instead of the second-order method to obtain a solution with a time discretization error below 10−9. Moreover, using the adaptive grid instead of the fixed grid resulted in a 64-fold reduction in the required number of grid points in the harmonic system and made it possible to simulate the He–H2 scattering for six times longer while maintaining reasonable accuracy. Applicability of the algorithm to high-dimensional quantum dynamics is demonstrated using the strongly anharmonic eight-dimensional Hénon–Heiles model.
Geometric integrators of the Schrödinger equation conserve exactly many invariants of the exact solution. Among these integrators, the split-operator algorithm is explicit and easy to implement but, unfortunately, is restricted to systems whose Hamiltonian is separable into kinetic and potential terms. Here, we describe several implicit geometric integrators applicable to both separable and nonseparable Hamiltonians and, in particular, to the nonadiabatic molecular Hamiltonian in the adiabatic representation. These integrators combine the dynamic Fourier method with the recursive symmetric composition of the trapezoidal rule or implicit midpoint method, which results in an arbitrary order of accuracy in the time step. Moreover, these integrators are exactly unitary, symplectic, symmetric, time-reversible, and stable and, in contrast to the split-operator algorithm, conserve energy exactly, regardless of the accuracy of the solution. The order of convergence and conservation of geometric properties are proven analytically and demonstrated numerically on a two-surface NaI model in the adiabatic representation. Although each step of the higher order integrators is more costly, these algorithms become the most efficient ones if higher accuracy is desired; a thousand-fold speedup compared to the second-order trapezoidal rule (the Crank-Nicolson method) was observed for a wavefunction convergence error of 10−10. In a companion paper [J. Roulet, S. Choi, and J. Vaníček, J. Chem. Phys. 150, 204113 (2019)], we discuss analogous, arbitrary-order compositions of the split-operator algorithm and apply both types of geometric integrators to a higher-dimensional system in the diabatic representation.
Abstract— The fundamental possibility of a numerical solution of the time-dependent Schrödinger equation for a wave function defined on a non-uniform coordinate grid by the spectral method using the fast Fourier transform algorithm is discussed. The method is based on reducing the non-uniform coordinate grid to a uniform one by a non-linear transformation of coordinates and approximating the obtained evolution operator using the Lie–Trotter–Suzuki product formula. Algorithms for the numerical solution of the first and second orders are developed.
wavepacket is a program to simulate the dynamics of a wave packet interacting with a time-dependent potential, with the time-dependent Schrödinger equation being solved on a one-, two-, or three-dimensional spatial grid using the split operator method. This new version fixes bugs present in the original program. Program Title: wavepacket CPC Library link to program files: https://doi.org/10.17632/gnxmkygxrw.1 Licensing provisions: GPLv3 Programming language: C (iso C99) Journal reference of previous version: Comput. Phys. Commun. 185 (2014) 407–414. Does the new version supersede the previous version?: Yes Reasons for the new version: Bug fixes. Summary of revisions: This new version fixes bugs in the program wavepacket [1]. In the original version, the function print_wf_bin, used (optionally) to write out the wave function as a binary file, contained a bug that would result in only part of the full wave function being written out when the program was run in parallel. Additionally, negative values of the parameters nt and/or nprint could lead to the program being stuck in an infinite loop. Additional checks have been added to ensure that both these parameters are positive (or zero). Nature of problem: Solves the time-dependent Schrödinger equation for a single particle interacting with a time-dependent potential. Solution method: The wave function is described by its value on a spatial grid and the evolution operator is approximated using the split-operator method [2,3], with the kinetic energy operator calculated using a Fast Fourier Transform. Additional comments including restrictions and unusual features: Simulation can be in one, two, or three dimensions. The fftw [4] implementation of the fast Fourier transform is used, so the corresponding library must be installed. Serial and parallel versions are compiled from the same source files. The parallel version is based on the Message-Passing Interface (mpi) [5]. C.M. Dion, A. Hashemloo, G. Rahali, Comput. Phys. Commun. 185 (2014) 407–414. M.D. Feit, J.A. Fleck, Jr., A. Steiger, Solution of the Schrödinger equation by a spectral method, J. Comput. Phys. 47 (1982) 412–433. M.D. Feit, J.A. Fleck, Jr., Solution of the Schrödinger equation by a spectral method II: Vibrational energy levels of triatomic molecules, J. Chem. Phys. 78 (1) (1983) 301–308. http://www.fftw.org http://www.mpi-forum.org
The split-operator algorithm combined with the dynamic Fourier method is a popular approach for quantum dynamics in systems with separable kinetic and potential energies because it is explicit, easy to implement, and preserves many geometric invariants of the exact evolution. Although the standard version has only second-order accuracy, it can be recursively symmetrically composed to obtain geometric integrators of arbitrary even order of accuracy in the time step. We have implemented such integrators for nonadiabatic quantum dynamics in the diabatic basis using the triple-jump, Suzuki-fractal, and several non-recursive optimal compositions of the second-order split-operator algorithm. Because the automatically generated splitting coefficients are redundant, we reduce the computational cost by pruning these coefficients and lower memory requirements by identifying unique coefficients. The composed integrators are exactly unitary, symplectic, symmetric, time-reversible, and stable regardless of the accuracy of the solution. The order of convergence and preservation of geometric properties are justified analytically and confirmed numerically on a two-surface NaI model in the diabatic representation. Efficiency of higher-order algorithms is sometimes questioned because the number of splitting steps grows exponentially with the order. Yet, we show that, if higher accuracy is desired, higher-order integrators become the most efficient; a 600-fold speedup compared to the second-order split-operator method and a 17000-fold speedup compared to the trapezoidal rule (the Crank-Nicolson method) was observed for convergence error of 10−10.
A matrix inverse free method to solve time-dependent Schrodinger equation is presented. The method is not subject to form of Hamiltonian and adopting real space grid system such as structured and unstructured grid, and achieves the order N algorithm even if we adopt unstructured grid systems to adjust complex regions. We have demonstrated some error evaluation problems that have exact solution and application in which it has non-Schrodinger type Hamiltonian, and we have successfully simulated the time-evolution of wave function for all of the case.
While scattering phenomena are discussed in most quantum mechanics textbooks, the analysis is usually performed using time-independent techniques. These techniques have certainly been successful in predicting the results of scattering events but fail to present an intuitive picture of what happens to the physical system during a collision. A remedy to this situation could be affected by the introduction of time-dependent techniques, utilizing the conceptually simple idea of a wavefunction evolving over time, but these methods have traditionally been computationally unmanageable. This paper describes a novel numerical technique, the split-operator Fourier transform method, which greatly reduces the amount of computational effort required to accurately solve the time-dependent Schrödinger equation. This method is applied to the Rutherford scattering problem in two dimensions and it is shown that the method is not only efficient enough to yield results in less than 15 min when implemented on a standard personal computer, but is also simple enough to be employed by an advanced undergraduate. The method also has the advantage of reinforcing the quantum-mechanical formalism in the student’s mind; both the evolution operator and the complementary relationship between operators and wavefunctions in coordinate and momentum space are exploited. The results of the calculations are presented as color graphs portraying the magnitude and phase of wave packets as they evolve through time.