Computing a vibrational spectrum numerically with a full-quantum treatment often involves expanding the unknown wavefunction as a linear combination of basis functions and solving the time-independent Schrödinger equation to obtain the expansion coefficients. This procedure yields a matrix eigenvalue problem that must be solved using linear algebra algorithms.

Collocation provides a formulation of the matrix eigenvalue problem by demanding that the Schrödinger equation be satisfied at collocation points in space. The collocation matrix eigenvalue problem is \[ \mathbf{B}^{-1}(\mathbf{T+VB})\mathbf{U} = \mathbf{UE}\] \[\mathbf{B}_{\alpha i} = b_i(\boldsymbol{q}_\alpha), \quad\mathbf{T}_{\alpha i}=T b_i(\boldsymbol{q})\Big|_{\boldsymbol{q} = \boldsymbol{q}_\alpha}, \quad\mathbf{V}_{\alpha\beta} = \delta_{\alpha\beta}V(\boldsymbol{q}_\alpha)~.\]

\(b_i(\boldsymbol{q})\) is a basis function, \(\boldsymbol{q}_\alpha\) a collocation point, \(T\) is the kinetic energy operator, and \(V(\boldsymbol{q})\) is the Born-Oppenheimer potential energy surface. The advantage of collocation over more common methods, such as variational methods, is greater flexibility. Collocation facilitates the use of optimized basis functions and coordinates. There is more leeway in the choice of points since no quadrature is needed.

In this work we show that collocation can be used with any choice of internal coordinates using a numerical kinetic energy operator, and can be combined with 1D contracted stretch, bend and torsion functions to efficiently compute the vibrational spectrum of flexible molecules.

The size of the basis set used can be large, thus iterative eigensolvers are preferred when computing a vibrational spectrum. They compute eigenvalues by means of matrix-vector products, they are more efficient than direct eigensolvers when only a portion of the spectrum is needed, and they avoid storing a prohibitively large matrix. Collocation yields non-symmetric matrices, and standard iterative methods for non-symmetric eigenvalue problems have unappealing computer memory requirements when a large number of eigenvalues is sought.

We present a new, low-storage iterative method based on the block Sakurai-Sugiura (SS) contour integration method to extract many eigenvalues from a single series of matrix-vector products. The method uses a multi-shift quasi-minimal residual (QMR) procedure to iteratively compute the quantities needed to use the SS method in any portion of the spectrum from a single series of matrix-vector products, enabling a systematic search of the eigenvalues.