Eckart ro-vibrational Hamiltonians via the gateway Hamilton operator: theory and practice

TL;DR

This paper discusses the development and application of Eckart ro-vibrational Hamiltonians using the gateway Hamiltonian method, demonstrating that rotation is unnecessary for constructing these Hamiltonians and highlighting the advantages of the gateway approach.

Contribution

It introduces a new application of the gateway Hamiltonian method to Eckart Hamiltonians, simplifying calculations by eliminating the need for rotation and providing a geometric interpretation of Eckart displacements.

Findings

Numerical examples confirm no need for rotation in Eckart Hamiltonians.

Gateway method simplifies calculations by avoiding differentiation of the rotating matrix.

A shortest Eckart displacement vector can be uniquely defined and is often shorter.

Abstract

Recently, a general expression for Eckart-frame Hamilton operators has been obtained by the gateway Hamiltonian method ({\it J. Chem. Phys.} {\bf 142}, 174107 (2015); {\it ibid.} {\bf 143}, 064104 (2015)). The kinetic energy operator in this general Hamiltonian is nearly identical with that of the Eckart-Watson operator even when curvilinear vibrational coordinates are employed. Its different realizations correspond to different methods of calculating Eckart displacements. There are at least two different methods for calculating such displacements: rotation and projection. In this communication the application of Eckart Hamiltonian operators constructed by rotation and projection, respectively, is numerically demonstrated in calculating vibrational energy levels. The numerical examples confirm that there is no need for rotation to construct an Eckart ro-vibrational Hamiltonian. The application of the gateway method is advantageous even when rotation is used, since it obviates the need for differentiation of the matrix rotating into the Eckart frame. Simple geometrical arguments explain that there are infinitely many different methods for calculating Eckart displacements. The geometrical picture also suggests that a unique Eckart displacement vector may be defined as the shortest (mass-weighted) Eckart displacement vector among Eckart displacement vectors corresponding to configurations related by rotation. Its length, as shown analytically and demonstrated by way of numerical examples, is equal to or less than that of the Eckart displacement vector one can obtain by rotation to the Eckart frame.

Tables

Table 1. Table 1: The rotational KEO

$2{\hat{T}}_{\text{rot}}$	${\hat{J}}_{\alpha }{\mu }_{\alpha \beta }{\hat{J}}_{\beta }$
${\mu }_{\alpha \beta }$	${\left({𝑰}^{{}_{}{}^{\prime }-1}\right)}_{\alpha \beta }$
$I_{\alpha \beta }^{{}^{\prime }}$	$I_{\alpha \beta }-\sum _{jkl}{\zeta }_{jl}^{\alpha }{\zeta }_{kl}^{\beta }{c}_j{c}_k=I_{\alpha \gamma }^{{}^{\prime \prime }}{\left({\left[{𝑰}^0\right]}^{-1}\right)}_{\gamma \delta }{I}_{\delta \beta }^{{}^{\prime \prime }}$
${\zeta }_{jl}^{\alpha }$	${ϵ}_{\alpha \beta \gamma }\sum _n{m}_n{d}_{\beta n}^{\text{E},j}{d}_{\gamma n}^{\text{E},l}$
$I_{\alpha \beta }^{{}^{\prime \prime }}$	$I_{\alpha \beta }^0+\frac{1}{2}\sum _k{a}_k^{\alpha \beta }{c}_k$
$a_k^{\alpha \beta }$	$2{ϵ}_{\alpha \gamma ϵ}{ϵ}_{\beta \delta ϵ}\sum _n{m}_n{a}_{\gamma n}^0{d}_{\delta n}^{\text{E},k}$

Table 2. Table 2: The vibrational KEO

$2{\hat{T}}_{\text{vib}}$	${\varrho }^2\sum _r{\hat{p}}_r{𝒢}_{rs}{\hat{p}}_s$
${𝒢}_{rs}$	$G_{rs}+{\varrho }^2{C}_{\alpha r}{\mu }_{\alpha \beta }{C}_{\beta s}$
$G_{rs}$	$\sum _i\frac{\partial s_r}{\partial c_i}\frac{\partial s_s}{\partial c_i}$
$C_{\alpha r}$	$\sum _{ij}{\zeta }_{ij}^{\alpha }{c}_i\frac{\partial s_r}{\partial c_j}$

Table 3. Table 3: The rotational-vibrational coupling KEO

$2{\hat{T}}_{\text{rot-vib}}$	$-{\varrho }^2\sum _r\left({\hat{J}}_{\alpha }{𝒞}_{\alpha r}{\hat{p}}_r+{\hat{p}}_r{𝒞}_{r\alpha }^T{\hat{J}}_{\alpha }\right)$
${𝒞}_{\alpha r}$	${\mu }_{\alpha \beta }{C}_{\beta r}$

Table 4. Table 4: The potential and the pseudo-potential energy

$V$	$V\left(\left\{a_{\beta n}^0+\sum _j{c}_j{d}_{\beta n}^{\text{E},j}\right\}\right)$
$V_{\text{ps}}$	$\tilde{V}-\frac{{\hslash }^2}{8}{\mu }_{\alpha \alpha }$
$\frac{8}{{\hslash }^2}\tilde{V}$	$-{\sum }_{ij}\left[\frac{\partial {\widetilde{𝒢}}_{ij}}{\partial c_i}\right]\frac{\partial \ln det𝑮}{\partial c_j}$ $+{\sum }_{ij}{\widetilde{𝒢}}_{ij}\left\{-\frac{{\partial }^2\ln det𝑮}{\partial c_i\partial c_j}+\frac{1}{4}\frac{\partial \ln det𝑮}{\partial c_i}\frac{\partial \ln det𝑮}{\partial c_j}\right\}$
	$-{\sum }_{ij}{\sum }_r\left[\frac{\partial }{\partial s_r}\frac{\partial s_r}{\partial c_i}\right]{\widetilde{𝒢}}_{ij}\frac{\partial \ln \left(det𝝁det𝑮\right)}{\partial c_j}$ $-{\sum }_{ijk}{\sum }_s\left[\frac{\partial }{\partial c_i}\frac{\partial s_s}{\partial c_j}\right]\frac{\partial c_k}{\partial s_s}\frac{\partial \ln \left(det𝝁det𝑮\right)}{\partial c_k}$
	$+\frac{1}{4}{\sum }_{ij}{\widetilde{𝒢}}_{ij}\left\{\frac{\partial \ln det𝝁}{\partial c_i}\frac{\partial \ln det𝑮}{\partial c_j}+\frac{\partial \ln det𝑮}{\partial c_i}\frac{\partial \ln det𝝁}{\partial c_j}\right\}$
${\widetilde{𝒢}}_{ij}$	${\delta }_{ij}+{\sum }_{kl}{c}_k{\zeta }_{ki}^{\alpha }{\mu }_{\alpha \beta }{\zeta }_{lj}^{\beta }{c}_l$

Table 5. Table 5: Numbering of the atoms, the atomic masses, and the coordinates b α ~ ​ n subscript 𝑏 ~ 𝛼 𝑛 b_{\tilde{\alpha}n} of the atoms in the initial system of axes.

Atom	Numbering of the atoms; $n$	Masses	X	Y	Z
H1	$1$	$1.00782522$	$r_1\cos \varphi$	$r_1\sin \varphi$	$0$
O	$2$	$15.99491502$	$0$	$0$	$0$
H2	$3$	$1.00782522$	$r_2$	$0$	$0$

Table 6. Table 6: Basis vectors for the vibrational space

$\alpha n$	${𝒅}^{\text{E},1}$	${𝒅}^{\text{E},2}$	${𝒅}^{\text{E},3}$
$x1$	-0.2409179258	-0.6669945045	0.0308844916
$y1$	0.9356174659	-0.1315179013	-0.0084552498
$z1$	0.0000000000	-0.0000000000	-0.0000000000
$x2$	0.0151800220	0.0420267243	0.0589871470
$y2$	-0.0589524156	0.0510523585	-0.0012189233
$z2$	0.0000000000	-0.0000000000	0.0000000000
$x3$	0.0000000000	0.0000000000	-0.9670531688
$y3$	0.0000000000	-0.6787199452	0.0278004437
$z3$	0.0000000000	0.0000000000	-0.0000000000

Table 7. Table 7: Derivatives of Eckart coordinates with respect to ϕ italic-ϕ \phi at ( r 1 = 1.358430 ​ Å , r 2 = 0.658430 ​ Å , ϕ = 1.322818 ​ rad ) formulae-sequence subscript 𝑟 1 1.358430 Å formulae-sequence subscript 𝑟 2 0.658430 Å italic-ϕ 1.322818 rad \left(r_{1}=1.358430\text{\AA},~{}r_{2}=0.658430\text{\AA},~{}\phi=1.322818\textrm{rad}\right) .

$h$	Approximate $\partial {𝒂}^{\text{E}}/\partial \varphi$	Exact $\partial {𝒂}^{\text{E}}/\partial \varphi$	Approx$-$Exact
0.1000000000	-0.4054459586	-0.4057377040	0.0002917454
0.1000000000	-0.0207341020	-0.0206473673	-0.0000867347
0.1000000000	0.0000000000	0.0000000000	0.0000000000
0.1000000000	0.0145783041	0.0145969975	-0.0000186934
0.1000000000	0.0263719843	0.0263829584	-0.0000109741
0.1000000000	0.0000000000	0.0000000000	0.0000000000
0.1000000000	0.1740777305	0.1740727981	0.0000049323
0.1000000000	-0.3978083597	-0.3980692607	0.0002609009
0.1000000000	0.0000000000	0.0000000000	0.0000000000
0.0000100000	-0.4057377040	-0.4057377040	0.0000000000
0.0000100000	-0.0206473675	-0.0206473673	-0.0000000002
0.0000100000	0.0000000000	0.0000000000	0.0000000000
0.0000100000	0.0145969975	0.0145969975	0.0000000000
0.0000100000	0.0263829584	0.0263829584	0.0000000000
0.0000100000	0.0000000000	0.0000000000	0.0000000000
0.0000100000	0.1740727982	0.1740727981	0.0000000000
0.0000100000	-0.3980692606	-0.3980692607	0.0000000001
0.0000100000	0.0000000000	0.0000000000	0.0000000000
0.0000000001	-0.4057600000	-0.4057377040	-0.0000222960
0.0000000001	-0.0206500000	-0.0206473673	-0.0000026327
0.0000000001	0.0000000000	0.0000000000	0.0000000000
0.0000000001	0.0145995000	0.0145969975	0.0000025025
0.0000000001	0.0263810000	0.0263829584	-0.0000019584
0.0000000001	0.0000000000	0.0000000000	0.0000000000
0.0000000001	0.1740700000	0.1740727981	-0.0000027981
0.0000000001	-0.3980600000	-0.3980692607	0.0000092607
0.0000000001	0.0000000000	0.0000000000	0.0000000000

Table 8. Table 8: Bending energy levels

( n1 n2 n3 )	BC	Gateway rot	Gateway proj	Gateway rot, no $V_{\text{ps}}$	Gateway proj, no $V_{\text{ps}}$
(0 0 0)	4630.3465	4630.295443	4630.295446	4650.23	4638.90
		$\Delta$	$\Delta$	$\Delta$	$\Delta$
(0 1 0)	1594.32	0.02	0.02	-0.78	-2.24
(0 2 0)	3152.01	0.04	0.03	-1.85	-4.58
(0 3 0)	4667.70	0.06	0.05	-3.4	-7.27
(0 4 0)	6134.11	0.78	0.07	-5.19	-10.85
(0 5 0)	7539.79	0.06	0.01	-10.18	-16.84
Basis size		$25\times 41\times 41$	$23\times 51\times 51$	$25\times 41\times 41$	$23\times 51\times 51$

Table 9. Table 9: Stretching energy levels

( n1 n2 n3 )	BC	Gateway rot	Gateway proj	Gateway rot, no $V_{\text{ps}}$	Gateway proj, no $V_{\text{ps}}$
		$\Delta$	$\Delta$	$\Delta$	$\Delta$
(1 0 0)	3656.49	0.04	0.04	0.43	0.21
(2 0 0)	7202.67	0.08	0.08	0.86	0.44
(3 0 0)	10602.76	0.1	0.1	1.3	0.7
(4 0 0)	13829.70	0.14	0.14	-0.56	0.91
(5 0 0)	16899.45	0.61	0.01	1.95	-0.22
(0 0 1)	3755.92	0.04	0.04	0.51	0.38
(0 0 2)	7444.93	0.08	0.08	-0.01	0.7
(0 0 3)	11034.09	0.12	0.12	1.47	1.02
(0 0 4)	14541.30	0.14	0.15	-2.57	-0.42
(0 0 5)	17954.91	0.18	0.18	2.36	1.53

Table 10. Table 10: Combinations

( n1 n2 n3 )	BC	Gateway rot	Gateway proj	Gateway rot, no $V_{\text{ps}}$	Gateway proj, no $V_{\text{ps}}$
		$\Delta$	$\Delta$	$\Delta$	$\Delta$
(1 1 0)	5234.29	0.06	0.06	-0.33	-1.97
(1 2 0)	6775.03	0.1	0.08	-1.33	-4.23
(1 3 0)	8273.24	0.4	0.07	-2.59	-6.91
(1 4 0)	9719.75	1.47	-0.4	-4.15	-11.54
(0 1 1)	5332.06	0.06	0.06	-0.19	-1.75
(0 1 2)	9002.14	0.1	0.1	0.35	-1.35
(0 1 4)	16057.58	0.17	0.18	1.38	-0.61
(0 1 5)	19449.25	0.2	0.19	-0.13	-0.28
(0 6 1)	12571.35	0.2	0.22	0.74	-0.88

Table 11. Table 11: Combinations

( n1 n2 n3 )	BC	Gateway rot	Gateway proj	Gateway rot, no $V_{\text{ps}}$	Gateway proj, no $V_{\text{ps}}$
		$\Delta$	$\Delta$	$\Delta$	$\Delta$
(1 1 1)	8809.59	0.09	0.09	0.23	-1.53
(1 2 1)	10332.40	0.12	0.11	-0.65	-3.61
(1 3 1)	11815.47	0.34	0.26	-1.76	-5.86
(2 1 1)	12156.52	0.13	0.12	0.67	-1.3
(3 1 1)	15355.27	0.16	0.15	1.04	-1.18
(1 1 2)	12408.42	-0.21	-0.4	-0.09	-1.09
(1 2 2)	13911.72	0.20	-0.47	-0.16	-3.27
(2 2 2)	17226.08	-0.08	0.79	-0.16	-3.7
(1 1 3)	15839.10	0.21	-0.24	1.22	-0.9

Table 12. Table 12: Analytical expressions of particular solutions, 𝒉 E , j superscript 𝒉 E 𝑗 \bm{h}^{\textrm{E},j} , of the Eckart conditions of an N 𝑁 N -atom non-linear molecule

(206)

Table 13. Table 13: Table  12 continued

(312)

Equations

$$u_{\overline{\alpha }n}=R_{\overline{\alpha }}+{\left[{\mathbf{S}}^{-1}\left(\theta ,\varphi ,\chi \right)\right]}_{\overline{\alpha }\beta }{a}_{\beta n}\left(\left\{s_r\right\}\right),$$

$$a_{\beta n}\left(\left\{s_r\right\}\right)=a_{\beta n}^0+d_{\beta n}\left(\left\{s_r\right\}\right),$$

$$\sum _n{m}_n{a}_{\beta n}^0=0,$$

$$\sum _n{m}_n{a}_{\beta n}\left(\left\{s_r\right\}\right)=0$$

$$\sum _n{m}_n{d}_{\beta n}\left(\left\{s_r\right\}\right)=0.$$

$$T$$

$$\sum _{n=1}^N{m}_n{\mathbf{a}}_n^0\times {\mathbf{d}}_n=0$$

$$T$$

$$\sum _{n=1}^N{m}_n{\mathbf{a}}_n=0,$$

$$\sum _{n=1}^N{m}_n{\mathbf{a}}_n^0\times {\mathbf{a}}_n=0,$$

$$\mathbf{E}\mathbf{d}=0,$$

$${\mathbf{d}}^T=\left({\mathbf{d}}_1^T,{\mathbf{d}}_2^T,\dots ,{\mathbf{d}}_N^T\right),$$

$${\mathbf{d}}^{\text{E}}=\sum _{j=1}^K{b}_j{\mathbf{h}}^{\text{E},j},$$

$$\displaystyle\textrm{diag}(\bm{m})=\left(\begin{array}[]{cccccccccc}m_{1}&m_{1}&m_{1}&m_{2}&m_{2}&m_{2}&\ldots&m_{N}&m_{N}&m_{N}\end{array}\right).$$

$${\left[{\mathbf{d}}^{\text{E},j}\right]}^T\mathbf{m}{\mathbf{d}}^{\text{E},i}={\delta }_{ji},$$

$${\mathbf{d}}_n^{\text{E}}=\sum _{j=1}^{3N-6}{c}_j{\mathbf{d}}_n^{\text{E},j}$$

$$u_{\overline{\alpha }n}=R_{\overline{\alpha }}+{\left[{\mathbf{S}}^{-1}\left(\theta ,\varphi ,\chi \right)\right]}_{\overline{\alpha }\beta }{\mathbf{a}}_{\beta n}^{\text{E}},$$

$${\mathbf{a}}_{\beta n}^{\text{E}}={\mathbf{a}}_{\beta n}^0+{\mathbf{d}}_{\beta n}^{\text{E}},$$

$${\mathbf{d}}_{\beta n}^{\text{E}}=\sum _{j=1}^{3N-6}{c}_j{d}_{\beta n}^{\text{E},j},$$

$$c_j={\left[{\mathbf{m}}^{1/2}{\mathbf{d}}^{\text{E},j}\right]}^T{\mathbf{m}}^{1/2}{\mathbf{d}}^{\text{E}}.$$

$$\mathbf{P}=\sum _j{\mathbf{m}}^{1/2}{\mathbf{d}}^{\text{E},j}{\left({\mathbf{m}}^{1/2}{\mathbf{d}}^{\text{E},j}\right)}^T.$$

$${\mathbf{m}}^{1/2}{\mathbf{d}}^{\text{E}}\left(\left\{s_r\right\}\right)=\mathbf{P}{\mathbf{m}}^{1/2}\mathbf{d}\left(\left\{s_r\right\}\right).$$

$$c_j\left(\left\{s_r\right\}\right)$$

$$c_j^0={\left[{\mathbf{m}}^{1/2}{\mathbf{d}}^{\text{E},j}\right]}^T{\mathbf{m}}^{1/2}{\mathbf{a}}^0.$$

$$\left[\frac{1}{2}\sum _{rs}{\hat{p}}_r{\mathcal{G}}_{rs}{\hat{p}}_s+V+V_{\text{ps}}\right]|\Psi ⟩=\lambda |\Psi ⟩,$$

$$\frac{\partial s_r}{\partial c_j}=\frac{\partial s_r\left(\mathbf{a}\right)}{\partial c_j}=\frac{\partial s_r\left({\mathbf{a}}^{\text{E}}\right)}{\partial c_j}=\frac{\partial s_r}{\partial {\mathbf{a}}^{\text{E}}}\frac{\partial {\mathbf{a}}^{\text{E}}}{\partial c_j}$$

$$={\frac{\partial s_r}{\partial \mathbf{a}}|}_{\mathbf{a}={\mathbf{a}}^{\text{E}}}{\mathbf{d}}^{\text{E},j}={\mathbf{s}}_r({\mathbf{a}}^{\text{E}})\cdot {\mathbf{d}}^{\text{E},j},$$

$$\frac{\partial c_j}{\partial s_r}={\left[{\mathbf{d}}^{\text{E},j}\right]}^T\mathbf{m}\frac{\partial \mathbf{a}\left(\left\{s_r\right\}\right)}{\partial s_r}.$$

$${\mathbf{a}}^{\text{E}}={\mathbf{a}}^0+\sum _{j=1}^{3N-6}{c}_j{\mathbf{d}}^{\text{E},j},$$

$$\frac{\partial {\mathbf{a}}^{\text{E}}}{\partial s_r}=\sum _{j=1}^{3N-6}\frac{\partial c_j}{\partial s_r}{\mathbf{d}}^{\text{E},j}.$$

$${\left(\frac{\partial {\mathbf{a}}_n^{\text{E}}}{\partial s_r}\right)}_{\mathbf{s}={\mathbf{s}}^0}$$

$$V_{\text{ps}}=-\frac{{\hslash }^2}{2}\left(\frac{\cos \varphi }{r_1{r}_2}+\frac{1}{4}{G}_{\varphi \varphi }\left(1+{\csc }^2\varphi \right)\right)$$

Full text

I Introduction

The importance of Eckart conditions and Eckart ro-vibrational Hamiltonians cannot be overestimated. Their importance is described succinctly by Sutcliffesutcliffe-hinchliffe

”It would, I believe, be widely agreed that the modern theory of molecular spectra began with publication by Carl Eckart in 1935 of his paper Some Studies Concerning Rotating Axes and Polyatomic Moleculeseckart . It would be probably also be widely agreed that the apogee of this work occurred in 1968 when James K. G. Watsonwatson , published Simplification of the molecular vibration-rotation hamiltonian which put Eckart’s classical mechanical form into a proper quantum mechanical one. This leads to the wave mechanical problem for molecular vibration-rotational motion specified by what we shall call the Eckart-Watson Hamiltonian.”

Nearly at the same time of Watson’s work, the exact, general form of the quantum mechanical rotational-vibrational kinetic energy operator (KEO) given in curvilinear vibrational coordinates has been established by Meyer and Günthardmeyer-gunthard , and Pickettpickett . However, accounting for Eckart’s rotating axes in this general KEO has turned out to be difficult and has remained the subject of ongoing research wei-carrington-1 ; wei-carrington-2 ; mardis-sibert ; wei-carrington-3 ; rey ; wei ; meremianin-1 ; pesonen ; mccoy-burleigh-sibert ; wang-sibert ; szidarovszky-fabri-csaszar ; wang-carrington ; fabri-matyus-csaszar ; sadri ; yachmenev-yurchenko ; lauvergnat-luis . Research on how to employ the general ro-vibrational KEO possibly by accounting for the Eckart conditions has gained practical motivation with the advance of experimental spectroscopic studies of highly excited vibrational sates and large-amplitude internal motions of molecules h5+ ; ch5+ ; ch4h2o-1 ; ch4h2o-2 .

If it is so difficult to use, one may wonder, why to use Eckart conditions at all. One reason is that they minimize Coriolis coupling, thereby they give justification of thinking about the nuclear motions of a molecule as rotational and vibrational motions. With the Coriolis terms minimized one might also expect that fewer basis functions may suffice to obtain converged ro-vibrational energy levels in a variational calculationwang-carrington ; yachmenev-yurchenko . The use of Eckart frame may be the best choice in evaluating dipole-moment matrix elementssueur . Furthermore, applications of the Eckart frame are not restricted to small and medium size molecules, but it has been employed in describing nuclear motion dynamics of biomolecules, e.g. proteins noguti-gon .

The current situation on using Eckart frame is well described by Wang and Carringtonwang-carrington

”…, it appears that there is a massive gulf between knowing that it would be best to use an Eckart frame and developing a theoretical/computational scheme for exploiting the Eckart advantage while at the same time using curvilinear coordinates that enable one to cope with large amplitude motion. When normal coordinates are used it is straightforward to use an Eckart frame, however, the Eckart KEO in internal coordinates is very complicated. Such KEOs have been derived for triatomic molecules, but never employed to compute spectra. This conundrum is resolved by computing G-matrix elements numerically.”

Although there have been interesting theoretical developments, such as the application of geometric algebra to derive Eckart KEOs pesonen , and a method of analytical differentiation of the rotation matrix transforming into the Eckart frame has been introducedwang-song ; changala ; changala-baraban , practical application of the Eckart frame has remained brute force numerical work considerably complicated by employing an energy operator with little resemblance of the Eckart-Watson Hamiltonian. All approaches to incorporate the Eckart conditions into the nuclear motion Hamiltonian have assumed explicitlywang-carrington ; fabri-matyus-csaszar ; sadri ; yachmenev-yurchenko ; lauvergnat-luis or implicitlypesonen a rotation matrix determined such that the Eckart-axis conditionsstrauss ; redding ; meyer ; redding2 ; kudin ; dierksen ; dymarsky ; galbraith ; krasnos , which should not be mistaken for the (rotational) Eckart conditions, be satisfied. However, studies on the gateway Hamiltonian methodszalay-eckart ; szalay-gateway ; szalay-aspects have shown another solution to this question: Projection.

The purpose of the present work is to simplify and advance the use of Eckart conditions in variational calculations of ro-vibrational energy levels by employing the gateway Hamiltonian method.

The KEO in the general Eckart Hamilton operator as obtained by the gateway method is of an expression nearly identical with that of the Eckart-Watson operator even when curvilinear vibrational coordinates are employed. Its use, as will be numerically demonstrated,

• •

obviates the need for differentiation of the matrix rotating into the Eckart frame,

• •

allows construction of an Eckart ro-vibrational KEO without rotation,

• •

and leads to defining optimal Eckart displacements.

The outline of the paper is as follows. In Section II the basic ideas and formulas of the gateway Hamiltonian method are reviewed. A method of solving the vibrational Schrödinger equation corresponding to the general Eckart Hamilton operator is described in Section III. Numerical examples are given in Section IV to demonstrate that: a) There is no need for differentiating the matrix rotating into the Eckart frame. b) The rotation and projection methods of calculating Eckart displacements lead to different Hamiltonians. c) Although they are different, these Hamilton operators have the same spectrum. By considering the geometry associated with the Eckart conditions it is explained why there are infinitely many different Eckart displacements (and Eckart Hamiltonians) and how optimal Eckart displacements may be defined (Section V). Numerical examples of calculating optimal Eckart displacements are also given. Section VI summarizes the results. Two important questions, but mainly mathematical in character, are detailed in Appendices. In Appendix A analytical solution of the Eckart conditions is presented, while Appendix B describes the method employed to calculate optimal Eckart displacements.

II The general-form Eckart ro-vibrational Hamilton operator

The derivation of the nuclear motion Hamilton operator of a molecule of $N$ atoms starts by relating the Cartesian laboratory system (LS) coordinates, $u_{\overline{\alpha}n}$ $\left(\overline{\alpha}=X,Y,Z;~{}~{}n=1,2,\ldots,N\right)$, of the atoms of the molecule to the translational, $R_{\overline{\alpha}}$, rotational, $\theta,\phi,\chi$, and vibrational, $s_{r},r=1,2,\ldots,3N-6$, coordinates of the molecule (assumed to be non-linear) by the equation

[TABLE]

where $\beta=x,y,z$ label the axes of a molecule fixed Cartesian system (MS), $\bm{S}\left(\theta,\phi,\chi\right)$ is a real rotation matrix parametrized by Euler angles,

[TABLE]

[TABLE]

and the displacements $d_{\beta n}\left(\left\{s_{r}\right\}\right)$ are as yet unknown functions of the coordinates $s_{r}$. Summation over repeated Greek indices is assumed here and throughout the paper. $a_{\beta n}^{0}$ are the coordinates of the atoms of a reference configuration. $m_{n}$ denotes the mass of the $n$th atom of the molecule. The origin of the MS is fixed to the center of mass of the reference configuration and its axes can be oriented and fixed in any convenient way to the reference configuration. Thus the MS has been completely determined and the atom coordinates, $a_{\beta n}\left(\left\{s_{r}\right\}\right)$, of any distorted configuration must be given with respect to this MS.

We aim at determining the $3N$ unknown functions $d_{\beta n}\left(\left\{s_{r}\right\}\right)$ ( or $a_{\beta n}\left(\left\{s_{r}\right\}\right)$ ).

This goal will be achieved in two steps: 1) A general form of displacements obeying the Eckart conditions will be established. 2) Then it is shown how dependence of these displacements on vibrational coordinates can be introduced.

II.1 The Eckart conditions

Together with the coordinates of overall translation and rotation there are $3N+6$ unknowns, but there are only $3N$ equations in Eq.  (1). Since one must be able to calculate the coordinates of nuclear motion from given LS coordinates and vica versa, we must complete Eq.  (1) with six additional equations. Three equations are obtained immediately by noting that Eqs.  (1) imply that $a_{\beta n}\left(\left\{s_{r}\right\}\right)$ are atom coordinates in a translation reduced configuration space, that is

[TABLE]

must hold. Then by considering Eqs. (2) and (3) one finds that the displacements must obey three equations

[TABLE]

These equation are called translational Eckart conditions. Three conditions are still missing. They determine the functions $a_{\beta n}\left(\left\{s_{r}\right\}\right)$ (and $d_{\beta n}\left(\left\{s_{r}\right\}\right)$ ), that is the coordinates of distorted configurations in the MS.

The simplest conditions are geometrical and somewhat ad hoc. For instance, having fixed a right-handed MS to the reference configuration such that the $x$ axis is parallel to the bond connecting the atoms $A$ and $B$ and the $y$ axis lies in the plane containing the atoms $A,~{}B$, and $C$ and points to the half-plane containing $C$, we may require these conditions to hold even for distorted configurations.

To obtain the missing equations in more reasonable way let us consider the expression of the classical kinetic energy szalay-eckart :

[TABLE]

where $\bm{\omega}$ is an angular velocity vector whose components are defined as elements of the skew symmetric matrix ${\bm{S}}\frac{d\bm{S}^{-1}}{dt}$. $\bm{d}_{n}$ is defined as a column vector, $\bm{d}_{n}=\left(d_{xn},d_{yn},d_{zn}\right)^{T}$, with superscript $T$ denoting transposition. Its dependence on the vibrational coordinates is, for simplicity, not indicated explicitly. Finally, $\bm{a}_{n}^{0}=\left(a_{xn}^{0},a_{yn}^{0},a_{zn}^{0}\right)^{T}$. The $\cdot$ means dot product. Note that by requiring satisfaction of the translational Eckart conditions and the three equations

[TABLE]

the expression of kinetic energy simplifies to

[TABLE]

That is the translational motion becomes completely separated from rotation and vibrations, whereas the rotational and vibrational motions are decoupled at the reference configuration (at zero displacements). Eqs. (7) are called the rotational Eckart conditions. The translational and rotational Eckart conditions together are referred to as Eckart conditions.

Displacements obeying the Eckart conditions will be called Eckart displacements and denoted as $\bm{d}_{n}^{\textrm{E}}$. Atom coordinates $\bm{a}_{n}=\left(a_{xn},a_{yn},a_{zn}\right)^{T}$ obeying equations similar to the Eckart conditions,

[TABLE]

will be called Eckart coordinates and denoted as $\bm{a}_{n}^{\textrm{E}}$.

Clearly, to derive an Eckart KEO one must calculate Eckart displacements (or Eckart coordinates).

II.2 Eckart displacements and Eckart coordinates

It should be clear from Subsection  II.1 that it is the Eckart conditions and not some of their derivatives, the Eckart-axis conditions for example, which enter directly into the derivation of the nuclear motion energy operator. Therefore, we must look for solutions of the Eckart conditions. By definition $\bm{d}_{n}^{\textrm{E}}$ are such solutions. In other words, they solve the system of equations, a homogeneous system of six linear equations with $3N$ unknowns, provided by the Eckart conditions:

[TABLE]

where

[TABLE]

and matrix $\bm{E}$ can be read out from the Eckart conditions.

The general solution of such a system of equationscarl-meyer is

[TABLE]

where $K=3N-6$ (since $\textrm{rank}(\bm{E})=6$), $b_{j}$ are free variables, and $\bm{h}^{\textrm{E},j}$ are particular solutions (that is $\bm{h}^{\textrm{E},j}$ obey the Eckart conditions).

In Appendix A analytical expressions of $\bm{h}^{\textrm{E},j}$ applicable to any molecule are derived. They depend only on the atomic masses and the atom coordinates of the reference configuration. Though it is less general than the method described in Appendix A, it is interesting to note that particular solutions of the Eckart conditions can be also obtained analytically by choosing $3N-6~{}(3N-5)$ internal coordinates and calculating the corresponding Wilson $\bm{s}_{j}$ vectorswilson . Due to the translational and rotational invariance of the $\bm{s}_{j}$ vectors calculated at the reference configuration the vectors $\bm{h}^{\textrm{E},j}=\bm{m}^{-1}\bm{s}_{j}$ are particular solutions, where $\bm{m}$ denotes a $3N$ by $3N$ diagonal matrix with

[TABLE]

The vectors $\bm{h}^{\textrm{E},j}$ span a subspace of the configuration space. It is called vibrational space. The vectors $\bm{h}^{\textrm{E},j}$ are, in general, not orthogonal. To simplify calculations we shall use in place of $\bm{h}^{\textrm{E},j}$ a set of vectors $\bm{d}^{\textrm{E},j}$ obeying the normalization condition

[TABLE]

and write the general solution of the Eckart conditions as

[TABLE]

with $c_{j}$ denoting free parameters.

Mathematically the $c_{j}$ are just free parameters. Physically, they are vibrational coordinates. The relations of the LS coordinates and nuclear motion coordinates assuming the Eckart conditions read as

[TABLE]

and now, one can derive the nuclear motion Hamiltonian corresponding to the Eckart conditions. The derivation gives the gateway Hamilton operator szalay-gateway . However, the question arises if one can use vibrational coordinates other than the $c_{j}$s.

II.3 Vibrational coordinate dependence of Eckart displacements and coordinates

The $c_{j}$ coefficients have geometrical significance. They are the coordinates of mass-weighted Eckart displacements, $\bm{m}^{1/2}\bm{d}^{\textrm{E}}$, in a coordinate system of $3N-6$ orthogonal axes in the space of $3N$-dimensional vectors satisfying the Eckart conditions. That is

[TABLE]

Yet, in other words, they are the coordinates of the projection of a mass-weighted displacement vector onto the vibrational space. The projection matrix is

[TABLE]

Therefore, the $c_{j}$ can be expressed in terms of geometrically defined coordinates.

In fact, $c_{j}$ can be expressed in terms of any sets of internal coordinates, $\left\{s_{r}\right\}$, as follows:

• •

Express the atom coordinates, $b_{\tilde{\alpha}n}$, where $\tilde{\alpha}=\texttt{X,Y,Z}$ label coordinates in an initial Cartesian system of axes, as functions of the chosen set of internal coordinates, i.e. $b_{\tilde{\alpha}n}=b_{\tilde{\alpha}n}\left(\left\{s_{r}\right\}\right)$. How this may be done has been studied in a number of papers thompson ; hilderbrandt ; lopata-kiss ; trigub ; essen . Then translate to the center of mass to obtain the coordinates $\bm{a}\left(\left\{s_{r}\right\}\right)$.

• •

Form the displacement vector $\bm{d}\left(\left\{s_{r}\right\}\right)=\bm{a}\left(\left\{s_{r}\right\}\right)-\bm{a}^{0}$.

• •

Project it to obtain the Eckart displacement vector $\bm{m}^{1/2}\bm{d}^{\textrm{E}}\left(\left\{s_{r}\right\}\right)$:

[TABLE]

• •

Then it follows, that

[TABLE]

where

[TABLE]

II.4 The general-form Eckart ro-vibrational Hamilton operator

The gateway Hamilton operator is an Eckart ro-vibrational Hamiltonian with exact KEO given in terms of the vibrational coordinates $c_{j}$. It is called gateway, since, due to the simple relation of the coordinates $c_{j}$ with other sets of vibrational coordinates, it can be easily transformed into Eckart Hamiltonians given in terms of other vibrational coordinates. The various terms in the general-form Eckart Hamilton operator given in terms of a general set of vibrational coordinates $s_{r}$, as obtained by transforming the gateway Hamiltonian, are summarized in Tables  1, 2, 3, and 4.

By examining the expressions one can see that no derivatives have to be calculated to obtain the rotational matrix $\bm{\mu}$, but one has to calculate the derivatives of the coordinates $s$ with respect to the coordinates $c$ to evaluate the vibrational matrix $\mathcal{G}$. One must also note that the vibrational matrix has contribution from the rotational matrix $\bm{\mu}$ and the matrix $\bm{C}$ appearing in the terms coupling rotational and vibrational momenta. Therefore, to check the validity and demonstrate the usefulness of the gateway method it suffices to consider the vibrational part of the general-form Eckart Hamiltonian.

In passing, it should be noted that the term $\tilde{V}$ was left out from the expression of the pseudo-potential in Ref.  szalay-gateway, . It was assumed that the pseudo-potential transformed as an ordinary scalar function under change of coordinates. But this is not a valid assumption littlejohn-reinsch . $\tilde{V}$ was also omitted in Ref.  szalay-aspects, , since it was assumed that Eq. (71) of Ref.  szalay-aspects, held. This equation does not hold in general, however.

III Solution of the vibrational Schrödinger equation

Numerical implementation of constructing the gateway vibrational Hamilton operator and a method of solving the corresponding Schrödinger equation is described below.

The vibrational Schrödinger equation reads as

[TABLE]

where $\hat{p}_{r}=-i\hbar\partial/\partial s_{r}$. When solving this equation we do not aim at extreme accuracy nor we want to obtain all vibrational energy levels. Our goal is just to show that the gateway method works fine. Then, if required, one can improve the results by employing more sophisticated methods to solve the eigenvalue equation. Therefore, we choose as simple a method of solution as possible:

• •

Direct product sinc discrete variable representation (DVR)colbert-miller ; szalay-smith ; littlejohn-cargo ; jerke , a special case of discretized continuous contracted Hermite distributed approximating functions szalay-smith , is employed.

• •

Approximate eigenvalues are obtained by Lanczos iterationlanczos ; cullum .

• •

The calculation of a Hamiltonian matrix vector product is carried out by partial summation bramley-carrington . It scales as $n^{3N-5}$ where $n$ is the number of basis functions applied to a single vibrational mode. (The Hamiltonian matrix is never calculated.)

The calculations consist of two main steps:

• •

The elements of the kinetic energy matrix and the potential energy are calculated on a grid.

• •

The Lanczos algorithm is employed to obtain approximate energy levels.

The first step is split. At first quantities independent of the grid points shown in Figure 1 are calculated.

Then the calculation continues by calculating the vibrational kinetic matrix and the potential and pseudo-potential values at the grid points as shown in Figure 2.

The difficult looking steps are the calculations of the derivatives $\frac{\partial s_{r}}{\partial c{j}}$ and the pseudo-potential. Even these can be done quite simply. Note that:

• •

When Eckart displacements are obtained by rotation,

[TABLE]

where $\bm{s}_{r}(\bm{a}^{\textrm{E}})$ is the Wilson $\bm{s}$-vector corresponding to the $r$th vibrational coordinate calculated at the $\bm{a}^{\textrm{E}}$ configuration. Note that the calculation of $\frac{\partial s_{r}}{\partial\bm{a}}$ is analytical, and clearly, there is no need for differentiating the matrix rotating into the Eckart frame.

• •

When Eckart displacements are obtained by projection the derivatives $\frac{\partial c_{j}}{\partial s_{r}}$ can be obtained analytically

[TABLE]

Then, the desired derivatives are obtained by inverting a matrix $\bm{D}$ whose elements are $D_{jr}=\frac{\partial c_{j}}{\partial s_{r}}$.

One might be interested in calculating the derivatives of Eckart coordinates. Numerical Eckart codeswang-carrington ; szidarovszky-fabri-csaszar ; sadri ; yachmenev-yurchenko do, in fact, calculate such derivatives. Since

[TABLE]

it follows that

[TABLE]

Therefore, when projection is employed to calculate Eckart displacements the calculation of the derivatives of Eckart coordinates with respect to the vibrational coordinates is trivial. When one employs rotation to generate Eckart displacements one can calculate the derivatives $\frac{\partial s_{r}}{\partial c_{j}}$ simply by the method described above. Then, the derivatives $\frac{\partial c_{j}}{\partial s_{r}}$ are obtained by inverting a matrix $\bm{A}$ whose elements are defined as $A_{rj}=\frac{\partial s_{r}}{\partial c_{j}}$. As the first of the numerical examples demonstrating the working of the gateway method we shall compare this method of differentiating Eckart coordinates with the simplest numerical differentiation scheme,

[TABLE]

where $h$ is small distortion along a vibrational coordinate, and $\bm{U}$ rotates the coordinates of a given distorted configuration into Eckart coordinates.

The higher order derivatives required to evaluating the pseudo-potential can be obtained similarly, without resorting to numerical differentiation.

IV Numerical examples

All results of numerical calculations to be presented refer to the H2O molecule. The equilibrium geometry and the potential energy surface obtained by Jensen by fitting to experimental data jensen are used in the calculations. Valence internal coordinates, the two bond lengths, $r_{1}$ and $r_{2}$, and the bond angle $\phi$ are used as vibrational coordinates.

The origin of the MS (Eckart-frame) is fixed to center of mass of the reference configuration and its axes are parallel to the axes of the coordinate system, the initial axis-system, shown in Figure 3.

$\bm{a}\left(r_{1},r_{2},\phi\right)$ are the coordinates of the atoms of distorted configurations whose center of mass coincides with the origin of the MS. They can be obtained from the atom coordinates summarized in Table  5 by translation. They are used along with the coordinates $\bm{a}^{0}=\bm{a}\left(r_{1}^{0},r_{2}^{0},\phi^{0}\right)$ of the reference configuration, where the equilibrium bond lengths, $r_{1}^{(e)}=r_{2}^{(e)}=0.95843~{}\text{\AA}$, and the equilibrium bond angle, $\phi^{(e)}=104.43976^{\circ}$, are taken as $r_{1}^{0},r_{2}^{0}$ and $\phi^{0}$, respectively, to calculate Eckart displacements and Eckart coordinates. To obtain Eckart displacements (and Eckart coordinates) we use the rotation and projection methods. The matrix $\bm{U}$ rotating to Eckart coordinates is obtained by using the method of Ref.  krasnos, . The basis vectors of the vibrational space given in Table 6 are obtained by employing the analytical formulas derived in Appendix A.

IV.1 Differentiation of Eckart coordinates

Table VII compares the derivatives of Eckart coordinates with respect to the bond angle at a distorted configuration. As expected in comparison of approximate numerical differentiation results with analytical ones, the deviations of the gateway and the numerical differentiation results are the smallest at an intermediate $h$. Thus, the results by the gateway method do qualify as ”analytical (exact)” results.

IV.2 The vibrational Hamiltonians

The potential energy surface and the elements of the vibrational matrix $\mathcal{G}$ are functions of the vibrational coordinates. Various 1D cuts of these surfaces, Figures  4, 5, 6, 7, and 8, show that the vibrational operators obtained by employing the rotation and projection methods, respectively, to construct Eckart displacements are different.

The vibrational $\mathcal{G}$ matrix should not depend on rotation to Eckart coordinates. Thus, it should be identical with Wilson’s G-matrix considered as function of the internal coordinates. Figures  4, 5, 6, and 7 show that this is indeed the case, thus providing another numerical evidence of the correctness of the gateway method.

To calculate the pseudo-potential we use the expression

[TABLE]

given in Ref.  johnson-reinhardt, when Eckart displacements are obtained by the rotation method. When using projection to obtain Eckart displacements, the pseudo-potential is calculated by employing the formulas given in Table  4. Even this complicated expression of the pseudo-potential can be calculated without resorting to numerical differentiation.

Tables  8, 9, 10, and 11 compare some of the converged vibrational energy levels as obtained by the gateway methods with the results of highly sophisticated DVR calculations by Bramley and Carrington bramley-carrington . They employed Radau coordinates along with the same geometry and potential data as the ones used here. The agreement between the gateway methods and Ref.  bramley-carrington, is very good. The agreement between the different gateway methods is very good either. Vibrational levels were also calculated without pseudo-potential. One can see that the energy levels involving only stretching are less influenced by removing or including the pseudo-potential than energy levels involving bending of the molecule. The contribution of the pseudo-potential is more significant in the case of the Eckart Hamiltonian calculated by the projection method.

By examining Figures  4, 5, 6, 7, and 8 one might wonder how it is that different curves correspond to the same coordinates. Recall that these coordinates are those of the configuration $\bm{a}$ but the function values correspond to those taken at Eckart configurations $\bm{a}^{\textrm{E}}$. The Eckart coordinates are obtained by rotation, $\bm{a}^{\textrm{E}}=\mathbb{U}\bm{a}$, and projection, $\bm{a}^{\textrm{E}}=\bm{P}\bm{a}$, respectively. $\mathbb{U}$ is a $3N\times 3N$ block diagonal matrix whose diagonal blocks are identical $3\times 3$ rotational matrices $\bm{U}(s)$ determined by the method of Ref.  krasnos, . Rotation does not change bond angles and bond lengths whereas projection does. Therefore, the values of internal coordinates corresponding to the different configurations are related as

[TABLE]

Thus, for instance, one has the relationships

[TABLE]

Clearly, the Eckart Hamilton operators obtained by projection and rotation are different. Nevertheless, they must have the same spectrum, as, in fact, suggested by the results of numerical calculations.

To show that analytically, we start by noting that with $\bm{P}\bm{a}(s)$ being Eckart coordinates corresponding to a configuration with internal coordinate values $\mathfrak{s}$, one can always find rotation $\mathbb{R}(\mathfrak{s})$ such that $\bm{P}\bm{a}(s)=\mathbb{R}\bm{a}(\mathfrak{s})$, where $\mathbb{R}(\mathfrak{s})$ is a $3N\times 3N$ block diagonal matrix whose diagonal blocks are identical $3\times 3$ rotational matrices $\bm{R}(\mathfrak{s})$ determined by the method of Ref.  krasnos, .

The Eckart Hamilton operators $\hat{H}^{\textrm{E}}\left(\mathbb{U}\bm{a}(s);s,\hat{p}_{s}\right)$ and $\hat{H}^{\textrm{E}}\left(\mathbb{R}\bm{a}(\mathfrak{s});\mathfrak{s},\hat{p}_{\mathfrak{s}}\right)$ one obtains by using the Eckart coordinates $\bm{a}^{\textrm{E}}(s)$ and $\bm{a}^{\textrm{E}}(\mathfrak{s})$, respectively, are identical since they differ only in notation. Since rotation of the molecule in the coordinate system of axes cannot change the spectrum of a ro-vibrational Hamilton operator, it is also true that the spectrum of the non-Eckart Hamilton operator $\hat{H}\left(\bm{a}(s);s,\hat{p}_{s}\right)$ is the same as that of the operator $\hat{H}^{\textrm{E}}\left(\mathbb{R}\bm{a}(\mathfrak{s});\mathfrak{s},\hat{p}_{\mathfrak{s}}\right)$. There remains to show that the spectrum of $\hat{H}^{\textrm{E}}\left(\mathbb{R}\bm{a}(\mathfrak{s});\mathfrak{s},\hat{p}_{\mathfrak{s}}\right)$ is, the same as that of the operator $\hat{H}^{\textrm{E}}\left(\bm{P}\bm{a}(s);s,\hat{p}_{s}\right)$. If it can be shown that these operators are related by invertible coordinate transformation, than, as a consequence they must have the same spectrum. A tentative proof is given below.

Let $\mathscr{D}_{s}\subseteq\mathscr{R}_{1}\times\mathscr{R}_{2}\times\mathscr{R}_{3N-6}$, a subset of $3N-6$–tuples of real numbers $\mathscr{R}$, be the domain accessible by the internal motions $s$. Projection of configurations $\bm{a}(s)$ implies the map $\varphi\left(\mathscr{D}_{s}\right)=\mathscr{D}_{\mathfrak{s}}$. Since in our case $\mathscr{D}_{\mathfrak{s}}$ is defined as the range of the map, the map is evidently onto. It is also one to one, if for every $\mathfrak{s}\in\mathscr{D}_{\mathfrak{s}}$ there is a unique $s\in\mathscr{D}_{s}$ such that $\mathfrak{s}=\varphi\left(s\right)$, that is $\varphi\left(s\right)=\varphi\left(s^{{}^{\prime}}\right)$ implies $s=s^{{}^{\prime}}$. Assume that $\varphi$ maps $s$ and $s^{{}^{\prime}}(\neq s)$ into $\mathfrak{s}$, that is $\bm{P}\bm{a}(s)=\mathbb{R}(\mathfrak{s})\bm{a}\left(\mathfrak{s}\right)$ and $\bm{P}\bm{a}(s^{{}^{\prime}})=\mathbb{R}^{{}^{\prime}}(\mathfrak{s})\bm{a}\left(\mathfrak{s}\right)$. By construction $\bm{P}\bm{a}(s)$ and $\bm{P}\bm{a}(s^{{}^{\prime}})$ are Eckart coordinates. Since the matrix transforming $\bm{a}\left(\mathfrak{s}\right)$ into Eckart coordinates is uniquely determined by the method of Ref.  krasnos, the matrices $\mathbb{R}(\mathfrak{s})$ and $\mathbb{R}^{{}^{\prime}}(\mathfrak{s})$ must be equal to this matrix. Then $\bm{P}\bm{a}(s)=\bm{P}\bm{a}(s^{{}^{\prime}})$ follows and we have the relationships

[TABLE]

and

[TABLE]

where $\bm{Q}$ is the projection matrix onto the translational-rotational subspace of the configuration space (i.e. $\bm{I}=\bm{P}+\bm{Q}$ with $\bm{I}$ denoting a $3N$ by $3N$ identity matrix). These equations show that changing the values of the internal coordinates from $s$ to $s^{{}^{\prime}}$ causes displacements having no components in the vibrational space. But this cannot happen since $s$ are genuine internal coordinates whose change should lead to non-zero displacement in the vibrational space as well. That is, $s^{{}^{\prime}}$ must be equal to $s$. Therefore, the map $\mathfrak{s}=\varphi(s)$ is one to one and onto, which means that it is invertibleszekeres .

In a variational calculation any Hamiltonian, $\hat{H}$ or $\hat{H}^{\textrm{E}}$, with exact KEO can be used. Nevertheless, the calculations may converge faster or more slowly depending on which Hamiltonian is employed. When one attempts to simplify a Hamiltonian by introducing approximations the quality of the approximate Hamiltonian derived may strongly depends on to which Hamiltonian are the approximations invoked. Eckart Hamiltonians may have advantage, whether the rate of convergence or the quality of the derived approximate Hamiltonians are considered. However, there are many different Eckart Hamiltonians (of which two have been considered in the numerical examples). Therefore, it may be useful to find a way of selecting a unique Eckart Hamiltonian out of the many different ones. This question can be reduced to finding Eckart displacements optimal in some appropriate sense. It is addressed in the next Section.

V Geometry of Eckart conditions: Optimal Eckart displacements

Pictures are often helpful in explaining and learning ideas. Figure  9 gives a simplified pictorial representation of the geometry of Eckart conditions. One can see immediately that to obtain Eckart coordinates (and displacements) one must find a map $\mathcal{M}$ mapping a general point of the translation reduced configuration space, $\bullet$, into a point, $\otimes$, of the vibrational space. One such map, namely, $\mathcal{M}=\bm{P}$, has been already considered in the previous Sections. Note that $\bm{P}$ solves the optimization problem $\min\limits_{\mathcal{M}}\lVert\bullet-\otimes\rVert^{2}$, and to each point $\bullet$ there corresponds a unique point $\otimes$ in vibrational space. ( $\lVert v\lVert$ denotes the length of a vector $v$.) Now imagine connecting the point $\oplus$ representing the origin of the MS with the point $\bullet$ by a straight line and drawing the hypersphere of radius $\lVert\bullet-\oplus\lVert$ in the (translation reduced) configuration space. Intersections of this sphere with the vibrational space give Eckart coordinates. In general, however, this does not give a unique image of the point $\bullet$. Therefore it is not a map schutz . To simplify, let us restrict to transformations of the form $\mathbb{U}=\mathrm{diag}\left(\bm{U}_{1},\bm{U}_{2},\ldots,\bm{U}_{3N-6}\right)$ with $\bm{U}_{1}=\bm{U}_{2}=\cdots=\bm{U}_{3N-6}=\bm{U}$ and $\bm{U}\in SO(3)$. Figure  10 shows such a transformation ${\textit{r}ed}\bullet=\mathbb{U}\bullet$ corresponding to mass weighted displacement $\bm{m}^{1/2}\bm{d}={\textit{r}ed}\bullet-{\textit{b}lue}\bullet=\mathbb{U}\bm{m}^{1/2}\bm{a}-\bm{m}^{1/2}\bm{a}^{0}$. The square of the length of this displacement, which may be called mass weighted squared displacement (MWSD), is $\mathrm{MWSD}=\lVert{\textit{r}ed}\bullet-{\textit{b}lue}\bullet\rVert^{2}=\left[\bm{d}\right]^{T}\bm{m}\bm{d}$. It has been shown in Ref.  krasnos, that the displacement vector calculated by $\mathbb{U}$ solving the problem of minimization $\min\limits_{\mathbb{U}}\left(\left[\bm{d}\right]^{T}\bm{m}\bm{d}\right)$ obeys the Eckart conditions. In other words, of all equivalent configurations, that is among all $\mathbb{U}\bm{a}$, the closest to the reference configuration is related to the reference configuration by Eckart displacements. Experience shows that it is uniquely determined.

Projection and rotation may be combined to generate Eckart coordinates as depicted in Figure  11. Observe that the points ${\textit{r}ed}\otimes=\bm{P}\mathbb{U}\bm{m}^{1/2}\bm{a}$ correspond to Eckart configurations for any choice of $\mathbb{U}$.

Therefore one can ask for an $\mathbb{U}$ leading to an Eckart configuration which is the closest to the reference configuration among all Eckart configurations which one can derive from a single distorted configuration by rotation. Translated into mathematical form it means minimization (with respect to $\mathbb{U}$) of the quantity

[TABLE]

Since, as explained in Ref.  szalay-eckart, , the projection matrix onto the vibrational space can be also expressed as

[TABLE]

where

[TABLE]

with $M$ denoting the mass of the molecule, and

[TABLE]

where $\varepsilon_{\delta\beta\gamma}$ is the Levi-Civita tensor and $\mathbf{I}^{0}$ stands for the rotational tensor of inertia of the reference configuration, one can obtain that

[TABLE]

In deriving Eq.  (40) use has been made of the fact the $\bm{a}^{0}_{n}$ and $\bm{a}_{n}$ obey the Eckart conditions and the translational Eckart conditions, respectively.

If $\bm{U}$ is chosen such that $\bm{U}\bm{a}_{n}$ are Eckart coordinates, then Eq.  (40) shows that ${\textit{b}lue}{\mathrm{MWSD}}=\min\limits_{\bm{U}}\mathrm{MWSD}$. Therefore, it follows that

[TABLE]

Results of numerical calculations shown in Figures  12 and 13 give numerical evidence. Some details of solving the optimization problem $\min\limits_{\bm{U}}{\textit{b}lue}{\mathrm{MWSD}}$ are described in Appendix B. One can see in Figure  12 that the optimization occasionally bogs down at local minima. At these points the optimization must be restarted by new initial parameters and the global minimum can be found. Nevertheless, the numerical results presented do confirm the inequality Eq.  (41).

VI Summary

The Eckart conditions treated as what they are, a homogeneous system of linear equations, has led to the gateway Hamiltonian method szalay-gateway ; szalay-aspects . In this approach the conditions define the space of vibrations and one uses them to determine the molecule fixed system coordinates of the atoms as functions of vibrational coordinates, and, eventually, to derive Eckart ro-vibrational Hamiltonians. The derivations assume that a basis for the vibrational space is already known. A prerequisite for practical use of the gateway method is, therefore, the availability of such a basis set. In the present work general analytical formulas of such a basis set have been derived and some of the practical advantages of the gateway Hamiltonian method have been numerically demonstrated. It has been pointed out that there is an infinite number of Eckart KEOs corresponding to any given set of curvilinear internal coordinates. It is suggested that one should use the KEO corresponding to in some sense optimally defined Eckart displacements. A possible definition of optimal Eckart displacements has been introduced and illustrated with numerical examples.

The results described may give food for thoughts:

• •

Construction and application of a Hamilton operator with optimal Eckart displacements.

• •

The geometrical interpretation of the Eckart conditions may elucidate and can help deriving the relationships of Eckart coordinates corresponding either to different electronic states duschinsky or isotopic speciesmeier of a molecule.

• •

Conditions other than the Eckart, e.g. those for using the instantaneous principal axis system, might be treated similarly by replacing linear algebra with the appropriate mathematical technique(s).

Acknowledgements.

The author is indebted to Ildikó Horváth for her support.

Appendix A Analytical solution of the Eckart conditions

Let $\bm{t}^{\alpha}$ and $\bm{r}^{\alpha}$ denote the row vectors of elements

[TABLE]

and

[TABLE]

By employing this notation the Eckart conditions,

[TABLE]

can be written as a system of linear equations

[TABLE]

with

[TABLE]

Since Eq.  (45) is a homogeneous system, its general solution, $\bm{d}^{\textrm{E}}$, can be written as

[TABLE]

where $K=3N-\textrm{rank}(\bm{E})$, $b_{j}$ are free variables, and $\bm{h}^{\textrm{E},j}$ are particular solutions (that is they obey the Eckart conditions). We shall determine $\bm{h}^{\textrm{E},j}$. Then, $\bm{h}^{\textrm{E},j}$ are employed to calculating a basis orthonormal with mass weighting and spanning the vibrational space.

Gaussian elimination carl-meyer applied to Eq.  (45) gives

[TABLE]

where

[TABLE]

and

[TABLE]

By introducing

[TABLE]

one can write $\bm{C}$ as

[TABLE]

Thus, we have the system of equations

[TABLE]

Then,

[TABLE]

With the dependent variables removed from the right-hand sides in Eq.  (79) one has

[TABLE]

By collecting the coefficients of the independent variables

[TABLE]

given in Eq.  (80) into columns one obtains the vectors

[TABLE]

respectively, given in Tables  12 and 13. One can check by simple analytical calculations employing Eqs.  (68,69) and (70) that these $\bm{h}^{\textrm{E},j}$ vectors do, indeed, obey the Eckart conditions Eq.  (44). For instance, by replacing the components of $\bm{d}$ with those of $\bm{h}^{\textrm{E},3}$ in Eq.  (44f) one obtains that

[TABLE]

One may have noticed that some terms can become singular. The singularities can be avoided simply by reorienting the coordinate system.

The vectors $\bm{h}^{\textrm{E},j}$ may not be orthonormal with mass weighting, i.e. $\left[\bm{h}^{\textrm{E},i}\right]^{T}\bm{m}\bm{h}^{\textrm{E},j}\neq\delta_{ij}$. One may use the Gram-Schmidt procedure to orthogonalize the vectors $\bm{m}^{1/2}\bm{h}^{\textrm{E},j}$ to obtain a set of orthonormal vectors $\bm{e}_{j}$:

[TABLE]

where $\|\bm{e}_{j}\|=\left(\left[\bm{l}_{j}\right]^{T}\bm{l}_{j}\right)^{1/2}$, and $\bm{I}$ is the $3N$ by $3N$ unit matrix. Then, basis vectors, $\bm{d}^{\textrm{E},j}$, spanning the vibrational space and orthonormal with respect to mass weighting can be calculated by the equations

[TABLE]

Appendix B Optimal Eckart displacements

For minimizing MWSD of Eq.  (40) a neater expression than Eq.  (40) arises by considering that

[TABLE]

Since a dot product $\bm{u}\cdot\bm{v}$ can be written as

[TABLE]

one obtains

[TABLE]

where use has been made of the equality

[TABLE]

To carry on recall that a rotation matrix $\bm{U}$ may be parametrized in terms of the scalar, $q_{0}$, and vector, $\bm{q}=\left(q_{1},q_{2},q_{3}\right)$, components of a quaternion $Q=\left[q_{0},\bm{q}\right]$ of unit norm, and a rotated vector $\bm{U}\bm{a}_{n}$ can be calculated as products of quaternions, since

[TABLE]

where $\left[0,\bm{U}\bm{a}_{n}\right]$ is the pure quaternion corresponding to the vector $\bm{U}\bm{a}_{n}$, and $A_{n}=\left[0,\bm{a}_{n}\right]$ denotes the pure quaternion corresponding to the vector $\bm{a}_{n}$.

Therefore we can rewrite Eq.  (315) as

[TABLE]

where $D^{\textrm{E},j}_{n}$ is the pure quaternion corresponding to the vector $\bm{d}_{n}^{\textrm{E},j}$, and the second equality follows from the fact that the norm of the product of quaternions is equal to the product of their norm. By expressing the quaternions in Eq.  (318) in terms of their components one can see that $c_{j}$ is a quadratic form

[TABLE]

with $\bm{q}^{T}=\left(q_{0},q_{1},q_{2},q_{3}\right)$ denoting a row vector containing the components of the quaternion $Q$.

Therefore,

[TABLE]

and it has to be minimized with respect to $\bm{q}$ under the normalization condition $\lVert\bm{q}\rVert=1$. The optimal $\bm{q}$ may be determined iteratively.

Note that

[TABLE]

Therefore,

[TABLE]

and we can find $\theta$ which minimizes MWSD at fixed $\bm{n}$. Next update $\theta$ and find the optimal $\bm{n}$ at this new, fixed $\theta$.

Note that

[TABLE]

Thus, with $\theta$ and $\phi$ fixed, we can find the optimal $\chi$ by minimizing

[TABLE]

Now having updated $\chi$, with $\theta$ and $\chi$ fixed we look for the optimal $\phi$. By using

[TABLE]

it can be derived that

[TABLE]

Minimization of this function of $\phi$ gives the optimal $\phi$. Then, update $\phi$ and iterate the procedure until convergence.

The 1D optimizations required are reduced to finding zeros of functions of single variable. To carry out calculations a MAPLEmaple code has been written.