Numerical solver of the time-dependent Schrodinger equation with Coulomb singularities

This paper addresses a very fundamental and important problem in the numerical analysis of atomic and molecular systems: How to discretize Hamiltonians with divergent potential terms, such as Coulomb singularities. At the point of a Coulomb singularity, the wave function cannot be described by a Taylor series expansion, which results in problems when standard discretization schemes are used. We propose using the known asymptotic form of the wave function near the singularity instead of the (nonexistent) Taylor series. This principle, namely discretization by asymptotic behavior correspondence (ABC), is employed in this paper for obtaining grid-discretizations for the Coulomb potential in Cartesian, cylindrical and spherical coordinate systems. We show that computations with the ABC discretization are faster and more precise than with a naive discretization by orders of magnitude. The ABC discretization is well suited for the standard numerical time propagators, such as the Crank-Nicholson, Peaceman-Rachford, and leapfrog schemes. We use the latter, since it is faster and has the same order of accuracy. The leapfrog scheme is generalized to allow absorbing potentials at the grid boundaries.