# Iterative Methods for Computing Vibrational Spectra

**F**) matrix. The second uses Smolyak quadrature and a pruned basis. The third uses a tensor rank reduction scheme.

Next Article in Journal

Next Article in Special Issue

Next Article in Special Issue

Previous Article in Journal

Chemistry Department, Queen’s University, Kingston, ON K7L 3N6, Canada

Received: 21 December 2017
/
Revised: 10 January 2018
/
Accepted: 11 January 2018
/
Published: 16 January 2018

(This article belongs to the Special Issue Computational Spectroscopy)

I review some computational methods for calculating vibrational spectra. They all use iterative eigensolvers to compute eigenvalues of a Hamiltonian matrix by evaluating matrix-vector products (MVPs). A direct-product basis can be used for molecules with five or fewer atoms. This is done by exploiting the structure of the basis and the structure of a direct product quadrature grid. I outline three methods that can be used for molecules with more than five atoms. The first uses contracted basis functions and an intermediate (**F**) matrix. The second uses Smolyak quadrature and a pruned basis. The third uses a tensor rank reduction scheme.

Effective numerical methods for solving the time-independent Schroedinger equation to compute vibrational spectra of polyatomic molecules have been developed in the last thirty years [1,2,3,4,5]. They are important when approximations, often based on perturbation theory, are not accurate enough. Almost all methods begin by choosing a basis in which to represent both the wavefunctions and the Hamiltonian and then solve a linear algebra problem. These two basic tasks are not independent: a basis with structure favours iterative linear algebra methods (vide infra). Computing vibrational spectra is useful because it helps experimentalists to assign and interpret measured spectra.

In this article, I present a subjective review of several methods for solving the time-independent Schroedinger equation to calculate vibrational spectra [4,6,7]. It is possible to generalize the methods I describe so that they can also be used to compute ro-vibrational spectra [8,9,10,11,12,13,14]. All of the methods presented here obtain solutions to the Schroedinger equation from a space built by evaluating matrix-vector products (MVPs) and are called iterative methods [15]. I shall ignore Multimode-type methods (MM) [5,16,17,18,19,20], that work when the potential energy surface (PES) is a sum of terms that depend on a subset of the coordinates [17,21] (denoted an MM representation) and work quite well for semi-rigid molecules for which normal coordinates are appropriate. Although widely used, I shall also ignore multiconfiguration time-dependent Hartree (MCTDH) methods [22,23]. They can be used with a block power method [24], an “improved relaxation” method, [25,26,27], or a block Lanczos method [28] to compute accurate vibrational energy levels. However, improved relaxation, the most popular MCTDH approach for calculating spectra, converges poorly if the density of states is high and therefore cannot be used to compute a large number of levels of a large molecule [29].

When using iterative methods, it is better not to calculate a Hamiltonian matrix. Many calculations are done with a basis so large that it would not be possible to store the Hamiltonian matrix in memory. Iterative methods require the evaluation of MVPs. How is it possible to compute MVPs without building a matrix representing the Hamiltonian? I shall first outline ideas that make it possible to use a product basis to evaluate MVPs without a Hamiltonian matrix. They exploit the structure of the basis, the quadrature grid, and the kinetic energy operator (KEO). For molecules with more than five atoms, vectors representing wavefunctions and the vector representing the PES on the quadrature grid [4,30,31] are so large that they require too much memory. In the rest of this chapter, I therefore review methods that obviate the need to store large vectors. The first method (Section 4) uses a contracted basis. To make the contracted basis method useful, it is essential that it be possible to evaluate MVPs in the contracted basis without transforming to a huge product grid. The second method (Section 5) uses a pruned basis and a pruned grid. Pruning significantly reduces the size of the largest vectors one must store. The third method (Section 6) builds a basis from MVPs by using tensor rank reduction.

When there are D vibrational coordinates, a direct product basis function is
where the indices $\left\{{n}_{k}\right\}$ are independent and ${n}_{c}=0,1,\dots ,{n}_{c}^{max}$. ${\varphi}_{{n}_{c}}\left({q}_{c}\right)$ is a 1D basis function for coordinate c. If ${n}_{c}^{max}=n\phantom{\rule{3.33333pt}{0ex}}\phantom{\rule{3.33333pt}{0ex}}\forall c$, then the direct product basis set has ${n}^{D}$ functions. The univariate functions are often ${\varphi}_{k}\left(x\right)={h}_{k}^{-1/2}{\left[w\left(x\right)\right]}^{1/2}{p}_{k}\left(z\right)$, where z is a function of x, ${p}_{k}\left(z\right)$ is a classical orthogonal polynomial, $w\left(x\right)$ is the corresponding weight function, and ${h}_{k}$ is a normalization factor. Such a basis is usually called a variational basis representation (VBR) [4,30].

$${\mathsf{\Phi}}_{{n}_{1},{n}_{2},\dots ,{n}_{D}}={\varphi}_{{n}_{1}}\left({q}_{1}\right){\varphi}_{{n}_{2}}\left({q}_{2}\right)\dots {\varphi}_{{n}_{D}}\left({q}_{D}\right),\phantom{\rule{3.33333pt}{0ex}}$$

Although there are problems for which a VBR basis is best, it is sometimes advantageous to use a discrete variable representation (DVR) basis [4,30,31,32]. In 1D, a standard DVR basis is a set of orthogonal but localized functions that spans the same space as a set of orthogonal de-localized functions, ${\varphi}_{k}\left(x\right)$. The 1D DVR Hamiltonian matrix eigenvalue problem is
where $\mathbf{K}$ is an exact kinetic matrix in a basis of ${\varphi}_{n}\left(q\right)$ (VBR) functions and ${\mathbf{V}}^{\mathbf{FBR}}$ is either a product or a quadrature approximation for the exact potential matrix [4]. One way to obtain the transformation matrix $\mathbf{T}\phantom{\rule{3.33333pt}{0ex}}$ is to diagonalize the matrix representing x in the VBR,
where $\mathbf{x}$ is the matrix representing x in the ${\varphi}_{n}\left(q\right)$ basis and $\mathbf{X}$ is a diagonal matrix whose nonzero values are eigenvalues [33]. Equation (2) can be written
where ${\mathbf{V}}^{\mathbf{diag}}$ is a diagonal matrix whose diagonal elements are values of the potential at the quadrature (DVR) points. A potential optimised DVR (PO-DVR) [34,35] is made from 1D basis functions that are solutions of 1D Schroedinger equations.

$${\mathbf{T}}^{\mathbf{T}}(\mathbf{K}+{\mathbf{V}}^{\mathbf{FBR}})\mathbf{T}\mathbf{U}=\mathbf{U}\mathbf{E}\phantom{\rule{3.33333pt}{0ex}},$$

$$\mathbf{x}\mathbf{T}=\mathbf{T}\mathbf{X}\phantom{\rule{3.33333pt}{0ex}},$$

$$({\mathbf{T}}^{\mathbf{T}}\mathbf{K}\mathbf{T}+{\mathbf{V}}^{\mathbf{diag}})\mathbf{U}=\mathbf{U}\mathbf{E}\phantom{\rule{3.33333pt}{0ex}},$$

Direct product DVR and VBR bases are popular [4,36,37,38,39,40]. Although direct product bases are huge, they can be used by exploiting the structure of the basis to efficiently evaluate the MVPs required to use iterative methods. The Lanczos and filter diagonalisation methods are popular iterative methods for solving the time-independent Schroedinger equation [40,41,42,43,44,45,46,47]. For a molecule with as many as five atoms, they make it possible, even with a direct product basis, to solve the vibrational Schroedinger equation with a general PES. The key ideas have been reviewed several times [4,5,36,37]. They are all based on doing sums sequentially [40,48].

In a direct product DVR, potential MVPs are trivial because the potential matrix is diagonal. When the KEO is a sum of products (SOPs) , with g terms each with D factors,
then kinetic MVPs can be efficiently evaluated by doing sums sequentially,
where ${h}_{{n}_{k}^{\prime},{n}_{k}}^{(k,l)}$ is an element of the $n\times n$ matrix representation of the factor ${\widehat{h}}^{(k,l)}\left({q}_{k}\right)$. Matrix elements of the full KEO are never computed.

$$\widehat{K}=\sum _{l=1}^{g}\prod _{k=1}^{D}{\widehat{h}}^{(k,l)}\left({q}_{k}\right),$$

$$\sum _{l=1}^{g}\sum _{{n}_{1}}{h}_{{n}_{1}^{\prime},{n}_{1}}^{(1,l)}\sum _{{n}_{2}}{h}_{{n}_{2}^{\prime},{n}_{2}}^{(2,l)}\dots \sum _{{n}_{D}}{h}_{{n}_{D}^{\prime},{n}_{D}}^{(f,l)}{u}_{{n}_{1},{n}_{2},\dots ,{n}_{D}}={u}_{{n}_{1}^{\prime},{n}_{2}^{\prime},\dots ,{n}_{D}^{\prime}}^{\prime}\phantom{\rule{3.33333pt}{0ex}},$$

If there are important singularities in the KEO, then a VBR basis is better than a DVR basis [49]. At an important singularity, the KEO is singular and vibrational wavefunctions have significant amplitude. In general, singularities occur whenever one coordinate takes a limiting value and another is undefined [50]. In a VBR basis, it is possible to evaluate the potential MVP by doing sums sequentially [48,49]. This enables one to avoid calculating potential matrix elements, which would require computing many-dimensional integrals. Consider a 2D example. The matrix–vector product is
where

$$\sum _{{n}_{1}}\sum _{{n}_{2}}{V}_{{n}_{1}^{\prime}{n}_{2}^{\prime},{n}_{1}{n}_{2}}{u}_{{n}_{1},{n}_{2}}={u}_{{n}_{1}^{\prime},{n}_{2}^{\prime}}^{\prime}\phantom{\rule{3.33333pt}{0ex}},$$

$${V}_{{n}_{1}^{\prime}{n}_{2}^{\prime},{n}_{1}{n}_{2}}=\int d{q}_{1}d{q}_{2}{\varphi}_{{n}_{1}^{\prime}}\left({q}_{1}\right){\varphi}_{{n}_{2}^{\prime}}\left({q}_{2}\right)V({q}_{1},{q}_{2}){\varphi}_{{n}_{1}}\left({q}_{1}\right){\varphi}_{{n}_{2}}\left({q}_{2}\right)\phantom{\rule{3.33333pt}{0ex}}.$$

In terms of $\mathbf{T}$ matrices (see Equation (3)) ,

$${V}_{{n}_{1}^{\prime}{n}_{2}^{\prime},{n}_{1}{n}_{2}}\approx \sum _{\alpha}\sum _{\beta}{\left(T\right)}_{{n}_{1}^{\prime},\alpha}{\left(T\right)}_{{n}_{2}^{\prime},\beta}V({\left({q}_{1}\right)}_{\alpha},{\left({q}_{2}\right)}_{\beta}){({T}^{\u2020})}_{\alpha ,{n}_{1}}{({T}^{\u2020})}_{\beta ,{n}_{2}}.$$

The matrix–vector product can be written,
and evaluated by doing sums sequentially,

$$\sum _{{n}_{1}}\sum _{{n}_{2}}\sum _{\alpha}\sum _{\beta}{\left(T\right)}_{{n}_{1}^{\prime},\alpha}{\left(T\right)}_{{n}_{2}^{\prime},\beta}V({\left({q}_{1}\right)}_{\alpha},{\left({q}_{2}\right)}_{\beta}){({T}^{\u2020})}_{\alpha ,{n}_{1}}{({T}^{\u2020})}_{\beta ,{n}_{2}}{u}_{{n}_{1},{n}_{2}}={u}_{{n}_{1}^{\prime},{n}_{2}^{\prime}}^{\prime}$$

$$\sum _{\alpha}{\left(T\right)}_{{n}_{1}^{\prime},\alpha}\sum _{\beta}{\left(T\right)}_{{n}_{2}^{\prime},\beta}V({\left({q}_{1}\right)}_{\alpha},{\left({q}_{2}\right)}_{\beta})\sum _{{n}_{1}}{({T}^{\u2020})}_{\alpha ,{n}_{1}}\sum _{{n}_{2}}{({T}^{\u2020})}_{\beta ,{n}_{2}}{u}_{{n}_{1},{n}_{2}}={u}_{{n}_{1}^{\prime},{n}_{2}^{\prime}}^{\prime}\phantom{\rule{3.33333pt}{0ex}}.$$

If there is an important singularity, good basis functions are always nondirect product functions which are products of functions of the coordinate which is undefined and the coordinate which takes a limiting value, with a shared index.

The ideas of this section make it possible to compute vibrational spectra without computing and storing a Hamiltonian matrix. However, for molecules with more than five atoms, the memory cost of storing vectors in a direct product basis is prohibitive. For molecules with more than five atoms, it is necessary to introduce other ideas to reduce the memory cost of calculations.

To include information about coupling in the basis functions, it is common to use basis functions that are products of factors that depend on more than one coordinate. I shall call the multi-dimensional factors contracted basis functions. It is important to devise good algorithms for evaluating MVPs in a contracted basis. Contracted bases are necessarily more complicated, i.e., they have less structure, and it is structure that is exploited to evaluate MVPS efficiently. An important advantage that contracted bases have is the reduced spectral range of the contracted-basis Hamiltonian matrix. Reducing the spectral range decreases the number of MVPs required to compute eigenvalues.

For molecules with more than three atoms, it is best to use contracted functions diagonalizing matrices that represent the Hamiltonian with one or more coordinates fixed. The basis functions are direct products of functions of different coordinates or groups of coordinates [6,51,52,53].

The most obvious way to evaluate MVPs is to transform from the contracted basis to a primitive basis, in which the contracted functions are determined. Computing matrix elements of the potential in the primitive basis requires storing the potential on a large (direct product) grid of points. Because the grid is huge, this is impractical for molecules with more than five atoms. An alternative is to store an intermediate matrix [6,54]. To explain how this is done, consider a ($J=0$) Hamiltonian in polyspherical coordinates [55,56,57]
with
$\theta $ represents all of the bend coordinates and r represents all of the stretch coordinates. The functions ${B}_{i}\left(r\right)$ and the operators ${T}_{b}^{\left(i\right)}\left(\theta \right)$ are known [55,56,58]. One constructs contracted bend functions from a Hamiltonian obtained by fixing the stretch coordinates at some reference geometry and contracted stretch functions from a Hamiltonian obtained by fixing all the bend coordinates at reference values. Products of the bend contracted functions and stretch contracted functions are the final basis functions.

$$H={T}_{ben}(\theta ,r)+{T}_{str}\left(r\right)+V(\theta ,r)$$

$$\begin{array}{ccc}\hfill {T}_{ben}(\theta ,r)& =& \sum _{i}{B}_{i}\left(r\right){T}_{b}^{\left(i\right)}\left(\theta \right)\hfill \\ \hfill {T}_{str}\left(r\right)& =& \sum _{i}\frac{-1}{2{\mu}_{i}}\frac{{\partial}^{2}}{\partial {r}_{i}^{2}}\phantom{\rule{3.33333pt}{0ex}}.\hfill \end{array}$$

The reduced-dimension Hamiltonian for the bend contraction is,

$${H}^{\left(b\right)}={T}_{ben}(\theta ,{r}_{e})+V(\theta ,{r}_{e}).$$

Its wavefunctions are denoted by
and the energies by ${E}_{b}$. The ${f}_{l}$ are primitive bend basis functions (l is a composite index) and the number of retained bend wavefunctions is denoted by ${n}_{b}$. Similarly, the reduced-dimension Hamiltonian for the stretch contraction is,
with the wavefunctions denoted by,
and the energies by ${E}_{s}$. The ${g}_{\alpha}$ are primitive DVR stretch basis functions ($\alpha $ is a composite index representing a multidimensional DVR function) and the number of retained stretch wavefunctions is denoted by ${n}_{s}$. ${\theta}_{e}$ and ${r}_{e}$ represent reference (often equilibrium) values of all the bend coordinates and all the stretch coordinates. The final basis is a product of the retained stretch and bend eigenfunctions

$${X}_{b}\left(\theta \right)=\sum _{l}{C}_{lb}{f}_{l}\left(\theta \right)$$

$${H}^{\left(s\right)}={T}_{str}\left(r\right)+V({\theta}_{e},r).$$

$${Y}_{s}\left(r\right)=\sum _{\alpha}{D}_{\alpha s}{g}_{\alpha}\left(r\right)$$

$$|bs\rangle =|{X}_{b}\rangle |{Y}_{s}\rangle \phantom{\rule{3.33333pt}{0ex}}.$$

The full Hamiltonian is
where
and
with

$$H={H}^{\left(b\right)}+{H}^{\left(s\right)}+\Delta T+\Delta V$$

$$\Delta V(\theta ,r)=V(\theta ,r)-V(\theta ,{r}_{e})-V({\theta}_{e},r)$$

$$\Delta T=\sum _{i}\Delta {B}_{i}\left(r\right){T}_{b}^{\left(i\right)}\left(\theta \right)$$

$$\Delta {B}_{i}\left(r\right)={B}_{i}\left(r\right)-{B}_{i}\left({r}_{e}\right)\phantom{\rule{3.33333pt}{0ex}}.$$

In the contracted basis, MVPs for $\Delta T$ and ${H}^{\left(b\right)}+{H}^{\left(s\right)}$ are easy [6].

If one uses an finite basis representation (FBR) primitive bend basis and a DVR primitive stretch basis, a matrix element of $\Delta V$ in the product contracted basis is,

$$\begin{array}{ccc}\hfill \langle {b}^{\prime}{s}^{\prime}|\Delta V(\theta ,r)|bs\rangle & =& \sum _{\genfrac{}{}{0pt}{}{{l}^{\prime}l}{\alpha}}{D}_{\alpha {s}^{\prime}}{C}_{{l}^{\prime}{b}^{\prime}}\langle {l}^{\prime}|\Delta V(\theta ,{r}_{\alpha})|l\rangle {C}_{lb}{D}_{\alpha s}\phantom{\rule{3.33333pt}{0ex}}.\hfill \end{array}$$

This may be re-written
where I have introduced an **F** matrix [6] defined by,

$$\begin{array}{c}\hfill \sum _{\alpha}{F}_{{b}^{\prime}b,\alpha}{D}_{\alpha {s}^{\prime}}{D}_{\alpha s}\phantom{\rule{3.33333pt}{0ex}},\end{array}$$

$${F}_{{b}^{\prime},b,\alpha}=<{b}^{\prime}|\Delta V(\theta ,{r}_{\alpha})|b>=\sum _{{l}^{\prime}l}{C}_{{l}^{\prime}{b}^{\prime}}{C}_{lb}\langle {l}^{\prime}|\Delta V(\theta ,{r}_{\alpha})|l\rangle \phantom{\rule{3.33333pt}{0ex}}.$$

The integral $\langle {l}^{\prime}|\Delta V(\theta ,{r}_{\alpha})|l\rangle $ is computed with quadrature. The ${F}_{{b}^{\prime}b,\alpha}$ elements are calculated (in parallel) and stored before matrix vector products are evaluated. The $\Delta V$ matrix–vector product
is done as follows:

$${u}_{{b}^{\prime}{s}^{\prime}}^{\prime}=\sum _{bs}\langle {b}^{\prime}{s}^{\prime}|\Delta V|bs\rangle {u}_{bs}\phantom{\rule{3.33333pt}{0ex}},$$

$$\begin{array}{ccc}\hfill {u}_{b\alpha}^{\left(1\right)}& =& \sum _{s}{D}_{\alpha s}{u}_{bs}\hfill \\ \hfill {u}_{{b}^{\prime}\alpha}^{\left(2\right)}& =& \sum _{b}{F}_{{b}^{\prime}b\alpha}{u}_{b\alpha}^{\left(1\right)}\hfill \\ \hfill {u}_{{b}^{\prime}{s}^{\prime}}^{\prime}& =& \sum _{\alpha}{D}_{\alpha {s}^{\prime}}{u}_{{b}^{\prime}\alpha}^{\left(2\right)}\phantom{\rule{3.33333pt}{0ex}}.\hfill \end{array}$$

The idea of reducing the memory cost of contracted-basis calculations by storing a matrix representation of $\Delta V$ was used in [54], where only the bend basis was contracted. The full power of the method is realized only when both stretch and bend bases are contracted [6,59,60,61,62,63,64]. Recently, similar ideas were used for Cl${}^{-}$ -H${}_{2}$O [65]. Yu has studied several molecules using similar ideas [66,67,68,69]. As in [54], he contracts only the bend part.

In this section, I present an alternative method for computing vibrational spectra of molecules with more than five atoms. In contrast to the idea of using contracted basis functions, it uses univariate basis functions, but uses only selected products, i.e., the basis of Equation (1) is pruned by removing functions that are deemed unimportant. This makes it possible to obviate the need to store large vectors.

It seems clear that one should discard basis functions that are not necessary. Many authors have implemented basis pruning strategies [16,18,19,53,62,63,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84,85,86]. Although pruning has the obvious advantage that it decreases the size of the vectors one must store and the spectral range of the Hamiltonian matrix, if one uses an iterative method, it complicates the evaluation of MVPs. A pruned basis necessarily has less structure than a direct product basis. In this section, I shall discuss how to evaluate MVPs when the pruning strategy retains some product structure.

Pruning is more efficient when used with a FBR and not a DVR [40]. The simplest VBR pruning condition is ${n}_{1}+\dots +{n}_{D}\le b$. The pruned basis is much smaller than the direct product basis. If ${n}_{c}=0,1,\dots ,b$ for $c=1,\dots ,D$ and $b=14$ then the size of the direct product is ∼6 $\times {10}^{11}$, for $D=10$; ∼4 $\times {10}^{17}$, for $D=15$; and ∼3 $\times {10}^{23}$, for $D=20$. On the other hand, if basis functions with ${n}_{1}+\dots +{n}_{D}>b=14$ are discarded, the basis size increase with D is less than linear: ∼2.0 $\times {10}^{6}$, for $D=10$, ∼7.7 $\times {10}^{7}$, for $D=15$, ∼1.4 $\times {10}^{9}$, for $D=20$. MVPs for the KEO in a pruned basis are straightforward. MVPs for the potential are only straightforward if a direct product quadrature grid is used [83], but if a direct product quadrature is used, one needs to store a potential vector about as large as the direct product vectors one avoids by pruning the basis. The most important advantage of pruning is therefore lost It is possible to find a nondirect product quadrature scheme that uses fewer points and to evaluate potential MVPs by doing sums sequentially. The ideas will be explained with the ${n}_{1}+\dots +{n}_{D}\le b$ pruning condition, but better pruning conditions will be briefly discussed at the end of this section.

The direct product quadrature is, in a sense, too good because it is so good that many matrix elements with basis functions removed by the pruning are also exact. A nondirect product Smolyak quadrature is better. It has far fewer points but maintains enough structure to allow efficient MVPs. For detail see [7,87,88]. To make a Smolyak quadrature, one needs a family of 1D quadrature rules for each coordinate. Quadratures are labelled by ${i}_{c}$, ${i}_{c}\phantom{\rule{3.33333pt}{0ex}}=\phantom{\rule{3.33333pt}{0ex}}1,2,\dots ,{i}_{c}^{max}$. The number of points in quadrature rule ${i}_{c}$ is ${m}_{c}\left({i}_{c}\right)$, where ${m}_{c}\left({i}_{c}\right)$ is a non-decreasing function of ${i}_{c}$. To evaluate MVPs efficiently, one needs to use points for which all points in rule ${i}_{c}-1$ are also in rule ${i}_{c}$. Such points are called nested. The standard way to write a Smolyak quadrature is $f({q}_{1},{q}_{2},{q}_{3},{q}_{4},{q}_{5},{q}_{6})({}^{1}w\left({q}_{1}\right)\phantom{\rule{3.33333pt}{0ex}}{}^{2}w\left({q}_{2}\right)\phantom{\rule{3.33333pt}{0ex}}{}^{3}w\left({q}_{3}\right)\phantom{\rule{3.33333pt}{0ex}}{}^{4}w\left({q}_{4}\right)\phantom{\rule{3.33333pt}{0ex}}{}^{5}w\left({q}_{5}\right)\phantom{\rule{3.33333pt}{0ex}}{}^{6}w\left({q}_{6}\right))$ is
where ${q}_{c}^{{k}_{c}}$ is a point in the quadrature labelled by ${i}_{c}$, ${}^{{i}_{c}}{w}_{{k}_{c}}$ is the corresponding weight and the 1D quadratures are designed to approximate
${C}_{{i}_{1},\dots {i}_{6}}^{smol}$ are coefficients; see [7]. H is increased until convergence is achieved. The union of the grids for which ${i}_{1}+\dots +{i}_{6}\le H$ is satisfied is called the Smolyak grid. The size of the Smolyak grid is orders of magnitude smaller than the direct product grid size.

$$\begin{array}{c}\hfill \sum _{{i}_{1}+\dots +{i}_{6}\le H}{C}_{{i}_{1},\dots {i}_{6}}^{smol}\sum _{{k}_{1}}^{{m}_{1}\left({i}_{1}\right)}\sum _{{k}_{2}}^{{m}_{2}\left({i}_{2}\right)}\sum _{{k}_{3}}^{{m}_{3}\left({i}_{3}\right)}\sum _{{k}_{4}}^{{m}_{4}\left({i}_{4}\right)}\sum _{{k}_{5}}^{{m}_{5}\left({i}_{5}\right)}\sum _{{k}_{6}}^{{m}_{6}\left({i}_{6}\right)}\\ \hfill {}^{{i}_{1}}{w}_{{k}_{1}}\phantom{\rule{3.33333pt}{0ex}}{}^{{i}_{2}}{w}_{{k}_{2}}\phantom{\rule{3.33333pt}{0ex}}{}^{{i}_{3}}{w}_{{k}_{3}}\phantom{\rule{3.33333pt}{0ex}}{}^{{i}_{4}}{w}_{{k}_{4}}\phantom{\rule{3.33333pt}{0ex}}{}^{{i}_{5}}{w}_{{k}_{5}}\phantom{\rule{3.33333pt}{0ex}}{}^{{i}_{6}}{w}_{{k}_{6}}f({q}_{1}^{{k}_{1}},{q}_{2}^{{k}_{2}},{q}_{3}^{{k}_{3}},{q}_{4}^{{k}_{4}},{q}_{5}^{{k}_{5}},{q}_{6}^{{k}_{6}})& \phantom{\rule{3.33333pt}{0ex}},\end{array}$$

$$\begin{array}{c}\hfill \int d{q}_{c}\phantom{\rule{3.33333pt}{0ex}}{}^{c}w\left({q}_{c}\right)f\left({z}_{c}\left({q}_{c}\right)\right)\phantom{\rule{3.33333pt}{0ex}}.\end{array}$$

It would be costly to use Equation (27) in the evaluation of MVPs because it would be necessary to evaluate the sum over ${i}_{1}+\dots +{i}_{6}\le H$ for each MVP. When the 1D quadrature rules are nested, one can [7] replace Equation (27) with
where
are “super weights” that are pre-computed [89]. ${N}_{c}$ is a maximum number of points for coordinate c [7]. ${N}_{c}$ depends on ${k}_{{c}^{\prime}}$ if $c>{c}^{\prime}$ and ${N}_{1}$ does not depend on ${k}_{1},\dots ,{k}_{D}$. Using the super weights, it is possible to evaluate a potential MVP by doing sums sequentially,
where ${T}_{nk}={h}_{k}^{-1/2}{p}_{k}\left(z\left({q}_{k}\right)\right)$. ${n}_{c}^{max}$ depends on ${n}_{{c}^{\prime}}$ if $c<{c}^{\prime}$.

$$\begin{array}{c}={\sum}_{{k}_{1}}^{{N}_{1}}{\sum}_{{k}_{2}}^{{N}_{2}}{\sum}_{{k}_{3}}^{{N}_{3}}{\sum}_{{k}_{4}}^{{N}_{4}}{\sum}_{{k}_{5}}^{{N}_{5}}{\sum}_{{k}_{6}}^{{N}_{6}}w({k}_{1},{k}_{2},{k}_{3},{k}_{4},{k}_{5},{k}_{6})\\ \times f({q}_{1}^{{k}_{1}},{q}_{2}^{{k}_{2}},{q}_{3}^{{k}_{3}},{q}_{4}^{{k}_{4}},{q}_{5}^{{k}_{5}},{q}_{6}^{{k}_{6}})& ,\end{array}$$

$$\begin{array}{c}\hfill w({k}_{1},\dots ,{k}_{6})=\sum _{{i}_{1}+\dots {i}_{6}\le H}{C}_{{i}_{1},\dots ,{i}_{6}}^{smol}{}^{{i}_{1}}{w}_{{k}_{1}}\dots {}^{{i}_{D}}{w}_{{k}_{6}},\end{array}$$

$$\begin{array}{c}\hfill {u}^{\prime}({n}_{6}^{\prime},{n}_{5}^{\prime},{n}_{4}^{\prime}{n}_{3}^{\prime},{n}_{2}^{\prime},{n}_{1}^{\prime})=\sum _{{k}_{1}=1}^{{N}_{1}}{T}_{{n}_{1}^{\prime}{k}_{1}}\sum _{{k}_{2}=1}^{{N}_{2}}{T}_{{n}_{2}^{\prime}{k}_{2}}\sum _{{k}_{3}=1}^{{N}_{3}}{T}_{{n}_{3}^{\prime}{k}_{3}}\\ \hfill \sum _{{k}_{4}=1}^{{N}_{4}}{T}_{{n}_{4}^{\prime}{k}_{4}}\sum _{{k}_{5}=1}^{{N}_{5}}{T}_{{n}_{5}^{\prime}{k}_{5}}\sum _{{k}_{6}=1}^{{N}_{6}}{T}_{{n}_{6}^{\prime}{k}_{6}}\\ \hfill w({k}_{1},{k}_{2},{k}_{3},{k}_{4},{k}_{5},{k}_{6})V({q}_{1}^{{k}_{1}},{q}_{2}^{{k}_{2}},{q}_{3}^{{k}_{3}},{q}_{4}^{{k}_{4}},{q}_{5}^{{k}_{5}},{q}_{6}^{{k}_{6}})\\ \hfill \sum _{{n}_{6}=0}^{{n}_{6}^{\mathrm{max}}}{T}_{{n}_{6}{k}_{6}}\sum _{{n}_{5}=0}^{{n}_{5}^{\mathrm{max}}}{T}_{{n}_{5}{k}_{5}}\sum _{{n}_{4}=0}^{{n}_{4}^{\mathrm{max}}}{T}_{{n}_{4}{k}_{4}}\sum _{{n}_{3}=0}^{{n}_{3}^{\mathrm{max}}}{T}_{{n}_{3}{k}_{3}}\sum _{{n}_{2}=0}^{{n}_{2}^{\mathrm{max}}}{T}_{{n}_{2}{k}_{2}}\sum _{{n}_{1}=0}^{{n}_{1}^{\mathrm{max}}}{T}_{{n}_{1}{k}_{1}}\\ \hfill u({n}_{6},{n}_{5},{n}_{4},{n}_{3},{n}_{2},{n}_{1})\phantom{\rule{3.33333pt}{0ex}},\end{array}$$

To use Equation (31), one first sums over ${n}_{1}$ to compute an intermediate vector $y{1}_{{k}_{1},{n}_{2},{n}_{3},{n}_{4},{n}_{5},{n}_{6}}$ and then sums over ${n}_{2}$ to compute an intermediate vector whose components are $y{2}_{{k}_{1},{k}_{2},{n}_{3},{n}_{4},{n}_{5},{n}_{6}}$ etc. At each step, the ${n}_{c}$ and ${k}_{c}$ indices are constrained among themselves. Everything is clearly explained in [90]. The Smolyak grid is a sum of smaller direct product grids and therefore has structure that makes it possible to evaluate MVPs by doing sums sequentially. It is also possible to do sums sequentially for any pruning condition of the form ${g}^{1}\left({n}_{1}\right)+\dots +{g}^{D}\left({n}_{D}\right)\le b$ [87]. Sometimes, ${g}^{c}\left({n}_{c}\right)={\alpha}_{c}{n}_{c}$ with ${\alpha}_{c}=\u230a\frac{{\omega}_{c}}{{\omega}_{lowest}}+0.5\u230b$ is a good choice. Often this choice can be improved. In general, basis functions with many non-zero indices for coordinates with large frequencies are unimportant. They can be pushed out of the basis by using ${G}^{c}\left({n}_{c}\right)>{\alpha}_{c}^{P}{n}_{c}$. In general, it is important to include in the basis product functions for which there are several non-zero indices for coordinates with small frequencies. Such functions are preferentially included in the basis by using ${G}^{c}\left({n}_{c}\right)<{\alpha}_{c}^{P}$. These are general guidelines which we have found useful, but they are not specific [87,88].

Contraction (Section 4) and pruning (Section 5) enable one to avoid storing vectors with as many elements as the direct product basis. For molecules with more than five atoms, this is essential. For example, if $D=12$ and $n=10$ then for a single vector one needs ∼8000 GB of memory to store a vector with ${n}^{D}$ components. In this section, I describe another approach for avoiding vectors with ${n}^{D}$ components. It $does$ use a direct product basis but exploits advantages of a SOPs PES [91,92,93]. The key idea is that in some cases, the ${n}^{D}$ coefficients, used to represent a function, can be computed from a much smaller set of numbers. For example, a product of functions of a single coordinate, ${\varphi}_{1}\left({q}_{1}\right){\varphi}_{2}\left({q}_{2}\right)\dots {\varphi}_{D}\left({q}_{D}\right)$, can be represented as
and it is only necessary to store $Dn$ numbers.

$$\sum _{{i}_{1}=1}^{n}{f}_{{i}_{1}}^{\left(1\right)}{\theta}_{{i}_{1}}^{1}\left({q}_{1}\right)\sum _{{i}_{2}=1}^{n}{f}_{{i}_{2}}^{\left(2\right)}{\theta}_{{i}_{2}}^{2}\left({q}_{2}\right)\dots \sum _{{i}_{D}=1}^{n}{f}_{{i}_{D}}^{\left(D\right)}{\theta}_{{i}_{D}}^{D}\left({q}_{D}\right)$$

We have developed computational methods that solve the Schroedinger equation by projecting into a basis of functions that are sums of products. Although the functions in this basis are SOPs, the basis is not a direct product basis. A single basis function with R terms is determined by only $RDn$ numbers. If D is large, this is much less than ${n}^{D}$. The key idea is to use basis functions that are sums of products of optimised factors. It works if the Hamiltonian is itself a SOPs. The SOPs basis functions are represented in a primitive product basis made from 1-D functions ${\theta}_{{i}_{j}}^{j}\left({q}_{j}\right)$ with ${i}_{j}=1,\dots ,{n}_{j}$ for each coordinate ${q}_{j}$. Any function can be expanded in this basis as

$$\mathsf{\Psi}({q}_{1},\dots ,{q}_{D})\simeq \sum _{{i}_{1}=1}^{{n}_{1}}\dots \sum _{{i}_{D}=1}^{{n}_{D}}{F}_{{i}_{1}{i}_{2}\dots {i}_{D}}\prod _{j=1}^{D}{\theta}_{{i}_{j}}^{j}\left({q}_{j}\right)\phantom{\rule{3.33333pt}{0ex}}.$$

The goal is to avoid explicitly introducing ${F}_{{i}_{1}{i}_{2}\dots {i}_{D}}$. This is possible if $\mathsf{\Psi}({q}_{1},\dots ,{q}_{D})$ is a SOPs. In that case,
where ${f}^{(\ell ,j)}$ is a one-dimensional vector associated with the ℓ-th term and coordinate j. The SOPs format for multidimensional functions is known as the canonical polyadic (CP) decomposition for tensors [94,95,96].

$${F}_{{i}_{1}{i}_{2}\dots {i}_{D}}=\sum _{\ell =1}^{R}\prod _{j=1}^{D}{f}_{{i}_{j}}^{(\ell ,j)}$$

Basis functions are made by applying the Hamiltonian. We have used a shifted block power method [91]. It is imperative that every basis vector be in the form of Equation (33). This is only the case if the Hamiltonian is of the form,
where ${h}_{kj}$ is a one-dimensional operator acting in a Hilbert space associated with coordinate ${q}_{j}$. PES can be forced into SOPs form by using, for example, potfit [5,97], multigrid potfit [98], or neural network methods [99,100,101].

$$H({q}_{1},\dots ,{q}_{D})=\sum _{k=1}^{T}\prod _{j=1}^{D}{h}_{kj}\left({q}_{j}\right),$$

When $\mathbf{H}$ is applied to a vector $\mathbf{F}$ to obtain a new vector ${\mathbf{F}}^{\prime}$, the number of terms in the vector increases. If there are T terms in $\mathbf{H}$, the rank (number of terms) of ${\mathbf{F}}^{\prime}$ is a factor of T larger than the rank of $\mathbf{F}$. All vectors have the form
where, for each term (ℓ) and each coordinate (j), ${\tilde{f}}_{{i}_{j}}^{(\ell ,j)}$ is a normalized 1-D vector, ${s}_{\ell}$ is a normalization coefficient, and ${n}_{j}$ is the number of basis functions for coordinate j. $\mathbf{H}$ can be applied to $\mathbf{F}$ by evaluating 1-D matrix–vector products,
where ${\left({\mathbf{h}}_{kj}\right)}_{{i}_{j},{i}_{j}^{\prime}}=\langle {\theta}_{{i}_{j}}^{j}|{h}_{kj}|{\theta}_{{i}_{j}^{\prime}}^{j}\rangle $.

$${F}_{{i}_{1}{i}_{2}\dots {i}_{D}}=\sum _{\ell =1}^{R}{s}_{\ell}\prod _{j=1}^{D}{\tilde{f}}_{{i}_{j}}^{(\ell ,j)}\phantom{\rule{4.pt}{0ex}}\mathrm{with}\phantom{\rule{4.pt}{0ex}}\sum _{{i}_{j}}^{{n}_{j}}{|{\tilde{f}}_{{i}_{j}}^{(\ell ,j)}|}^{2}=1\phantom{\rule{3.33333pt}{0ex}},$$

$${\left(\mathbf{H}F\right)}_{{i}_{1}^{\prime}\dots {i}_{D}^{\prime}}=\sum _{{i}_{1},{i}_{2},\dots ,{i}_{D}}\sum _{k=1}^{T}\prod _{{j}^{\prime}=1}^{D}{\left({\mathbf{h}}_{k{j}^{\prime}}\right)}_{{i}_{{j}^{\prime}}^{\prime}{i}_{{j}^{\prime}}}\sum _{\ell =1}^{R}\prod _{j=1}^{D}{s}_{\ell}{\tilde{f}}_{{i}_{j}}^{(\ell ,j)}.$$

$$=\sum _{k=1}^{T}\sum _{\ell =1}^{R}\prod _{j=1}^{D}\sum _{{i}_{j}}{\left({\mathbf{h}}_{kj}\right)}_{{i}_{j}^{\prime}{i}_{j}}{s}_{\ell}{\tilde{f}}_{{i}_{j}}^{(\ell ,j)},$$

To avoid having vectors with many terms, we must reduce the rank. To do this, we replace ${F}_{{i}_{1}{i}_{2}\dots {i}_{D}}^{\mathrm{old}}$
where ${R}_{\mathrm{new}}<{R}_{\mathrm{old}}$ and choose ${}^{\mathrm{new}}{\tilde{f}}_{{i}_{j}}^{(\ell ,j)}$ to minimize $\parallel {F}^{\mathrm{new}}-{F}^{\mathrm{old}}\parallel $. We use the same ${R}_{\mathrm{new}}$ for all reductions. An alternating least squares (ALS) algorithm described in [96] is used to carry out the reduction.

$$\begin{array}{ccc}\hfill {F}_{{i}_{1}{i}_{2}\dots {i}_{D}}^{\mathrm{old}}& =\hfill & \sum _{\ell =1}^{{R}_{\mathrm{old}}}{s}_{\ell}\prod _{j=1}^{D}{{}^{\mathrm{old}}\tilde{f}}_{{i}_{j}}^{(\ell ,j)}\hfill \\ & \u27f9\hfill & {F}_{{i}_{1}{i}_{2}\dots {i}_{D}}^{\mathrm{new}}=\sum _{\ell =1}^{{R}_{\mathrm{new}}}{s}_{\ell}\prod _{j=1}^{D}{{}^{\mathrm{new}}\tilde{f}}_{{i}_{j}}^{(\ell ,j)}\phantom{\rule{3.33333pt}{0ex}},\hfill \end{array}$$

In [91], it was demonstrated that these ideas work for a 20D Hamiltonian of coupled oscillators. A rank of only 20 (i.e., 20 terms in each of the basis functions) was sufficient to converge about 40 states of a 20D Hamiltonian. Related ideas were used successfully for molecules with as many as 10 atoms [91,92,93].

Iterative eigensolvers make it possible to calculate vibrational spectra without storing a Hamiltonian matrix. The Lanczos algorithm, filter diagonalization, and a re-started Lanczos or Arnoldi method available as ARPACK [102] are common iterative eigensolvers.

It is easiest to use an iterative eigensolver when the Hamiltonian is a SOP. A Taylor series potential is a SOPs. In many cases, the vibrational KEO is a SOPs. In normal coordinates, it is only a SOPs if one expands elements of the effective moment of inertia tensor [103] or sets them to zero (approximates the KEO). If the PES is not a SOPs, it can be massaged into SOPs form [101,104]. It is harder to use an iterative eigensolver when the Hamiltonian is not an SOPs. In this case, quadrature (or collocation) is used and the Hamiltonian matrix is usually not sparse. In a product basis, the cost of Hamiltonian MVPs scales as ${n}^{D+1}$ [40,48,49]. This favourable scaling is obtained by exploiting structure. Product basis/product grid methods are methods of first resort for molecules with four or five atoms. They also work extremely well for Van der Waals molecules when only the intermonomer coordinates are treated explicitly [8,11,12,105,106].

In this review article, I describe three methods that obviate the need to store vectors with as many components as the product basis. The first method uses basis functions that are eigenfunctions of a Hamiltonian obtained by setting a subset of the coordinates equal to reference values. In Section 4, a method is described for evaluating MVPs in a product contracted basis. It does not require storing vectors as large as the primitive product basis set. The key idea is to store an intermediate matrix, called the F matrix. The second method uses basis functions that are products of univariate functions. Functions deemed unimportant are removed from a direct-product basis by imposing a pruning condition. The pruning condition is chosen so that the pruned basis has structure. When used with a general PES, it is necessary to use the pruned basis in conjunction with a nondirect product quadrature grid that satisfies two requirements. It must have fewer points than the product quadrature and it must have sufficient structure to be able to evaluate MVPs by doing sums sequentially. A Smolyak quadrature satisfies both requirements. However, to use it, the quadrature must be written not as a sum over quadrature levels, as is usually the case, but as a sum over points (i.e., not Equation (27) but Equation (30)). The third method uses SOPs basis functions. It works only if the Hamiltonian is a SOPs. The basis functions are determined by reducing the rank of the vectors obtained from MVPs. With rank reduction methods, it is possible to compute vibrational spectra for molecules with more than a dozen atoms [107].

The research described in this paper was supported by the Canadian Natural Sciences and Engineering Research Council. Many excellent students and postdocs made important contributions to the development and implementation of the ideas described here. I am especially grateful to Gustavo Avila, Matthew Bramley, Arnaud Leclerc, and Xiao-Gang Wang.

The author declares no conflicts of interest.

- Schinke, R. Photodissociation Dynamics; Cambridge University Press: Cambridge, UK, 1993. [Google Scholar]
- Tannor, D.J. Introduction to Quantum Mechanics: A Time- Dependent Perspective; University Science Books: Sausalito, CA, USA, 2007. [Google Scholar]
- Kosloff, R. Time-dependent quantum-mechanical methods for molecular dynamics. J. Phys. Chem.
**1988**, 92, 2087–2100. [Google Scholar] [CrossRef] - Light, J.C.; Carrington, T., Jr. Discrete-variable representations and their utilization. Adv. Chem. Phys.
**2000**, 114, 263–310. [Google Scholar] - Bowman, J.M.; Carrington, T.; Meyer, H.-D. Variational quantum approaches for computing vibrational energies of polyatomic molecules. Mol. Phys.
**2008**, 106, 2145–2182. [Google Scholar] [CrossRef] - Wang, X.-G.; Carrington, T. New ideas for using contracted basis functions with a Lanczos eigensolver for computing vibrational spectra of molecules with four or more atoms. J. Chem. Phys.
**2002**, 117, 6923–6934. [Google Scholar] [CrossRef] - Avila, G.; Carrington, T., Jr. Nonproduct quadrature grids for solving the vibrational Schrödinger equation. J. Chem. Phys.
**2009**, 131, 174103. [Google Scholar] [CrossRef] [PubMed] - Dawes, R.; Wang, X.-G.; Jasper, A.; Carrington, T. Nitrous oxide dimer: a new potential energy surface and ro-vibrational spectrum of the non-polar isomer. J. Chem. Phys.
**2010**, 133, 134304:1–134304:14. [Google Scholar] [CrossRef] [PubMed] - Sarkar, P.; Poulin, N.; Carrington, T., Jr. Calculating rovibrational energy levels of a triatomic molecule with a simple Lanczos method. J. Chem. Phys.
**1999**, 110, 10269–10274. [Google Scholar] [CrossRef] - Szidarovszky, T.; Fabri, C.; Csaszar, A.G. The role of axis embedding on rigid rotor decomposition analysis of variational rovibrational wave functions. J. Chem. Phys.
**2012**, 136, 174112. [Google Scholar] [CrossRef] [PubMed] - Leforestier, C.; Braly, L.B.; Liu, K.; Elroy, M.J.; Saykally, R.J. Fully coupled six-dimensional calculations of the water dimer vibration-rotation-tunneling states with a split Wigner pseudo spectral approach. J. Chem. Phys.
**1997**, 106, 8527. [Google Scholar] [CrossRef] - Brown, J.; Wang, X.-G.; Dawes, R.; Carrington, T. Computational study of the rovibrational spectrum of (OCS)2. J. Chem. Phys.
**2012**, 136, 134306:1–134306:12. [Google Scholar] [CrossRef] [PubMed] - Wang, X.-G.; Carrington, T. Computing ro-vibrational levels of methane with internal vibrational coordinates and an Eckart frame. J. Chem. Phys.
**2013**, 138, 104106:1–104106:20. [Google Scholar] [CrossRef] [PubMed] - Wang, X.-G.; Carrington, T.; McKellar, A.R.W. Theoretical and experimental study of the rovibrational spectrum of He2-CO. J. Phys. Chem. A
**2009**, 113, 13331–13341, (Robert Field Festschrift). [Google Scholar] [CrossRef] [PubMed] - Golub, G.H.; van Loan, C.F. Matrix Computations; JHU Press: Baltimore, MD, USA, 2012; Volume 3. [Google Scholar]
- Carter, S.; Bowman, J.M.; Handy, N.C. Extensions and tests of “multimode”: A code to obtain accurate vibration/rotation energies of many-mode molecules. Theor. Chim. Acta
**1998**, 100, 191–198. [Google Scholar] [CrossRef] - Carter, S.; Culik, S.J.; Bowman, J.M. Vibrational self-consistent field method for many-mode systems: A new approach and application to the vibrations of CO adsorbed on Cu (100). J. Chem. Phys.
**1997**, 107, 10458–10469. [Google Scholar] [CrossRef] - Benoit, D.M. Fast vibrational self-consistent field calculations through a reduced mode–mode coupling scheme. J. Chem. Phys.
**2004**, 120, 562–573. [Google Scholar] [CrossRef] [PubMed] - Meier, P.; Neff, M.; Rauhut, G. Accurate Vibrational Frequencies of Borane and Its Isotopologues. J. Chem. Theory Comput.
**2011**, 7, 148–152. [Google Scholar] [CrossRef] [PubMed] - Rauhut, G. Configuration selection as a route towards efficient vibrational configuration interaction calculations. J. Chem. Phys.
**2007**, 127, 184109. [Google Scholar] [CrossRef] [PubMed] - Alis, O.; Rabitz, H. General Foundations of High Dimensional Model Representations. J. Math. Chem.
**1999**, 25, 197–233. [Google Scholar] - Beck, M.H.; Jäckle, A.; Worth, G.A.; Meyer, H.D. The multiconfiguration time-dependent Hartree (MCTDH) method: A highly efficient algorithm for propagating wavepackets. Phys. Rep.
**2000**, 324, 1–105. [Google Scholar] [CrossRef] - Manthe, U.; Meyer, H.-D.; Cederbaum, L. Wave-packet dynamics within the multiconfiguration Hartree framework: General aspects and application to NOCl. J. Chem. Phys.
**1992**, 97, 3199–3213. [Google Scholar] [CrossRef] - Manthe, U. The state averaged multiconfigurational time-dependent Hartree approach: Vibrational state and reaction rate calculations. J. Chem. Phys.
**2008**, 128, 064108. [Google Scholar] [CrossRef] [PubMed] - Meyer, H.-D.; le Quéré, F.; Léonard, C.; Gatti, F. Calculation and selective population of vibrational levels with the Multiconfiguration Time-Dependent Hartree (MCTDH) algorithm. Chem. Phys.
**2006**, 329, 179–192. [Google Scholar] [CrossRef] - Richter, F.; Gatti, F.; Léonard, C.; le Quéré, F.; Meyer, H.-D. Time-dependent wave packet study on trans-cis isomerization of HONO driven by an external field. J. Chem. Phys.
**2007**, 127, 164315. [Google Scholar] [CrossRef] [PubMed] - Doriol, L.J.; Gatti, F.; Iung, C.; Meyer, H.-D. Computation of vibrational energy levels and eigenstates of fluoroform using the multiconfiguration time-dependent Hartree method. J. Chem. Phys.
**2008**, 129, 224109. [Google Scholar] [CrossRef] [PubMed] - Wodraszka, R.; Manthe, U. Iterative Diagonalization in the Multiconfigurational Time-Dependent Hartree Approach: Ro-vibrational Eigenstates. J. Phys. Chem. A
**2013**, 117, 7246–7255. [Google Scholar] [CrossRef] [PubMed] - Vendrell, O.; Gatti, F.; Meyer, H.D. Full dimensional (15-dimensional) quantum-dynamical simulation of the protonated water dimer. II. Infrared spectrum and vibrational dynamics. J. Chem. Phys.
**2007**, 127, 184303. [Google Scholar] [CrossRef] [PubMed] - Light, J.C.; Hamilton, I.P.; Lill, J.V. Generalized discrete variable approximation in quantum-mechanics. J. Chem. Phys.
**1985**, 82, 1400–1409. [Google Scholar] [CrossRef] - Bac̆ić, Z.; Light, J.C. Theoretical Methods for Rovibrational States of Floppy Molecules. Annu. Rev. Phys. Chem.
**1989**, 40, 469–498. [Google Scholar] [CrossRef] - Wei, H.; Carrington, T. Discrete variable representations of complicated kinetic energy operators. J. Chem. Phys.
**1994**, 101, 1343–1360. [Google Scholar] [CrossRef] - Harris, D.O.; Engerholm, G.G.; Gwinn, W.D. Calculation of Matrix Elements for One-Dimensional Quantum-Mechanical Problems and the Application to Anharmonic Oscillators. J. Chem. Phys.
**1965**, 43, 1515–1517. [Google Scholar] [CrossRef] - Echave, J.; Clary, D.C. Potential optimized discrete variable representation. Chem. Phys. Lett.
**1992**, 190, 225–230. [Google Scholar] - Wei, H.; Carrington, T. The discrete variable representation for a triatomic Hamiltonian in bond length-bond angle coordinates. J. Chem. Phys.
**1992**, 97, 3029–3037. [Google Scholar] [CrossRef] - Carrington, T. Methods for calculating vibrational energy levels. Can. J. Chem.
**2004**, 82, 900–914. [Google Scholar] [CrossRef] - Csaszar, A.G.; Fabri, C.; Szidarovszky, T.; Matyus, E.; Furtenbacher, T.; Czako, G. The fourth age of quantum chemistry: molecules in motion. Phys. Chem. Chem. Phys.
**2012**, 14, 1085–1106. [Google Scholar] [CrossRef] [PubMed] - Matyus, E.; Czako, G.; Sutcliffe, B.T.; Csaszar, A.G. Vibrational energy levels with arbitrary potentials using the Eckart-Watson Hamiltonians and the discrete variable representation. J. Chem. Phys.
**2007**, 127, 084102. [Google Scholar] [CrossRef] [PubMed] - Yu, H.; Muckerman, J.T. A General Variational Algorithm to Calculate Vibrational Energy Levels of Tetraatomic Molecules. J. Mol. Spectrosc.
**2002**, 214, 11–20. [Google Scholar] [CrossRef] - Bramley, M.J.; Carrington, T. A general discrete variable method to calculate vibrational energy levels of three-and four-atom molecules. J. Chem. Phys.
**1993**, 99, 8519–8541. [Google Scholar] [CrossRef] - Wall, M.R.; Neuhauser, D. Extraction, through filter-diagonalization, of general quantum eigenvalues or classical normal mode frequencies from a small number of residues or a short-time segment of a signal. I. Theory and application to a quantum-dynamics model. J. Chem. Phys.
**1995**, 102, 8011–8022. [Google Scholar] [CrossRef] - Mandelshtam, V.A.; Taylor, H.S. Harmonic inversion of time signals and its applications. J. Chem. Phys.
**1997**, 107, 6756–6769. [Google Scholar] [CrossRef] - Iung, C.; Leforestier, C. Direct calculation of overtones: Application to the CD3H molecule. J. Chem. Phys.
**1995**, 102, 8453–8461. [Google Scholar] [CrossRef] - Le Quéré, F.; Leforestier, C. Quantum exact three-dimensional study of the photodissociation of the ozone molecule. J. Chem. Phys.
**1990**, 92, 247–253. [Google Scholar] [CrossRef] - McNichols, A.; Carrington, T. Vibrational energy levels of formaldehyde calculated from an internal coordinate Hamiltonian using the Lanczos algorithm. Chem. Phys. Lett.
**1993**, 202, 464–470. [Google Scholar] [CrossRef] - Huang, S.-W.; Carrington, T. A comparison of filter diagonalisation methods with the Lanczos method for calculating vibrational energy levels. Chem. Phys. Lett.
**1999**, 312, 311–318. [Google Scholar] [CrossRef] - Lee, S.; Chung, J.S.; Felker, P.M.; Cacheiro, J.L.; Fernández, B.; Bondo Pedersen, T.; Koch, H. Computational and experimental investigation of intermolecular states and forces in the benzene–helium van der Waals complex. J. Chem. Phys.
**2003**, 119, 12956. [Google Scholar] [CrossRef] - Manthe, U.; Koeppel, H. New method for calculating wave packet dynamics: Strongly coupled surfaces and the adiabatic basis. J. Chem. Phys.
**1990**, 93, 345–356. [Google Scholar] [CrossRef] - Bramley, M.J.; Tromp, J.W.; Carrington, T., Jr.; Corey, C.G. Efficient calculation of highly excited vibrational energy levels of floppy molecules: The band origins of H+ 3 up to 35,000 cm
^{−1}. J. Chem. Phys.**1994**, 100, 6175–6194. [Google Scholar] [CrossRef] - Sutcliffe, B.T. Coordinate Systems and Transformations. In Handbook of Molecular Physics and Quantum Chemistry; Wilson, S., Ed.; Wiley: Chichester, UK, 2003; Volume 1, Part 6, Chapter 31; pp. 485–500. [Google Scholar]
- Carter, S.; Handy, N.C. A variational method for the determination of the vibrational (J = 0) energy levels of acetylene, using a Hamiltonian in internal coordinates. Comput. Phys. Commun.
**1988**, 51, 49–58. [Google Scholar] [CrossRef] - Bramley, M.J.; Handy, N.C. Efficient calculation of rovibrational eigenstates of sequentially bonded four-atom molecules. J. Chem. Phys.
**1993**, 98, 1378. [Google Scholar] [CrossRef] - Yu, H.-G. An exact variational method to calculate vibrational energies of five atom molecules beyond the normal mode approach. J. Chem. Phys.
**2002**, 117, 2030. [Google Scholar] [CrossRef] - Bramley, M.J.; Carrington, T., Jr. Calculation of triatomic vibrational eigenstates: Product or contracted basis sets, Lanczos or conventional eigensolvers? What is the most efficient combination? J. Chem. Phys.
**1994**, 101, 8494. [Google Scholar] [CrossRef] - Gatti, F.; Iung, C.; Menou, M.; Justum, Y.; Nauts, A.; Chapuisat, X. Vector parametrization of the N-atom problem in quantum mechanics. I. Jacobi vectors. J. Chem. Phys.
**1998**, 108, 8804. [Google Scholar] [CrossRef] - Mladenović, M. Rovibrational Hamiltonians for general polyatomic molecules in spherical polar parametrization. I. Orthogonal representations. J. Chem. Phys.
**2000**, 112, 1070–1081. [Google Scholar] [CrossRef] - Chapuisat, X.; Belfhal, A.; Nauts, A. N-body quantum-mechanical Hamiltonians: Extrapotential terms. J. Chem. Phys.
**1991**, 149, 274–304. [Google Scholar] [CrossRef] - Gatti, F.; Iung, C. Exact and constrained kinetic energy operators for polyatomic molecules: The polyspherical approach. Phys. Rep.
**2009**, 484, 1–69. [Google Scholar] [CrossRef] - Wang, X.-G.; Carrington, T., Jr. A contracted basis-Lanczos calculation of vibrational levels of methane: Solving the Schrödinger equation in nine dimensions. J. Chem. Phys.
**2003**, 119, 101–117. [Google Scholar] [CrossRef] - Tremblay, J.C.; Carrington, T., Jr. Calculating vibrational energies and wave functions of vinylidene using a contracted basis with a locally reorthogonalized coupled two-term Lanczos eigensolver. J. Chem. Phys.
**2006**, 125, 094311. [Google Scholar] [CrossRef] [PubMed] - Wang, X.-G.; Carrington, T., Jr. Vibrational energy levels of CH5+. J. Chem. Phys.
**2008**, 129, 234102. [Google Scholar] [CrossRef] [PubMed] - Lee, H.-S.; Light, J.C. Molecular vibrations: Iterative solution with energy selected bases. J. Chem. Phys.
**2003**, 118, 3458. [Google Scholar] [CrossRef] - Lee, H.-S.; Light, J.C. Iterative solutions with energy selected bases for highly excited vibrations of tetra-atomic molecules. J. Chem. Phys.
**2004**, 120, 4626. [Google Scholar] [CrossRef] [PubMed] - Wang, X.-G.; Carrington, T., Jr. Contracted basis Lanczos methods for computing numerically exact rovibrational levels of methane. J. Chem. Phys.
**2004**, 121, 2937–2954. [Google Scholar] [CrossRef] [PubMed] - Wang, X.-G.; Carrington, T. Using monomer vibrational wavefunctions as contracted basis functions to compute rovibrational levels of an H
_{2}O-atom complex in full dimensionality. J. Chem. Phys.**2017**, 146, 104105:1–104105:15. [Google Scholar] [CrossRef] [PubMed] - Yu, H.-G. Two-layer Lanczos iteration approach to molecular spectroscopic calculation. J. Chem. Phys.
**2002**, 117, 8190–8196. [Google Scholar] [CrossRef] - Yu, H.-G. Full-dimensional quantum calculations of vibrational spectra of six-atom molecules. I. Theory and numerical results. J. Chem. Phys.
**2004**, 120, 2270–2284. [Google Scholar] [CrossRef] [PubMed] - Yu, H.-G. Converged quantum dynamics calculations of vibrational energies of CH
_{4}and CH_{3}D using an ab initio potential. J. Chem. Phys.**2004**, 121, 6334. [Google Scholar] [CrossRef] [PubMed] - Yu, H.-G. A rigorous full-dimensional quantum dynamics calculation of the vibrational energies of H3O-2. J. Chem. Phys.
**2006**, 125, 204306. [Google Scholar] [CrossRef] [PubMed] - Yurchenko, S.N.; Thiel, W.; Jensen, P. Theoretical ROVibrational Energies (TROVE): A robust numerical approach to the calculation of rovibrational energies for polyatomic molecules. J. Mol. Spectrosc.
**2007**, 245, 126–140. [Google Scholar] [CrossRef] - Dawes, R.; Carrington, T. How to choose 1-D basis functions so that a very efficient multidimensional basis may be extracted from a direct product of the 1-D functions: Energy levels of coupled systems with as many as 16 coordinates. J. Chem. Phys.
**2005**, 122, 134101:1–134101:14. [Google Scholar] [CrossRef] [PubMed] - Colbert, D.T.; Miller, W.H. A novel discrete variable representation for quantum mechanical reactive scattering via the S-matrix Kohn method. J. Chem. Phys.
**1992**, 96, 1982. [Google Scholar] [CrossRef] - Shimshovitz, A.; Tannor, D.J. Phase-Space Approach to Solving the Time-Independent Schrodinger Equation. Phys. Rev. Lett.
**2012**, 109, 070402. [Google Scholar] [CrossRef] [PubMed] - Carter, S.; Handy, N.C. The variational method for the calculation of ro-vibrational energy levels. Comput. Phys. Rep.
**1986**, 5, 117–171. [Google Scholar] [CrossRef] - Halonen, L.; Noid, D.W.; Child, M.S. Local mode predictions for excited stretching vibrational states of HCCD and H
^{12}C^{13}CH. J. Chem. Phys.**1983**, 78, 2803. [Google Scholar] [CrossRef] - Halonen, L.; Child, M.S. Local mode theory for C
_{3v}molecules: CH_{3}D, CHD_{3}, SiH_{3}D, and SiHD_{3}. J. Chem. Phys.**1983**, 79, 4355. [Google Scholar] [CrossRef] - Maynard, A.; Wyatt, R.E.; Iung, C. A quantum dynamical study of CH overtones in fluoroform. II. Eigenstate analysis of the v
_{CH}=1 and v_{CH}=2 regions. J. Chem. Phys.**1997**, 106, 9483. [Google Scholar] [CrossRef] - Maynard, A.T.; Wyatt, R.E.; Iung, C. A quantum dynamical study of CH overtones in fluoroform. I. A nine-dimensional ab initio surface, vibrational spectra and dynamics. J. Chem. Phys.
**1995**, 103, 8372. [Google Scholar] [CrossRef] - Iung, C.; Leforestier, C.; Wyatt, R.E. Wave operator and artificial intelligence contraction algorithms in quantum dynamics: Application to CD
_{3}H and C_{6}H_{6}. J. Chem. Phys.**1993**, 98, 6722. [Google Scholar] [CrossRef] - Poirier, B. Using wavelets to extend quantum dynamics calculations to ten or more degrees of freedom. J. Theor. Comput. Chem.
**2003**, 2, 65. [Google Scholar] [CrossRef] - Poirier, B.; Salam, A. Quantum dynamics calculations using symmetrized, orthogonal Weyl-Heisenberg wavelets with a phase space truncation scheme. III. Representations and calculations. J. Chem. Phys.
**2004**, 121, 1704–1724. [Google Scholar] [CrossRef] [PubMed] - Poirier, B.; Salam, A. Quantum dynamics calculations using symmetrized, orthogonal Weyl-Heisenberg wavelets with a phase space truncation scheme. II. Construction and optimization. J. Chem. Phys.
**2004**, 121, 1690–1703. [Google Scholar] [CrossRef] [PubMed] - Wang, X.-G.; Carrington, T. The utility of constraining basis function indices when using the lanczos algorithm to calculate vibrational energy levels. J. Phys. Chem. A
**2001**, 105, 2575–2581. [Google Scholar] [CrossRef] - Halverson, T.; Poirier, B. One Million Quantum States of Benzene. J. Phys. Chem. A
**2015**, 119, 12417–12433. [Google Scholar] [CrossRef] [PubMed] - Brown, J.; Carrington, T. Using an expanding nondirect product harmonic basis with an iterative eigensolver to compute vibrational energy levels with as many as seven atoms. J. Chem. Phys.
**2016**, 145, 144104:1–144104:10. [Google Scholar] [CrossRef] [PubMed] - Brown, J.; Carrington, T. Assessing the utility of phase-space-localized basis functions: Exploiting direct product structure and a new basis function selection procedure. J. Chem. Phys.
**2016**, 144, 244115:1–244115:10. [Google Scholar] [CrossRef] [PubMed] - Avila, G.; Carrington, T. Solving the vibrational Schrödinger equation using bases pruned to include strongly coupled functions and compatible quadratures. J. Chem. Phys.
**2012**, 137, 174108. [Google Scholar] [CrossRef] [PubMed] - Avila, G.; Carrington, T., Jr. Pruned bases that are compatible with iterative eigensolvers and general potentials: New results for CH
_{3}CN. Chem. Phys.**2017**, 482, 3–8. [Google Scholar] [CrossRef] - Petras, K. Fast calculation of coefficients in the Smolyak algorithm. Numer. Algorithms
**2001**, 26, 93–109. [Google Scholar] [CrossRef] - Avila, G.; Carrington, T., Jr. Using nonproduct quadrature grids to solve the vibrational Schrödinger equation in 12D. J. Chem. Phys.
**2011**, 134, 054126. [Google Scholar] [CrossRef] [PubMed] - Leclerc, A.; Carrington, T. Calculating vibrational spectra with sum of product basis functions without storing full-dimensional vectors or matrices. J. Chem. Phys.
**2014**, 140, 174111:1–174111:13. [Google Scholar] [CrossRef] [PubMed] - Thomas, P.S.; Carrington, T., Jr. Using nested contractions and a hierarchical tensor format to compute vibrational spectra of molecules with seven atoms. J. Phys. Chem. A
**2015**, 119, 13074–13091. [Google Scholar] [CrossRef] [PubMed] - Thomas, P.S.; Carrington, T., Jr. An intertwined method for making low-rank, sum-of-product basis functions that makes it possible to compute vibrational spectra of molecules with more than 10 atoms. J. Chem. Phys.
**2017**, 146, 204110:1–204110:15. [Google Scholar] [CrossRef] [PubMed] - Zhang, T.; Golub, G.H. Rank-One Approximation to High Order Tensors. Matrix Anal. Appl.
**2001**, 23, 534–550. [Google Scholar] [CrossRef] - Beylkin, G.; Mohlenkamp, M.J. Numerical operator calculus in higher dimensions. Proc. Natl. Acad. Sci. USA
**2002**, 99, 10246–10251. [Google Scholar] [CrossRef] [PubMed] - Beylkin, G.; Mohlenkamp, M.J. Algorithms for Numerical Analysis in High Dimensions. Sci. Comput.
**2005**, 26, 2133–2159. [Google Scholar] [CrossRef] - Meyer, H.D.; Gatti, F.; Worth, G.A. (Eds.) Multidimensional Quantum Dynamics: MCTDH Theory and Applications; Wiley-VCH: Weinheim, Germany, 2009. [Google Scholar]
- Pelaez, D.; Meyer, H.-D. The multigrid POTFIT (MGPF) method: Grid representations of potentials for quantum dynamics of large systems. J. Chem. Phys.
**2013**, 138, 014108. [Google Scholar] [CrossRef] [PubMed] - Manzhos, S.; Carrington, T. Using neural networks to represent potential surfaces as sums of products. J. Chem. Phys.
**2006**, 125, 194105. [Google Scholar] [CrossRef] [PubMed] - Manzhos, S.; Carrington, T. Using redundant coordinates to represent potential energy surfaces with lower-dimensional functions. J. Chem. Phys.
**2007**, 127, 014103. [Google Scholar] [CrossRef] [PubMed] - Manzhos, S.; Carrington, T. Using neural networks, optimized coordinates, and high-dimensional model representations to obtain a vinyl bromide potential surface. J. Chem. Phys.
**2008**, 129, 224104. [Google Scholar] [CrossRef] [PubMed] - Lehoucq, R.B.; Sorensen, D.C.; Yang, C. ARPACK Users’ Guide: Solution of Large-Scale Eigenvalue Problems With Implicitly Restarted Arnoldi Methods SIAM, Philadelphia, PA, 1998. Available online: http://www.caam.rice.edu/software/ARPACK (accessed on 21 December 2017).
- Watson, J.K.G. Simplification of the molecular vibration-rotation hamiltonian. Mol. Phys.
**1968**, 15, 479–490. [Google Scholar] [CrossRef] - Jaeckle, A.; Meyer, H.-D. Product representation of potential energy surfaces. J. Chem. Phys.
**1996**, 104, 7974. [Google Scholar] [CrossRef] - Chen, H.; Light, J.C. Vibrations of the carbon dioxide dimer. J. Chem. Phys.
**2000**, 112, 5070–5080. [Google Scholar] [CrossRef] - Li, H.; Roy, P.-N.; Le Roy, R.J. Analytic Morse/long-range potential energy surfaces and predicted infrared spectra for CO
_{2}-H_{2}. J. Chem. Phys.**2010**, 132, 214309. [Google Scholar] [CrossRef] [PubMed] - Thomas, P.S.; Carrington, T. Unpublished work. J. Chem. Phys.
**2018**.

© 2018 by the author. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).