The rovibrational Hamiltonian for ammonia-like molecules

The rovibrational Hamiltonian for ammonia-like molecules

Spectrochimica Acta Part A 58 (2002) 601– 628 www.elsevier.com/locate/saa The rovibrational Hamiltonian for ammonia-like molecules Jan Makarewicz *, ...

294KB Sizes 1 Downloads 52 Views

Spectrochimica Acta Part A 58 (2002) 601– 628 www.elsevier.com/locate/saa

The rovibrational Hamiltonian for ammonia-like molecules Jan Makarewicz *, Alexander Skalozub 1 Faculty of Chemistry, A. Mickiewicz Uni6ersity, Grunwaldzka 6, PL 60 -780 Poznan, Poland Received 18 July 2001; accepted 17 September 2001

Abstract A new exact quantum mechanical rovibrational Hamiltonian operator for ammonia-like molecules is derived. The Hamiltonian is constructed in a molecular system of axes, such that its z% axis makes a trisection of the pyramidal angle formed by three bond vectors with the vertex on the central atom. The introduced set of the internal rovibrational coordinates is adapted to facilitate a convenient description of the inversion motion. These internal coordinates and the molecular axis system have a remarkable property, namely, the internal vibrational angular momentum of the molecule equals zero. This property significantly reduces the Coriolis coupling and simplifies the form of the Hamiltonian. The correctness of this Hamiltonian is proved by a numerical procedure. The orthogonal Radau vectors allowing us to define a similar molecular axis system and the internal coordinates are considered. The Hamiltonian for the Radau parameterization takes a form simple enough to carry out effectively variational calculations of the molecular rovibrational states. Under the appropriate choice of the variational basis functions, the Hamiltonian matrix elements are fully factorizable and do not have any singularities. A convenient method of symmetrization of the basis functions is proposed. © 2002 Elsevier Science B.V. All rights reserved. Keywords: Rovibrational Hamiltonian; Inversion motion; Umbrella-like coordinates; Vibrational angular momentum; Coriolis coupling

1. Introduction Theoretical prediction of the molecular rovibrational spectra from the first principles is based on the solution of the differential Schro¨dinger equation describing nuclear molecular motions. Similarly as in the classical mechanics, ‘‘the main difficulty in integrating a given differential equation lies in introducing convenient variables, * Corresponding author. Fax: + 48-61-8658008. E-mail address: [email protected] (J. Makarewicz). 1 On leave from: Ukrainian State University of Chemical Technology, pr. Gagarina 8, 320005 Dnepropetrovsk, Ukraine.

which there is no rule for finding’’ [1]. This statement of Jacobi is valid up-to-date since the choice of convenient dynamical variables is an essential point in modern methods developed to solve the molecular rovibrational problem. At present, the most widely used accurate computational techniques are based on variational and discrete representation methods implemented mainly for triatomic [2–12] and tetraatomic [13 –26] molecules. In both methods the Hamiltonian operator defining the molecular Schro¨dinger equation must be represented by the Hamiltonian matrix in some functional basis set. The computer requirements of these methods are determined mainly by the size of the required functional basis

1386-1425/02/$ - see front matter © 2002 Elsevier Science B.V. All rights reserved. PII: S 1 3 8 6 - 1 4 2 5 ( 0 1 ) 0 0 6 6 0 - 6

602

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

set defined in a space of the selected coordinates describing the molecular motions. Thus, the choice of these coordinates affects greatly the efficiency of the variational calculations. The optimal coordinate set that would minimize the computational effort of calculating the eigenvalues and eigenfunctions of the Hamiltonian matrix is desirable. Such a set should satisfy a few criteria. The coordinates should span the whole configuration space accessible for the nuclei and the boundary conditions for this space should be simple enough to enable easy computation of the Hamiltonian matrix elements. The essential singular terms in the Hamiltonian should not appear since they need a special treatment to avoid the problems in evaluation of the Hamiltonian matrix elements. The coordinates should be compatible with the symmetry inherent in the Hamiltonian in such a way that the symmetry operators can be executed on the basis functions. Then, it is possible to split the Hamiltonian matrix into smaller sub-matrices representing each symmetry type. The reduction of the Hamiltonian matrix dimension considerably enhances the efficiency of the matrix diagonalization. The interaction terms in the kinetic energy and the nuclear potential energy function should be minimized in a large domain of the coordinate space. The Hamiltonian in such coordinates becomes nearly separable, i.e. it can be approximately expressed as a sum of uncoupled terms depending on separate subsets of the whole coordinate space. The approximate separability enables a construction model Hamiltonians of lower dimensions. The eigenfunctions of such Hamiltonians provide efficient contracted basis sets for accurate variational solution in full dimension. In optimal coordinates the Hamiltonian should be factorizable. This means that each its term should be represented as a product of operators depending on nonintersecting subsets of coordinate space. The most desired is full factorization where each term is a product of one-dimensional operators. The matrix elements of a fully factorizable Hamiltonian can be easily determined by one-dimensional numerical or analytical integrations resulting in substantial reduction of computational effort.

The obvious coordinates that are precursors of any rovibrational curvilinear coordinates are 3N Cartesian components of N nuclear vectors measured relative to the laboratory system of axes (LSA). These coordinates fulfil most of the optimum criteria mentioned above, except for the separability and factorizability of the potential energy function. This function depends only on relative distances between the nuclei, not on their positions, so it cannot be presented in a separable and factorizable form. More convenient coordinates adapted to describe the potential function are obtained by separating out three LSA components of the center of mass vector and choosing three rotational and 3N− 6 remaining coordinates describing the internal nuclear motions. In this way, a partial separation is achieved in the potential energy because it depends only on 3N− 6 translationally and rotationally invariant internal coordinates. The rotational coordinates are defined by introducing the molecular system of axes (MSA) attached to moving nuclei [27]. The orientation of MSA with respect to LSA is usually described by three Euler angles, which play the role of the rotational coordinates. The choice of MSA affects the Coriolis interaction between the rotational and the internal molecule motions. The physical origin of this interaction is the internal angular momentum generated by the nuclei moving in the rotating molecule. If MSA is defined in such a way, that the internal angular momenta generated by vibrational modes are small, then the rotational motion is well separated from the internal motion. It is well known that MSA satisfying the Eckart conditions [27] cancels the Coriolis coupling at the equilibrium configuration of a molecule executing small amplitude vibrations. Thus, the Eckart axes are optimal for defining the rotational motion of semirigid molecules. The Eckart idea has been extended to molecules performing large amplitude motions [28–39]. In almost all variational calculations on the rovibrational states of triatomic and tetraatomic molecules the Eckart axes have not been used. Instead, MSA has been embedded to two internuclear vectors because of the simplicity of the Coriolis terms. However, these terms have not been minimized.

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

The problem of the optimal MSA has been studied in detail for a triatomic molecule [35]. It has been proved that the optimal MSA overlaps with the system of internal axes which cancel the Coriolis coupling at the equilibrium or at an averaged molecular configuration. The rovibrational Hamiltonian obtained with the minimized rotation –vibration interactions conserved a simple factorizable form. For a symmetric triatomic molecule, the Coriolis interaction has been reduced using a bisector-embedded MSA [36]. It is difficult to find the optimal internal coordinates even for simple molecules. Familiar types of the internal coordinates do not satisfy all the optimum criteria. The rectilinear normal coordinates [40] widely used for semirigid molecules generate the unfactorizable pseudopotential term in the Watson Hamiltonian [41]. This term becomes singular at linear configurations and causes serious problems. The disadvantages of these coordinates in the treatment of large amplitude motions are commonly known, so other curvilinear coordinates are exploited in this problem. The valence coordinates defined by the bond lengths, interbond angles and dihedral angles have useful properties. They are the local spherical polar coordinates that parametrize a chosen set of N − 1 internuclear vectors. For typical molecules, the potential energy function is approximately separable if the internuclear vectors are selected according to the pattern of the chemical bonds. However, it is not always possible to define these vectors in accordance with the bonding pattern. For example, a cyclic molecule with N nuclei has N bonds but only N −1 internuclear vectors are allowed to describe its geometry. In such a case the valence coordinates are not optimal. The disadvantages of these coordinates are compensated by the factorizable form of the kinetic energy operator T. , which is especially simple when the internuclear vectors are chosen as some orthogonal vectors. This operator was first derived for three body systems using the valence coordinates based on the bond vectors [42– 45], orthogonal Jacobi vectors [46,47] and more general internuclear vectors [36,48,49]. The general exact operator T. for polyatomic molecules ex-

603

pressed by the bond-angle coordinates based on the Jacobi vectors has been derived by several authors [20,23,50]. The first correct rovibrational operator T. in valence coordinates based on the nonorthogonal internuclear vectors has been derived for tetraatomic molecules [51,52] by direct transformation of the Laplacian operator from the LSA Cartesian coordinates to the rotational and internal coordinates [53]. The computer algebra has been employed to overcome the complexity of the derivation. A pure vibrational part of the operator T. for some polyatomic molecules has been obtained using the same derivation method [54– 56]. Recently the Cartesian operator T. of a polyatomic molecule has been transformed to the full rovibrational operator T. expressed by the internal valence coordinates and the Euler angles defined for a simple two-vector embedded MSA [24]. All terms depending on the Euler angles have been replaced by the angular momentum operators in a clever way at the end of the derivation. The complexity of a direct transformation method appeared due to rotational noninvariance of the Euler angles used as the rotational coordinates. This method has been significantly simplified by introducing an infinitesimal rotational coordinate system yielding T. directly in terms of the angular momentum operators, instead of the Euler angles [57,58]. An invariant form of the classical rovibrational kinetic energy T was obtained a long time ago by Eckart [27]. He transformed T in the classical Lagrangian form, expressed by the LSA Cartesian velocities, to the Lagrangian written as a quadratic form of the vibrational and angular velocities, specified by the kinetic energy matrix K. The Hamiltonian form of T expressed by a quadratic form G of generalized vibrational momenta and angular momenta was obtained by inverting the matrix K. It is well known that the classical kinetic energy T can be transformed to the correct operator T. using the Podolsky rule [59]. Thus, the approach based on the classical mechanics dispenses from the cumbersome Euler angles. The general formulas for the elements of the vibrational G matrix for the nonorthogonal valence coordinates have

604

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

been derived by Decius [60]. With the known G matrix, the exact vibrational quantum mechanical operator T. for polyatomic molecules has been determined by the Podolsky method only recently [61,62]. The rotational invariance is also apparent in the direct transformation method developed by Louck et al. [63,64]. Their method has not been widely used although it offers an elegant, and general route [65,66]. In this paper we derive and discuss interesting features of the exact rovibrational operator T. for ammonia-like molecules in the internal coordinates adapted to the inversion motion. These coordinates, further called the ‘umbrella-like’ (UL) coordinates, have been used recently by Handy et al. [22] and independently by the authors [67]. In Section 2 we reveal advantages of the UL coordinates and introduce the MSA appropriate to describe the rotational motion. In Section 3 we focus on the derivation of the full rovibrational operator T. in the MSA rotational coordinates and the internal UL coordinates. We have tried first the classical route by constructing and inverting the matrix K. The determination of K was not too difficult, but our attempt to invert this matrix failed, although we used symbolic algebra programs. For this reason we exploit here the method given in Ref. [64] allowing a direct transformation of the linear momentum operators from the LSA to MSA Cartesian coordinates. As a result, the operator T. is further expressed by the internal coordinates. To determine the rotational part of T. we construct in Section 4 the MSA basis in terms of the nuclear bond vectors. The dynamical properties of this basis are characterized in Section 5 by the parameters describing a sensitivity of the rotating MSA to the variations of the moving nuclear bond vectors observed in the LSA. The explicit forms of the vibrational linear and angular momenta are given in Sections 6 and 7. They allow a determination, in Sections 8 and 9, of all terms contributing to the rovibrational operator T. . In Section 10 we simplify T. by introducing the Radau nuclear vectors and propose a convenient method for the symmetrization of this operator.

2. Internal coordinates and molecular system of axes for ammonia-like molecules In applications it is not always clear how to construct good internal coordinates explicitly. This difficulty arises in the case of a simple highly symmetric NH3 molecule for which the exact rovibrational Hamiltonian is known only in valence coordinates [51,52]. However, these coordinates determined from three equivalent NH bond vectors are not compatible with the permutation symmetry S3 of this molecule. The set of the angular valence coordinates contains only two equivalent bond angles q12 and q13, and one dihedral angle ~, i.e. the ‘book’ angle between the two planes NH1H2 and NH1H3. This set is not invariant under the cyclic permutation (123) of the protons. One can try to use the invariant set of three equivalent bond angles (q12, q13, q23). However, such a coordinate set is not complete because it does not span the whole configuration space. In this set, a given nuclear configuration cannot be distinguished from the inverted one. Recently, the three bond angles have been proposed to derive the Hamiltonian for NH3-like molecules without inversion motion [57]. Although the final expressions have not been presented, the complexity of this Hamiltonian is apparent because the Jacobian of the transformation from the Cartesian to the bond angle coordinates is unfactorizable. Also, the chosen MSA with the z% axis parallel to the sum of the bond vectors r1 + r2 + r3 generates unfactorizable Coriolis and rotation singular terms in this Hamiltonian. Additional difficulty arises from the coupled integration ranges for the bond angles. The integration range for q23 depends on the values of the angles q12 and q13, so it would be difficult to calculate the Hamiltonian matrix elements with the angular basis functions. The total rovibrational problem for the ammonia molecule can eventually be treated using the valence coordinates with the book angle ~ as the inversion coordinate. The calculations could be performed ignoring the cyclic permutations and exploiting the lower permutation symmetry group S2. However, it has been recently found that such a coordinate set is inappropriate to describe the

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

inversion motion in the framework of flexible models [68]. Such models are unavoidable for larger molecules containing NH2 group, which executes the inversion motion. In the flexible nonrigid inverter model developed by Papousˇek and S& pirko [69–73], the small amplitude vibrations adjust adiabatically to the varying inversion coordinate. During the inversion motion the flexible geometry of the molecule varies along the minimum energy path along the inversion coordinate. This path can be determined by employing the standard ab initio methods. The flexible model for the ammonia molecule has been studied for the book angle ~ as the inversion coordinate and the remaining valence coordinates as flexible parameters. It has been found that the flexible model can be defined only in a limited range of this inversion coordinate [68]. At the value of ~ equal 60°, a bifurcation appears in the behavior of the optimized internal coordinates. As a result, the molecular geometry optimized along the minimum energy path exhibits a broken symmetry lower than C36. At ~:45.5° the flexible model cannot be defined at all, since the optimization process becomes divergent owing to a catastrophe in the potential energy function. This fact motivated us to search for another convenient inversion coordinate [67]. We have found, that the bifurcation can be shifted to a very high energy region if the potential energy function can be represented by the UL coordinates. Thus, they are proper coordinates for a description of the inversion motion even in highly excited states.

Fig. 1. The local molecular system of axes and the internal angular coordinates for the ammonia-like molecule.

605

The UL coordinates are presented in Fig. 1. The angle h as the inversion coordinate is defined in conjunction with the local system of axes. Its origin is placed on the central N atom. The z% axis of this system makes the same angle h with all three NH bonds. The x% axis lies in the plane formed by the first bond and the z% axis. The y% axis runs perpendicular to the x%z% plane forming a right-handed x%y%z% system of coordinates. Thus, the z% axis is a trisector of the pyramid of the ammonia-like molecule. The angle ~2 is the polar angle of the second NH bond, defined as a dihedral angle between the x%z% plane and the plane containing the z% axis and the second NH bond. The angle ~3 is defined in an analogous way as the polar angle of the third NH bond. The angles ~2 and ~3 are constant during the inversion motion of the ammonia molecule. Thus, the flexible model describing this motion is much simpler than those defined in terms of the bond angles. We performed the ab initio calculations of the minimum energy path along the inversion coordinate h in several simple amines like NH2OH, NH2CH3, etc. and we found that ~2 and ~3 change a little with varying h. They behave as typical harmonic modes. The UL variables define the Cartesian coordinates of the three nuclear bond vectors r%ii, where i= 1, 2, 3 labels the vectors and i=x, y, z labels their components in the selected local system of axes. These components are given by r%1x = r1 sin h

r%2x = r2 sin h cos ~2

r%3x = r3 sin h cos ~3

r%1y = 0

r%2y = r2 sin h sin ~2

r%3y = r3 sin h sin ~3

r%1z = r1 cos h

r%2z = r2 cos h

r%3z = r3 cos h

The MSA is obtained by shifting the local system of axes to the nuclear center of mass, without changing its orientation. The presented internal coordinates have already been used by Handy et al. [22]. Their angles (i, q3, q2) correspond to ours (h, ~2, 2p− ~3). They considered the allowed ranges of the angular coordinates and concluded that the proper ranges are 05 q2, q3 5 2p, 05 i5p ‘‘bearing in mind that the potential will forbid the atoms from crossing’’ [22]. Naturally, this is a reasonable ap-

606

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

proximation, however, the rigorous constraints for the angular coordinates are 05~2 5~3 52p, 0 5 h 5p, irrespective of the potential properties. Although ~2 varies independent of ~3, it cannot be larger than ~3. At fixed bond lengths, the coordinates (h, ~2, ~3) and (p − a, 2p −t2, 2p − t3) describe the same molecular configuration. Thus, the constraint ~2 5~3 provides the one to one mapping between the Cartesian and the UL coordinates. The UL internal coordinates introduced have a form similar to that of the standard spherical polar ones, however, the inversion angle h is a collective coordinate for the three bond vectors. We use this similarity to derive the Jacobian of the transformation C“ UL from the Cartesian LSA coordinates to our rotational and the UL internal coordinates. As an intermediate step, let us consider first the transformation C“ V from the Cartesian to the rotational and valence coordinates, which are simply the spherical polar coordinates under the corresponding choice of the local coordinate axes. So, let us choose the auxiliary local axis system with the z¦ axis directed along the bond vector r1, the x¦ axis in the plane formed by the vectors r1 and r2, and the y¦ axis perpendicular to this plane. The Jacobian of the transformation C“ V is well known [61,62]. It is proportional to the invariant volume element dVs in the spherical polar coordinates

The relation between the valence and the UL coordinates follows from the standard spherical geometry formulas t1i = cos2 h+ sin2 h cos ~i and cos ~ sin q12 sin q13 = cos q23 − cosq12 cos q13

(1) where q12 and q13 are the bond angles and ~ is the dihedral (torsional) angle between the planes containing the vectors r1 and r2, and the vectors r1 and r3. The volume element for the internal UL coordinates, can be written as dV =r 21r 22r 23Jdr1dr2dr3dhd~2d~3,

(2)

where t12 cos q12, t13 cos q13 and J = #(t12,t13,~)/#(h,~3,~2) is the Jacobian of the transformation from the internal valence to the UL coordinates. Thus, the Jacobian of the full transformation C“UL is sin r%r 21r 22r 23J, where the standard factor sin r% appears from the MSA Euler angles.

(3)

The form of the Jacobian matrix defines the relations between the derivatives Æ Ç Æ Ã # à Ã#t12 à #h à à #h à # à Ã#t à à = à 12 Ã#~3 à à #~3 à # à Ã#t12 Ã#~ à à #~ È 2É È 2

Æc11 Ã = Ãc21 Ã Èc31

#t13 #h #t13 #~3 #t13 #~2

c12 c22 c32

dVs = (r 21r 22r 23 sin q12 sin q13) dr1 dr2 dr3 dq12 dq13 d~,

for i= 2,3

ÇÆ Ç #~ à à # à #h à Ã#t12 à #~ à à # à Ãà à #~3 à Ã#t13 à #~ à à # à #~2 à à #~ à ÉÈ É Æ Ç Ã # à Ã#t à c13 Ç Ã 12 à à # c23 à à Ã. à Ã#t13 à c33 É Ã # à à #~ Ã È É

(4)

They have to be written exactly in the above order (allowing the cyclic permutations) to conserve the right-handed axis system. After calculating the necessary derivatives we obtain J= det c =4 sin3 h sin

~1 ~2 ~3 sin sin , 2 2 2

(5)

where we introduce also an angle ~1 = ~3 − ~2 to achieve a symmetric description of the molecular geometry. The full Jacobian sin r%r 21r 22r 23J presented here agrees with that obtained earlier by Handy et al. [22] who used a different derivation method. Let us compare the configuration space volumes for the valence angular coordinates

&

1

dt12

&

−1

&

1

dt13

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628 y

d~ =8y

(6)

−y

−1

and for the UL coordinates

&

J dh d~2 d~3

&

y

=4

sin3 h dh

&

0

2y

sin

0

~3 d~3 2

&

~3

sin

0

~1 ~2 sin d~2 2 2

= 8y.

(7)

Notice that when the last integral over ~2 in Eq. (7) is calculated in the interval (0, 2p), instead of (0, ~3), then we obtain zero which is a wrong result. This fact supports our correct choice of the allowed ranges for the angular UL coordinates. In the following, we will derive vibrational momenta for which we need the elements of the inverse matrix c − 1. They are

1 c− 11 =

cot h cos



~3 2



,

~2 ~ ~ sin 2 sin 1 4 1− sin h sin 2 2 2 ~3 2 cot h sin 2 −1 c 12 = − , ~ ~ ~ 2 1−sin2 h sin2 2 sin 2 sin 1 2 2 2 2

2





~3 ~ ~ − sin2 h sin 2 cos 1 2 2 2 1 c− , 13 = − ~ ~ ~1 2 2 2 2 2 2 sin h 1− sin h sin sin sin 2 2 2 ~2 cot h cos 2 −1 c 21 = − , ~ ~ ~ 4 1−sin2 h sin2 3 sin 3 sin 1 2 2 2 sin









~ ~ ~2 − sin2 h sin 3 cos 1 2 2 2 1 c− , 22 = ~ ~ ~1 3 2 2 2 3 2 sin h 1− sin h sin sin sin 2 2 2 ~2 2 cot h sin 2 −1 c 23 = , ~ ~ ~ 2 1− sin2 h sin2 3 sin 3 sin 1 2 2 2 sin









~2 ~3 cos 2 2 1 c− , 31 = − sin h ~1 sin 2 ~2 ~3 sin cos 2 2 1 c− , 32 = 2 cos h ~1 sin 2 ~3 ~2 sin cos 2 2 1 c− . 33 = 2 cos h ~1 sin 2

607

cos

(8)

3. Transformation of the laboratory system to the molecular system of axes The kinetic energy operator T. of a molecule composed of N nuclei with masses mi takes the simplest form in the LSA Cartesian coordinates N 1 1 #2 T. L = − ' 2 % , 2 i = 1 mi #R 2i

(9)

where Ri are the nuclear position vectors measured relative to LSA. The translational motion of the molecule as a whole can be separated from the relative motions of the nuclei by introducing a new set of N−1 internuclear vectors rk as linear combinations N

ri = % Aij Rj

(10)

j=1

and the nuclear center of mass vector rCM in the form of the Nth vector in this set: N

rCM = rN = M − 1 % mjRj,

(11)

j=1

where M=N i = 1mi is the total nuclear mass. Under the linear transformation (Ri )“ (ri ) the kinetic energy of the center of mass separates from the total kinetic energy. What remains is the kinetic energy of the internal motions about the center of mass 1 N − 1 N − 1 1 #2 T. = − ' 2 % % 2 i = 1 j = 1 mij #ri #rj with the effective masses given by

(12)

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

608 N

1 AikAjk = % mij k = 1 mk

(13)

The choice of the vectors (ri ) is an intermediate step in the construction of 3N −6 curvilinear coordinates describing the internal dynamics of the molecule. For the ammonia-like molecules, natural are the internuclear vectors ri along the chemical bonds formed between the central atom (i=0) and the peripheral atoms (i= 1, 2, 3). These nuclear bond vectors ri =Ri −R0 define a simple matrix A, which gives the following effective masses [74] 1 −1 1 + m− and mij =m0 m− ii =m i 0

for i "j. (14)

− i'



# #Rik # #r%ji # = − i' +% #rih (ria #Rik j #ria #r%ji

N−1

J. = − i' % ri × i=1

# . #ri

(15)

Now, let us consider the transformation LSA“ MSA. The LSA Cartesian coordinates rih =eh ·ri transform to the MSA coordinates r%ih =e%h ·ri according to rih = Rhir%ii,

(16)

where we apply Einstein’s summation convention only to the coordinate labels denoted by the Greek letters. The rotationally invariant set of the MSA Cartesian coordinates r%ih with h = x, y, z is a precursor of any curvilinear internal coordinate set {qw }. For example, the molecular bond lengths or bond angles are defined as scalar products of the ri vectors that are expressed by r%ih. The R matrix depends on the Euler angles r, €,  and its elements Rhi (r,€,) =eh ·e%i are direction cosines specifying the orientation of MSA relative to LSA. They can be treated as the rotational variables involved in the molecular wavefunction [63,64]. Using this fact, we can apply the chain rule to express the linear momentum operators pˆih = − i'#/#ria by the coordinates r%ia and Rhi :

(17)

After evaluating the derivatives #r%ji /#rih we come to −i'

# = #rih − i'



#Rik # # #R # +Rih + % ki rjk #rih #Rik #r%ii j (rih #r%ji



(18) Now, let us consider the LSA components J. h = eh ·J. of the angular momentum of Eq. (15) J. h = − i'mhik %rii i

The transformation A also separates the center of mass contribution to the total angular momentum. We shall ignore this contribution and consider only the angular momentum associated with the N −1 vectors



# , #rik

(19)

where mhik is the antisymmetric Levi–Civita tensor. Substituting Eq. (18) into Eq. (19) we obtain J. h = (J. hRik )

N−1 # # + % (J. hr%ii ) #Rik i = 1 #r%ii

(20)

where J. h act within the brackets. The angular momentum operators J. h play the role of the generators of the molecular axes rotations. They can be expressed in terms of the derivatives over the Euler angles. It is clear that r%ih, as the scalar products e%h ·ri, are independent of the Euler angles since they are rotationally invariant, so (J. hr%ii )= 0. Using the standard formulae for J. h and Rhi (r,€,) in the parametrization of the Euler angles we have (J. hRik )= i'mhilRlk

(21)

This result allows us to express J. h by the elements Rhi of the rotation matrix J. h = − i'mhilRik

# #Rlk

(22)

Now, we can see that the representation of J. h through Rhi enables one to treat the differentiation over Rhi, treated as independent variables. This leads to the standard commutation relation [J. h,J. i ]= i'mhikJ. k that follows directly from Eq. (22). The angular momentum operator can be projected onto the MSA which gives the components

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

J. %i = e%i ·J. = e%i ·ehJ. h. They can also be expressed by Rhi : J. %i = − i'mhzlRhiRzk

# #Rlk

# #Rkz

(25)

(26)

The above definition of the angular momentum MSA components J%i helps us to see the physical meaning of the terms in Eq. (18) defining the linear momenta. First, we evaluate the derivatives #Rki /#rih from the definition Rki =ek ·e%i : #Rki #e% =ek · i. #rih #rih

(27)

Since the basis orts are orthonormal, i.e. e%i ·e%k = lik, they satisfy the equation #e%i #e% · e%k + e%i · k = 0. #rih #rih

(28)

It is well known that the solution to this equation can be written as the vector product of e%i and some vector d ih that characterizes the response of the rotating molecular axes to the variation of the LSA components of the nuclear vectors: #e%i =e%i × d ih #rih

(29)

# # = − i'Rkh #r%ik #rih



− i'd ihzmzli Rkl



# # + %r%jl , #Rki #r% j ji (32)

where the presence of the angular momentum components is apparent. Now, we define the intrinsic vibrational angular momentum operators by L. i = − i'mikh %r%jk j

# . #r%jh

(33)

They do not satisfy the usual commutation relations because they are defined by the MSA components r%ih that are mutually dependent. The transformation law for the linear momenta takes a clear final form −i'





# # = Rkh − i' − d ihi(J%. i − L. i ). #rih #r%ik

(34)

analogous to that in Ref. [64] where the transformation of vectors Ri –rCM has been considered. The above law is general and convenient in applications, since all complications inherent in using explicitly the Euler angles disappear. Additionally, no matrix inversions are involved in the derivation procedure. The main problem consists in determining the coefficients d ihi. They can be calculated using a systematic procedure, which needs only the MSA basis orts as known functions of the nuclear bond vectors. Under the transformation of the linear momenta, the operator T. of Eq. (12) takes a simple and physically meaningful form T. =T. V + T. R + T. C

This relationship applied to Eq. (27) gives #Rki =Rklmlizd ihz. #rih

− i'

(24)

The operators J%i fulfill the ‘anomalous’ commutation rules [J. %h,J. %i ]= − i'mhikJ. %k

(31)

With this result we transform Eq. (18) to the form

Surprisingly, we have not found this interesting equation in the accessible literature, so we derive it in Appendix A. It can also be directly proved using the explicit form of the R matrix. The relation (24) reduces Eq. (23) to J. %i = i'milzRkl

1 #Rik d ihz = mzlkRil . 2 #rih

(23)

This equation can be simplified using the relation mhikRilRkz = Rhimilz

609

(35)

where (30)

Using the orthogonality of the rotation matrix RikRih =lkh and the value of the convolution of the Levi– Civita tensor we invert Eq. (29) relative to d ihz and obtain

1 1 T. V = % % pˆ pˆ + % bhipˆ ih 2 i j mij ih jh i

(36)

represents a pure vibrational contribution to the kinetic energy, generated by the vibrational momenta pˆ ih = − i'#/#r%ih. The operator T. R given by

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

610

1 T. R = vhiN. hN. i +vhN. h 2

(37)

describes the energy of a rotor characterized by the effective reciprocal tensor of inertia vhi =% % i

j

1 i d dj mij kh ki

(38)

defined by the mutually orthogonal sets of the nuclear vectors, it is sufficient to derive the simpler operator T. A describing only the frame. The second part of Eq. (43) does not cause difficulties. This result is valid regardless of the internal coordinates chosen for the frame and the way the MSA embedded to the frame is defined.

This energy is created by ‘reduced’ angular momentum components N. h =J. %h −L. h.The operator 1 T. C = (M. hN. h + N. hM. h ) 2

4. Construction of the MSA basis (39)

describes the Coriolis interaction through the operators M. h = − % % i

j

1 R d j pˆ mij ik ih ik

(40)

The linear terms 1 1 j bhi = % (N. kd ik )Rih 2 j mij

e%3 = l1r1 × r2 + l2r1 × r3 + l3r2 × r3 (41)

and 1 1 j j vh = % % [d ikz(N. zd kh ) −Rik (pˆ ikd ih )] 2 i j mij

(42)

have pure quantum mechanical origin, because they occur due to noncommuting operators. The operators T. R and T. C contain also the vibrational contributions from the operators L. h. The derivation of the rovibrational operator T. for a given molecule can be substantially simplified. Notice, that if some nuclear vector ri is not involved in the construction of the MSA, then the corresponding tensor d i is zero. Let us consider the set A of the vectors ri, called the frame, which are involved in the definition of the MSA, i.e. i A. The set B of the remaining vectors is called the top. If these two sets are orthogonal in the 1 −1 metric m − ij , i.e. m ij =0 if i A and j  B, then the rovibrational operator T. is decomposed into two parts T. = T. A (pˆAM. A,N. )+

1 1 % % pˆ pˆ 2 j  B k  B mjk jh kh

To derive the vectors d ih we have to define the basis of the MSA. As the vectors ri of the chemical bonds determine the molecular geometry in a natural way, we use them to construct the basis vectors. The unit vector e%3 along the z% axis can be chosen according to

(43)

where the Coliolis operators M. A h are defined only through the momenta pˆ i from the set A. Thus, in order to derive the operator T. for a complex molecule, composed of the frame and the top

This can be done even for the planar molecular configuration since for three coplanar vectors the above definition determines the direction along the z% axis in a unique way. This definition is correct also for two collinear vectors, in contrast to the bases used for the valence coordinates. In our case the singularity appears only for three collinear vectors but such linear configurations are inaccessible in the bound states of typical ammonia-like molecules. To find the coefficients li we use the fact that the z% axis makes the same angles with the three vectors ri, i.e. it makes trisection of the pyramidal angle with the vertex on the central atom e%3·r1 = r1 cos h= l3r1·(r2 × r3), e%3·r2 = r2 cos h= l2r2·(r1 × r3)= − l2r1·(r2 × r3), e%3·r3 = r3 cos h= l1r3·(r1 × r2)= l1r1·(r2 × r3). From the product r1·(r2 × r3)= r1r2r3 J cot h we obtain the coefficients li : l1 =

sin h , r1r2 J

l2 = −

sin h , r1r3 J

Thus, the basis ort e%3 is e%3 =



l3 =



sin h r1 r2 r1 r3 r2 r3 × − × + × . J r1 r2 r1 r3 r2 r3

sin h . r2r3 J

(44)



J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

According to our definition of the MSA, the basis ort e%1 lies in the plane formed by r1 and e%3 so r1 = x%1e%1 + z%1e%3 with x%1 =r1·e%1 =r1 sin h and z%1 = r1·e%3 = r1 cos he%3. Now, we have e%1 =

1 r1 −cot he%3. sin h r1

(45)

Since e%2 can be defined as the vector product e%2 = e%3 × e%1, we can easily express e%2 as 1 r e%2 = e%3 × 1. sin h r1

a12 =



1 ~2 ~1 cot cot − 1 , 2r2 sin h 2 2

a13 = −



611



1 ~3 ~1 cot cot + 1 . 2r3 sin h 2 2

Summing up, the basis of the MSA orts is defined as follows 3

e%i = % aiiri

for i=1,2 and e%3 = (e%1 × e%2).

i=1

(46)

(50)

Substituting Eq. (44) for e%3 into Eq. (46) and calculating the two-fold vector products we obtain

The coefficients aii are structural elements of the theory, as they come into the main terms in the Hamiltonian. For further considerations it is convenient to introduce coefficients Sii related to aii

e%2 = −

1 [r (l r r cos q12 +l2r1r3 cos q13) r1 sin h 1 1 1 2

+ r2(l3r1r3 cos q13 −l1r 21)

where we use the scalar products r1·r2 = r1·r2 cos q12 and r1·r3 =r1r3 cos q13. The relations between the valence and the UL coordinates permit writing e%2 = a21r1 +a22r2 + a23r3,



(47)



~ 1 ~ a21 = − cot 2 + cot 3 , 2r1 sin h 2 2

 

 

3

3

i=1

i=1

% aii ri = % Sii = 0

a1i a2j − a2i a1j =

(48)

J%= 4 sin

Here the first index refers to the basis ort and the second one to the bond vector. To represent the basis ort e%1 in the analogous form, let us write it also as the vector product e%1 = e%2 × e%3. After substituting here e%2 from Eq. (47) and e%3 from Eq. (44) we have





mijk sin h , r i rj J

(53)

mijk , J%

(54)

(49)

~ ~ ~1 sin 2 sin 3. 2 2 2

(55)

To prove Eq. (52) it is enough to form the dot product of e%i and e%3. To prove Eq. (53) we substitute e%i into the definition of e%3 and compare the result with Eq. (44). Moreover, it can be shown that: a11 = cot

~ 2 + ~3 a21, 2

a13 = cot

~2 a . 2 23

e%1 = a11r1 +a12r2 + a13r3, 1 ~ ~ cot 2 cot 3 − 1 , 2r1 sin h 2 2

(52)

where we introduce the reduced Jacobian

1 ~ ~ a23 = cot 3 − cot 1 . 2r3 sin h 2 2

a11 = −

for i= 1, 2

and the specific ‘commutation’ relations

S1i S2j − S2i S1j =

~ 1 ~ a22 = cot 2 + cot 1 , 2r2 sin h 2 2

where

(51)

They satisfy the following relations

− r3(l2r 21 +l3r1r2 cos q12)]

where

Sii = aii ri sin h.

a12 = cot

~3 a , 2 22

The analogous relations between Sli and S2i will be used in further considerations.

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

612

5. Determination of the ” coefficients All components d ihi will be subsequently calculated from the definition (31). For the index i =1 corresponding to the component x we have d ih1 =





#Ri3 #Ri2 #Ri2 1 Ri2 −Ri3 = − Ri3 . 2 #rih #rih #rih

The same procedure applied for the component with i=2 d ih2 = Ri3 yields

The scalar product Ri2 =ei · e%2 =a21r1i + a22r2i +a23r3i can be written in terms of the structural elements S2j

#Ri1 , #rih



d ih2 = a1i Rh3 − Rhi =





r%iir%i3 r 2i



S1i r%ii Rh3 − Rhi cos h . ri sin h ri

(59)

3

Ri2 =

r 1 % S2j ji. rj sin h j = 1

The calculation of the component with i=3

The differentiation of Ri2 gives

 

 

#Ri2 # rji r # S2j 1 = % S + % ji . #rih sin h j = 1 2j #rih rj j = 1 rj #rih sin h 3

3

The property 3j = 1 S2j =0, see Eq. (52), allows us to write the derivatives of Ri2 as follows

       

3 #Ri2 1 # rji = % S2j #rih rj #rih sin h j = 1

 

r3i r1i # S23 − . r3 r1 #rih sin h

(56)

After multiplying this expression by Ri3 the last two terms disappear, as they are proportional to differences of the z% MSA components of the vectors rj /rj. Now one obtains d ih1 = −

 

Ri3 3 # rji % S2j . sin h j = 1 #rih rj

(57)

Since the Cartesian LSA components of the vectors ri are independent, the derivatives in Eq. (57) can be calculated by the standard way

 





# rji # rji l l = = ij lhi − 3ij rjirjh. #rih rj #rih rjkrjk rj rj Using this result one has d = −a2i i h1





r% r% Rh3 −Rhi ii2 i3 ri







(60)

is more difficult as in this case the contribution from the derivatives of Sii over the Cartesian coordinates does not vanish. The elements S22 and S23 which contribute to #Ri2/#rih, specified in Eq. (56), can be written as

'

1− cos q13 2(1− cos q12)(1− cos q23)

S23 =− sin h

'

and

1− cos q12 . 2(1− cos q13)(1− cos q23)

Now, we need derivatives of cos qjk over rih. Because cos qjk =

rj rk · , rj rk

it is easy to show that # cos qjk lhi = (l r + likrji ) #rih rj rk ij ki −





l rjirki lij r + ik r . rj rk r 2j jh r 2k kh

(61)

With the above results, after lengthy algebraic manipulations we obtain d ih3 = −

S2i r% =− Rh3 −Rhi cos h ii . ri sin h ri



#Ri2 #R #R 1 R − Ri2 i1 = Ri1 i2 2 i1 #rih #rih #rih

S22 = sin h

r r # S22 + 2i − 1i r2 r1 #rih sin h +

d ih3 =

(58)

li1Rh2 + cot hd ih1. r1 sin h

(62)

The analogous coefficients can be derived much easier for the valence coordinates and the MSA based on the two vectors r1 and r2. They are equal

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

d0 ih1 =

li1Rh2 , r1

d0 ih3 = −

d0 ih2 = −

li1Rh1 , r1

# cos q13 r%3h r%1h = − cos q13 2 , #r%1h r1r3 r1

li2Rh2 − cot q12d0 ih1, r2 sin q12

6. Vibrational momenta To express the kinetic energy operator in terms of the internal UL coordinates it is necessary to transform the momentum operators pˆ ih = − i'#/ #r%ih from the Cartesian MSA coordinates to the internal coordinates. Here, this transformation is executed in two steps. In the first step C“ V, we pass from the Cartesian MSA coordinates to the valence coordinates. After the second step V“ UL, using the connection between the valence and the UL coordinates, we finally arrive at the UL coordinates. The valence coordinates are constructed in terms of the MSA Cartesian coordinates treated as independent variables. In such a case, the transformation can be simply performed because we can use the chain rule

The calculation of #~/#r%ih is more difficult. We use Eq. (3) where the angle ~ is related to the angles qij. After differentiation of both sides of this equation over r%ih the term sin ~ sin q12 sin q13 appears on its left-hand side. We can represent this term by the internal coordinates as follows. Let us differentiate Eq. (3) over ~ under fixed q12 and q13. As a result we obtain sin ~ sin q12 sin q13 = − =−

# # cos q13 #~ # . + #r%ih # cos q13 #r%ih #~

# cos q23 #~ # (cos2 h+ sin2 h sin ~1). #~

Since #h 1 = c− 31 , #~

#~3 1 = c− 32 , #~

− J cot h

#~ = (cos ~ cot q12 sin q13 − cos q13) #r%ih # cos q12 + (cos ~ cot q13 sin q12 #r%ih − cos q12)

(63)

#~2 1 = c− 33 , #~

we have sin ~ sin q12 sin q13 = J cot h. Using this result we obtain the equation

3 # #rj # # cos q12 # =% + #r%ih j = 1 #r%ih #rj #r%ih # cos q12

# cos q13 # cos q23 + . #r%ih #r%ih

This equation can be rewritten in terms of the UL coordinates #~ #r%ih ~ ~ ~ ~ sin 1 cos 3 + sin2 h sin 1 sin 2 2 2 2 2 # cos q12 = ~ ~ #r%ih sin 2 1− sin2 h sin2 2 2 2 ~1 ~2 ~1 ~3 sin cos − sin2 h sin sin 2 2 2 2 # cos q13 − ~3 ~ #r%ih sin 1− sin2 h sin2 3 2 2 # cos q23 + . #r%ih − J cot h

Since rj = r%jir%ji

# cos q13 = 0, #r%2h

# cos q13 r%1h r%3h = − cos q13 2 . #r%3h r1r3 r3

where the ‘tilde’ distinguishes this variant from the previous one for the UL coordinates. These coefficients allow us to write immediately the Hamiltonian of Ref. [24].

+

613

and

cos qjk =

r%ji r%ki rj rk

we can immediately apply Eq. (61) to #cos qjk /#r%ih by replacing the vector LSA components by the MSA components. As a result, we obtain # cos q12 r%2h r% = −cos q12 1h , #r%1h r1r2 r 21 # cos q12 r%1h r% = −cos q12 22h, #r%2h r1r2 r2

# cos q12 =0, #r%3h

















J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

614

Finally, we need the operators #/# cos qli and #/#~ which can be determined from Eq. (4) by inverting the matrix c. The elements of the matrix c − 1 are given in Eq. (8). Thus, we have all necessary quantities for determining the needed derivatives and after cumbersome algebraic manipulations we have

pˆ 3y = sin h sin ~3pˆ r 3 +

cos h sin ~1 pˆ h pˆ 1x = sin hpˆ r 1 + r1 J%

pˆ 3z = cos hpˆ r 3 −







cos2 h ~3 ~2 cot pˆ ~ 2 +cot pˆ ~ 3 , r1 sin h 2 2

+



sin h sin ~1 pˆ h r1 J%







cos h sin ~3 cos ~2 pˆ h r2 J%

 



~2 ~1 pˆ , +cot 2 2 ~3

pˆ 2y =sin hsin ~2pˆ r 2 −

cos h sin ~3 sin ~2 pˆ h J% r2



cos ~2 cos h ~1 + pˆ + sin ~2 cot pˆ ~ 2 r2 sin h ~ 2 r2 sin h 2 2



+ cot

 

~2 ~ + cot 1 pˆ ~ 3 , 2 2







 

cos h ~ ~ ~ cot 1 pˆ ~ 2 + cot 2 + cot 1 pˆ ~ 3 , r2 2 2 2

pˆ 3x = sin h cos ~3pˆ r 3 +

cos h sin ~2 cos ~3 pˆ h r3 J%



sin ~3 cos2 h ~ pˆ ~ 3 − cos ~3 cot 1 pˆ ~ 3 r3 sin h r3 sin h 2



− cot

 

~3 ~ −cot 1 pˆ ~ 2 , 2 2



sin h sin ~2 pˆ h r3 J%



 

cos h ~1 ~3 ~1 cot pˆ ~ 3 − cot − cot pˆ ~ 2 , r3 2 2 2

The total vibrational angular momentum contains the contributions from individual angular momenta created by the varying bond vectors ri. The definition of the corresponding individual operators is compatible with Eq. (33): l. ji = mikhr%jkpˆ jh

(65)

Substituting the MSA Cartesian coordinates and the linear momentum operators expressed by the internal coordinates into the definition (65) we subsequently obtain l. 1x = − i' cot h l. 1y = − i'







# # + , #~2 #~3

l. 1z = i'

(66)

sin ~1 # J% #h



−cot h cot

sin h sin ~3 pˆ 2z =cos hpˆ r 2 + pˆ r2 J% h −

 

~3 ~1 − cot pˆ , 2 2 ~2

7. Vibrational angular momenta

sin ~2 cos2 h ~ pˆ ~ 2 + cos ~2 cot 1 pˆ ~ 2 r2 sin h r2 sin h 2

+ cot



where pˆ q − i'#/#q for q= ri, h and ~i, and J% is the reduced Jacobian of Eq. (55).

cos h ~ ~ cot 3 pˆ ~ 2 +cot 2 pˆ ~ 3 , r1 2 2

pˆ 2x = sin h cos ~2pˆ r 2 −

+



cos ~3 cos2 h ~1 pˆ ~ 3 − sin ~3 cot pˆ ~ 3 r3 sin h r3 sin h 2

− cot

(64)

1 pˆ 1y = − (pˆ + pˆ ~ 3), r1 sin h ~ 2 pˆ 1z =cos hpˆ r 1 −

+

cos h sin ~2 sin ~3 pˆ h r3 J%



~3 # ~ # + cot 2 2 #~2 2 #~3



n

,

# # + , #~2 #~3

l. 2x = −i'

!

sin ~2 sin ~3 # # − cot h cos ~2 #~2 J% #h



− cot h sin ~2 cot



+ cot

~1 # 2 #~2

 n"

~2 ~ # + cot 1 2 2 #~3

,

!

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

l. 2y = −i' −

cos ~2 sin ~3 # # − cot h sin ~2 #~2 J% #h



+cot h cos ~2 cot



+ cot l. 2z = − i'

~1 # 2 #~2

pˆ 1y =

 n"

~2 ~1 # + cot 2 2 #~3

!

sin ~2 sin ~3 # # −cot h cos ~3 J% #h #~3



+ cot h sin ~3 cot



− cot

~1 # 2 #~3

!

,



−cot h cos ~3 cot



− cot

pˆ 2y = sin h sin ~2pˆ r 2 −

cos h sin h l. + cos ~2l. 2z, r2 2x r2

sin h (sin ~2l. 2x − cos ~2l. 2y ), r2



cos h sin ~2 cos t3 sin t3 pˆ a − pˆ r r3 J% r3 sin a 3

 n"

pˆ 3z = cos hpˆ r 3 +

,



 

t cos2 h t t cos t3 cot 1 pˆ ~ 3 − cot 3 − cot 1 pˆ ~ 2 , r3 sin a 2 2 2

pˆ 3y = sin h sin ~3pˆ r 3 −

cos h sin h l. 3x + cos ~3l. 3z, r3 r3

sin h (sin ~3l. 3x − cos ~3l. 3y ) r3

The above formulae can be written in an extremely compact form

# . l. 3z = − i' #~3 Now, we can see the extraordinary character of the trisection of the pyramidal angle formed by the vectors of the three chemical bonds. From the above formulas it immediately follows that L. x =L. y =L. z 0.

cos h sin h l. 2y − sin ~2l. 2z, r2 r2



~1 # 2 #~3

~ # ~3 − cot 1 2 2 #~2

sin h l. , r1 1y

pˆ 3x = sin h cos ~3pˆ r 3 +

sin ~2 cos ~3 # # −i' −cot h sin ~3 #~3 J% #h

(68)

pˆ 2x = sin hcos ~2pˆ r 2 +

pˆ 2z = cos hpˆ r 2 +

 n"

~3 ~1 # − cot 2 2 #~2

cos h l. 1y, r1

l. 1z l. 1x =− , r1 sin h r1 cos h

pˆ 1z = cos hpˆ r 1 −

,

# , #~2

l. 3x = − i' −

l. 3y =

pˆ 1x = sin hpˆ r 1 +

615

(67)

Thus, the physical meaning of the trisection is that the MSA adjusts immediately to the vibrational nuclear motions in such a way that the components of the vibrational angular momenta of the chemical bond vectors compensate each other. The resulting components of the total vibrational angular momentum are exactly zero. This property reduces the significant part of the Coriolis interactions generated by the vibrational motion. Moreover, this simplifies considerably the rovibrational Hamiltonian. The vibrational angular momenta permit the following representation of the linear vibrational momenta

pˆ ih =

r%ih m r% l. pˆ ri − hik 2ii ik. ri ri

(69)

Obviously, the vibrational angular momenta are constrained not only by 3i = 1l. ih 0. They satisfy also the equations r%ihl. ih 0, which allow us to write the square of the linear momenta in an elegant form typical for the spherical coordinates l. 2 pˆ 2i = pˆ 2ri + i2. ri

(70)

8. The operators of the rotational and Coriolis interaction energy The rotational energy of Eq. (37) is determined by the symmetric v tensor. Its components are defined in Eq. (38) only by the coefficients d ihi,

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

616

derived in Section 5. As a result, we can represent all components vhi in terms of the structural elements Sih S 221 S2 S2 + 22 2 + 23 2 2 m11r 1 m22r 2 m33r 3

vxx =

+



2 S21S22 f(~2) S21S23 f(~3) + r1r2 r1r3 m0



S22S23 f(~1) + , r2r3

(71)

vxy =



3

+ +



G rx1 =



n

f(~1) (S S +S12S23), m0r2r3 22 13



G rz1 = cot hG rx1,

 

G rx2 =

vyz =cot h vxy +

sin ~2 sin ~3 S12 + S13 , m0r1r2 m0r1r3

G ry2 =



 

sin h cos h 1 −cos ~2 1− cos ~3 S12 + S13 , m0 r2 r3



sin ~2 sin ~3 S22 − S , m0r1r2 m0r1r3 23

+ cot2 h vxx −



sin h cos h 1− cos ~2 1− cos ~3 S22 + S23 , m0 r2 r3

G ry1 =

f(~3) (S11S23 +S21S13) m0r1r3

1 2 sin hm11r 21

(73)

(74)

f(~2) (S12S21 +S11S22) m0r1r2

vxz =cot h vxx −

vzz =

1 j Rikd ih s wik. mij

j Substitution of d ih and of s wik into this expression yields the following components G wh:

S12S13 f(~1) , r2r3

 

3

i=1 j=1

S11S21 S12S22 S13S23 + + m11r 21 m22r 22 m33r 23

+

(72)

w

G wh = − % %

2 S11S12 f(~2) S11S13 f(~3) + + m0 r1r2 r 1r3 +

T. C = % (G whpˆ w + pˆ wG wh)J. %h,

where the elements G wh of the general rovibration G matrix, have the form

S 211 S2 S2 + 12 2 + 13 2 vyy = 2 m11r 1 m22r 2 m33r 3



internal coordinates in the following order r1, r2, r3, h, ~2, ~3. The explicit form of s wih can be directly found from Eq. (64). The coefficients s wih are the generalized s vectors components [75]. Substitution of pˆ ik into M. h gives M. h = G whpˆ w which leads to the standard form of the Coriolis interaction energy operator







where f(~i ) =sin2 h +cos2 h cos ~i. The linear terms vh are reported in Appendix B. The components M. h of the Coriolis interaction operator are defined in Eq. (40). They can be easily expressed by the internal coordinates if the Cartesian components of the momentum operators pˆ ih are written in the form pˆ ih =s wihpˆ w, where w labels the



sin h cos h 1− cos ~2 1−cos ~1 S11 + S13 , m0 r1 r3

− 2sin ~2 2sin ~3 S22 − S23 . m0r1r2 m0r1r3



sin h cos h 1− cos ~2 1− cos ~1 S21 + S23 , m0 r1 r3

G rz2 =

sin ~2 + cot hG rx2, m0r1

G rx3 =

sin h cos h 1− cos ~3 1− cos ~1 S21 + S22 , m0 r1 r2



G ry3 =







sin h cos h 1−cos ~3 1− cos ~1 S11 + S12 , m0 r1 r2



G rz3 =

sin ~3 + cot hG rx3. m0r1

G hx = −v12 − +

617

1 sin ~2 sin ~3 sin ~1 = − + , I0 m0r1r2 m0r1r3 m0r2r3

(75)

1 f(~2) S21 + S22 m11r 21 m0r1r2



f(~3) S , m0r1r3 23



G hy = − v22 + +



J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

1 1 f(~2) f(~3) f(~1) =− + + − , I1 m11r 21 m0r1r2 m0r1r3 m0r2r3

1 f(~2) S + S m11r 21 11 m0r1r2 12

1 1 f(~2) f(~3) f(~1) = − + − , I2 m22r 22 m0r1r2 m0r1r3 m0r2r3



f(~3) S , m0r1r3 13



G hz =cot h G hx − G ~x2 = − cot h





G ~y2 =cot h



+ I

−1 23





sin ~2 sin ~3 S12 − S13 , m0r1r2 m0r1r3

1 I− cot 3

1 + I− 23 cot



1 1 1 f(~2) = + −2 , m0r1r2 I12 m11r 21 m22r 22



~1 ~3 1 − I− S 13 cot 2 2 21

n

~1 ~ 1 1 cot 3 S22 +I − −I − 3 0 S23 , 2 2

1 I− cot 3



~1 ~3 1 − I− S 13 cot 2 2 11



n

~ ~ 1 1 cot 1 −I − cot 3 S12 +I − 3 0 S13 , 2 2

1 cos ~2 + sin2 hm11r 21 m0r1r2

G ~z2 = −

1 1 f(~2) f(~3) f(~1) = + − − , I3 m33r 23 m0r1r2 m0r1r3 m0r2r3



 

cot2 h sin ~2 sin ~3 S22 + S23 +cot hG ~x2, m0 r1r2 r1r3

+

G ~x3 =cot h





1 I− cot 2

1 + I− 23 cot ~3 y

G = − cot h









~1 ~2 1 + I− S 12 cot 2 2 11



n

~1 ~ 1 1 +I − cot 2 S13 +I − 2 0 S12 , 2 2

1 cos ~3 + sin2 hm11r 21 m0r1r3



n

~1 ~ 1 1 cot 2 S23 +I − +I − 2 0 S22 , 2 2

1 I− cot 2

1 + I− 23 cot

G ~z3 = −

~1 ~2 1 S21 + I− 12 cot 2 2

sin ~2 sin ~3 + cot2 h S22 + S m0r1r2 m0r1r3 23



+ cot hG ~x3, where we introduce the effective moments of inertia defined by

1 1 f(~3) 1 = + −2 , m0r1r3 I13 m11r 21 m33r 23 1 1 1 f(~1) = + −2 . I23 m22r 22 m33r 23 m0r2r3

9. Vibrational G matrix The vibrational operator T. V given in Eq. (36) contains only the linear momenta pˆ ih. It can be converted by the substitution pˆ ih = s wihpˆ w to the standard form





1 T. V = % %G vwpˆ wpˆ v + h vpˆ v , 2 v w

(76)

where the elements of the matrix G are expressed by the s vectors 3

3

1 v w s ihs jh i = 1 j = 1 mij

G vw = % %

(77)

The coefficients h v standing at the linear terms pˆ v will be considered in Appendix C. The known s vectors yield the following expressions for G vw cos2 h+ sin2 h cos ~2 , m0

G r1,r1 =

1 , m11

G r1,r3 =

cos2 h+ sin2 h cos ~3 , m0

G r2,r2 =

1 , m22

G r3,r3 =

1 , m33

G r1,r2 =

G r2,r3 =

cos2 h+ sin2 h cos ~1 , m0

(78)



618

G r1,h =

sin h cos h 1 −cos ~2 1 − cos ~3 S12 + S13 m0 r2 r3



, G

r2,h

+

G r3,h =

G r1,~2 =



sin ~2 , m0r2



− 2 cos2 h 1 −cos ~2 1 −cos ~3 S22 + S23 r2 r3 m0 sin ~3 , m0r3

− sin2 h G r2,~2 =



− 2 cos2 h 1 −cos ~2 1 −cos ~1 S21 + S23 m0 r1 r3 − (1+ cos2 h)

 

sin ~2 , m 0r1

 



G h,h =

1



2 11 1

m r −2

G h,~2 =

 

sin ~3 , m 0r1

!

2 11 1

m r

S11 +

n

f(~2) f(~3) S12 + S , m0r1r2 m0r1r3 13

cot h sin ~1 ~ ~ sin ~3 − cot 3 − cot 1 J% m11r 21 2 m22r 22 2 sin ~2 ~3 ~1 + cot − cot m33r 23 2 2





G h,~3 =

!

cot h sin ~1 ~ − cot 2 J% m11r 21 2





~3 2

n



   

~1 2

n"

,





sin ~3 ~ sin ~2 ~ ~ cot 2 + cot 1 − cot 1 2 2 m22r 2 2 2 m33r 3 2

+

1 ~ sin ~3 sin2 h sin ~2 + cot 2 m0r1r2 2

+ f(~2) 1+ cos ~1 + sin ~1 cot −

~2 2

n



1 [sin ~3(sin ~3 + J%) m0r1r3

+ f(~3)(2+ cos ~2 + cos ~1)]

 



1 ~ sin ~3 sin2 h sin ~1 + cot 1 m0r2r3 2



1 G ~2,~2 = 2 sin h

>

1+ cos2 h cot2

~3 2

m11r 21 ~ 1+cos2 h cot2 1 2 + m22r 22 ~ ~ cos2 h cot 3 − cot 1 2 2 + 2 m33r 3



+vyy 1

 

+ f(~1) 1+ cos ~2 + sin ~2 cot

− 2cos2 h 1 −cos ~3 1 −cos ~1 S21 + S22 m0 r1 r2 − (1+ cos2 h)



1 ~1 sin ~2 sin2 h sin ~1 + cot m0r2r3 2

− f(~1) 1+ cos ~3 − sin ~3 cot

+

− sin2 h sin ~3 sin ~1 = − , r1 r2 m0

G r3,~3 =



− sin h sin ~2 sin ~1 + , m0 r1 r3 2

G r2,~3 =



− 2 cos2 h 1 −cos ~2 1 −cos ~3 S22 + S23 m0 r2 r3 − sin2 h

G r1,~3 =

+

1−cos ~1 S12 , r2



 

+ f(~3) 1+ cos ~1 − sin ~1 cot

sin h cos h 1 −cos ~3 (S11 −1) m0 r1 +

r3,~2



1−cos ~1 S13 , r3



1 [sin ~2(sin ~2 − J%) m0r1r2

+ f(~2)(2+ cos ~3 + cos ~1)] 1 ~3 − sin ~2 sin2 h sin ~3 + cot m0r1r3 2

sin h cos h 1 −cos ~2 = (S11 −1) m0 r1 +

G



J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

+

~1 2

n"



2



2 sin2 h(1− cos ~2 − f(~2)) m0r1r2

+cos2 h(1− cos ~2 + f(~2)) cot



n

~2 ~ ~ cot 3 − cot 1 2 2 2

,



2 cos2 h(1− cos ~3 +f(~3)) m0r1r3

cot +



J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628





n

~3 ~3 ~1 cot − cot 2 2 2



cot



~1 ~ ~ cot 3 − cot 1 2 2 2

G ~2,~3 =

1 sin2 h

>

+

n"

1+ cos2 h cot

~2 ~3 cot 2 2



 



~ ~2 cot 1 2 2



1 sin2 h(1− cos ~3 −f(~3)) m0r1r3



+2 cos2 h(1−cos ~1 +f(~1))cot2

>

1 + cos2 h cot2



+

+



~2 ~ + cot 1 2 2 m22r 22

1 + cos2 h cot2 m33r 23

~1 2

n"

~1 2

,





n

~3 ~ ~ cot 2 + cot 1 2 2 2

2 cos2 h(1 −cos ~1 + f(~1)) m0r2r3



~1 ~ ~ cot 2 +cot 1 2 2 2

n"

.

In order to check the correctness of our Hamiltonian we compared the derived vibrational (G vw), Coriolis (G wh), and rotational (G hi vhi ) components of the tensor G with those derived by using the Eckart approach founded on the classical kinetic energy [27]. This energy in the classical Lagrangian form is a quadratic form of the vibrational and angular velocities, defined by the matrix K. In the Eckart method, the tensor G is obtained by inverting K. We determined and inverted the matrix K by employing a numerical algorithm of Ref. [68]. We have found a full agreement between the values of K − 1 and their analytical counterparts in G. Our vibrational matrix G 6v agrees with that given by Handy et al. [22], except for their g 54, g 65 and g 66 reported with misprints. In g 54, s2s1 cot2 i in the first term should be replaced by − s2c1 cot2 i. In g 65, (2c2c12 + s1s3) in the sixth term should be replaced by (− 2c2c12 + s1s3). In g 66, 4vr 22 in the third term should be replaced by 4vr 21 and − c2c3 in the sixth term should be replaced by c2c3.

~2 2

m11r 21

cos2 h cot

n

n

~ ~3 cot 1 2 2

1 sin2 h(1− cos ~1 −f(~1)) m0r2r3

1 G ~3,~3 = 2 sin h



n

2 sin2 h(1 −cos ~3 − f(~3)) m0r1r3

cot

1 sin2 h(1− cos ~2 −f(~2)) m0r1r2

+2 cos2 h(1−cos ~3 +f(~3))cot





m11r 21 ~ ~ ~ cos2 h cot 1 cot 2 + cot 1 2 2 2 + 2 m22r 2 ~ ~ ~ cos2 h cot 1 cot 3 − cot 1 2 2 2 − m33r 23

−2 cos2 h(1−cos ~2 +f(~2)) cot +



~ ~2 ~ cot 2 +cot 1 2 2 2

+ cos2 h(1 −cos ~3 + f(~3)) cot

,



+

2 cos2 h(1 −cos ~2 + f(~2)) m0r1r2

cot

2 cos2 h(1− cos ~1 +f(~1)) m0r2r3



619

2

10. Symmetrized Hamiltonian in the Radau coordinates The bond vectors used in the definition of the UL coordinates create the coupling terms in the kinetic energy operator of Eq. (12). To eliminate such terms it is desirable to introduce some orthogonal nuclear vectors diagonalizing this opera-

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

620

tor. The transformation from the nonorthogonal vectors to various systems of the orthogonal nuclear vectors has been widely discussed in the literature [76–80]. For our purposes the most convenient are the orthogonal vectors that are geometrically and physically close to the bond vectors. So we employ the Radau coordinates for the following reasons. (i) They imitate the usual bond vectors because the canonical point of the Radau coordinates is placed close to the central atom in ammonia-like molecules [24]. Thus, they are well adapted to the description of the potential energy function. (ii) They are symmetrically equivalent under permutations of identical terminal nuclei, so the symmetry of the molecule can be exploited easily. (iii) They conserve our choice of the MSA and of the internal coordinates. They adequately describe the inversion motion in a natural way because the canonical point for the planar configuration of the molecule also lies in the molecular plane. The transformation from the bond vectors to the Radau ones can be formally realized by the transition to the limit m0 “ . In fact, the effective masses in the kinetic energy operator are 1 −1 1 1 −1 m− + m− and m − for the bond ii =m i 0 ij =m 0 −1 vectors, and mii =mi and m ij =0 for the Radau vectors. Obviously, for m0 “ we reach the orthogonal representation for the kinetic energy operator. As a result we can write T. Rd =





1 3 1 2( −i') 1 % pˆ 2ri + pˆ ri + G hhpˆ 2h 2 i = 1 mii ri 2 1 1 + G ~2~2pˆ 2~ 2 + G ~3~3pˆ 2~ 3 +G ~2~3pˆ ~ 2pˆ ~ 3 2 2 +G h~2pˆ hpˆ ~ 2 +G h~3pˆ hpˆ ~ 3 +h hpˆ h +h ~2pˆ ~ 2 1 +h ~3pˆ ~ 3 +M. hJ. h + vhiJ. hJ. i +vhJ. h, 2

(79)

where the G matrix and the other coefficients at the operators can be easily obtained from those 1 given in Sections 8 and 9 by taking m − 0 =0. The vibrational part of the operator T. Rd can be written in an extremely compact and elegant form if

we take into account the definition of the vibrational angular momenta: T. Rd =





1 3 1 − 2i' l. 2i % pˆ 2ri + pˆ ri + 2 + h hpˆ h ri 2 i = 1 mii ri + h ~2pˆ ~ 2 + h ~3pˆ ~ 3.

It is important that the matrix elements of T. Rd can be calculated without any problems although some terms in T. Rd proportional to (sin(~3 − ~2)/ 2) − 2 are unfactorizable. To avoid a direct evaluation of the two-dimensional integrals containing such terms it is enough to take the variational basis functions in the form C= sin

~1 ~ ~ sin 2 sin 3b. 2 2 2

(80)

Under such a choice the matrix elements are fully factorized if b is a product of one-dimensional functions depending only on a single coordinate. In addition, these elements do not contain any essential singularities as they are cancelled by the Jacobian and the weight factor in Eq. (80). Let us consider a part of the function b depending on ~2 and ~3 in the product form b(~2,~3)= „2(~2)„3(~3). For the symmetric ammonia-like molecule of C36 symmetry point group the functions „2(~2) and „3(~3) should be localized around the points ~2e = 2/3p and ~3e = 4/3p, correspondingly. A more convenient would be the set of functions localized around the same point ~e = 2/3p. To achieve this let us pass from the angle ~3 to the angle ~%3 = 2y −~3 for which ~%3e = ~e. It is clear that b(~2,~3) cannot be explicitly symmetrized according to the full permutation symmetry group S3 containing the cyclic permutation (123), because of the lack of the function „1(~1) in the product. However, we can extend the basis set by including „1(~1) into the product b(~2,~%3)= „1(~1)„2(~2)„3(~%3) and treating the angle ~1 = ~3 − ~2 = 2p− (~%3 + ~2) as a variable on equal footing with ~2 and ~%3. This allows us to exploit fully the symmetry of ammonia-like molecules. The set of one-dimensional basis functions „1(~1) will also be localized about ~e = 2/3p. In this way only one basis set of one-dimensional functions can be employed for all three coordinates. This





 

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

significantly reduces the number of the matrix elements necessary to build the total Hamiltonian matrix. To execute the differential operators on the basis functions we have to transform the operator T. Rd to the form containing the derivatives over ~1, ~2, and ~%3. The part of the wave function depending on these coordinates can be written as

G h~1 =

Then, operating by #/#~2 on the function C1(~2,~3,~3 − ~2) is equivalent to operating by #/ #~2 −#/#~1 on the function c2(~2,~%3,~1), if ~1 is treated as an independent variable:



#c1(~2,~3,~3 − ~2) # #

− c (~ ,~% ,~ ). #~2 #~2 #~1 2 2 3 1 Similarly, for ~3 one has



After executing these operators on the basis functions we must return to the old variables ~2 and ~%3, and then calculate the necessary matrix elements as the corresponding integrals over ~2 and ~%3. The totally symmetrized form of the operator T. sRd is



1 3 1 −2i' 1 % pˆ 2ri + pˆ ri + G hhpˆ 2h 2 i = 1 mii ri 2 1 1 1 + G ~1~1pˆ 2~ 1 + G ~2~2pˆ 2~ 2 + G ~%3~%3pˆ 2~ %3 2 2 2 −G ~1~2pˆ ~ 1pˆ ~ 2 − G ~1~%3pˆ ~ 1pˆ ~ %3 −G ~2~%3pˆ ~ 2pˆ ~ %3 +G h~1pˆ hpˆ ~ 1 +G h~2pˆ hpˆ ~ 2 +G h~%3pˆ hpˆ ~ %3

(81)

+h hpˆ h +h ~1pˆ ~ 1 +h ~2pˆ ~ 2 +h ~%3pˆ ~ %3 1 +M. hJ. h + vhiJ. hJ. i +vhJ. h. 2 where the explicit forms of the G matrix elements for the angular coordinates are

(82)





,



cot h sin ~1 ~2 sin ~%3 ~2 ~1 cot − cot + cot J% m1r 21 2 m2r 22 2 2





sin ~2 ~ cot 1 , 2 m3r 3 2





~% 2 ~2 + cot 3 1 2 2 G ~1~1 = 2 sin h m1r 21 ~ ~% 1+cos2 h cot2 2 1+ cos2 h cot2 3 2 2 + + , m2r 22 m3r 23 cos2 h cot

pˆ ~ 3 “ −pˆ ~ %3 +pˆ ~ 1.





sin ~2 ~% ~ cot 3 +cot 1 2 m3r 3 2 2

<



Now it is clear that the transformation of T. Rd reduces to a simple replacement

T. sRd =



+

#C1(~2,~3,~3 − ~2) # #

− + C (~ ,~% ,~ ). #~3 #~%3 #~1 2 2 3 1

pˆ ~ 2 “pˆ ~ 2 −pˆ ~ 1,

G h~%3 =

sin ~%3 ~2 sin ~2 ~%3 cot + cot , 2 2 m22r 2 2 m33r 3 2

cot h sin ~1 ~%3 sin ~%3 ~1 cot + cot 2 2 J% m1r 1 2 m2r 2 2 −

C = C(~2,~3)=C1(~2,~3,~3 −~2) = C2(~2,~%3,~1).



cot h sin ~1 ~2 ~%3 − cot + cot J% m11r 21 2 2 +

G h~2 =

621

G ~2~2 =

1 sin2 h

<

1+ cos2 h cot2

~%3 2

m1r 21 ~1 1+cos2 h cot2 2 + 2 m2r 2 ~ ~% cos2 h cot 3 + cot 1 2 2 + m3r 23

1 G ~%3~%3 = 2 sin h

<

=



1+ cos2 h cot2



~2 2



m1r 21 ~2 ~1 cos2 h cot + cot 2 2 + m2r 22 ~1 1+cos2 h cot2 2 + , 2 m3r 3

=

2

2

,

=

622

>

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628





~%3 ~2 ~%3 cot + cot 2 2 2 2 m1r 1 ~ ~1 2 1− cos2 h cot cot 2 2 + m2r 22 ~%3 ~% ~ cot 3 + cot 1 cos2 h cot 2 2 2 − , 2 m3r 3 ~2 ~2 ~%3 cos2 h cot cot + cot 1 2 2 2 G ~1~%3 = 2 sin h m1r 21 ~2 ~ ~ cot 2 + cot 1 cos2 h cot 2 2 2 + 2 m 2r 2 ~ ~% 1− cos2 h cot 3 cot 1 2 2 + , m3r 23

1 G ~1~2 = 2 sin h

cos2 h cot

?



>





G ~2~%3 =

+

1 sin ~2 − sin ~1(1−cos ~2) m2r 22 s12

+

1 sin ~2 , m3r 23 s13





+

1 sin ~%3 m2r 22 s12

+

1 sin~%3 − sin ~1(1−cos ~%3) ,where sij s13 m3r 23

The Coriolis operators M. h are given by



M. x = −

~2 ~% cot 3 2 2







m1r 21 ~ ~ ~ cos2 h cot 1 cot 2 + cot 1 2 2 2 + 2 m 2r 2 ~ ~% ~ 1 cot 3 + cot 1 cos2 h cot 2 2 2 + m3r 23



?

,

hh

The element G is the same as for the original Radau coordinates. The linear terms h w are



−i' 1 sin ~1 2 m1r 21 s23

M. y =



1 sin ~1 − sin ~2(1 − cos ~1) m2r 22 s12

+

1 sin ~1 − sin ~%3(1 − cos ~1) , m3r 23 s13

n

S23 ~% cot 3 pˆ ~ 1 m3r 23 2 S22 S21 ~%3 ~1 cot − cot 2 2 m1r 1 2 m 2r 2 2

S23 ~ ~% cot 3 +cot 1 m3r 23 2 2 S21 ~ cot 2 2 m1r 1 2



pˆ ~ 2





S22 ~ S ~ ~ cot 2 + cot 1 − 232 cot 1 pˆ ~ %3, m2r 22 2 2 m 3r 3 2



+ cot h −









S11 ~ ~% cot 2 + cot 3 m1r 21 2 2



S12 ~ S ~% cot 2 + 132 cot 3 pˆ ~ 1 2 m 2r 2 2 m3r 3 2

+ cot h (83)

(84)

S11(1−S11) S 212 S2 − − 132 pˆ h 2 2 m1r 1 m2r 2 m3r 3

+

+



~ S21 ~% S ~ cot 3 +cot 2 − 222 cot 2 2 m 1r 1 2 2 m2r 2 2

+ cot h − +

h ~1 =

      

+ cot h − +



S21(1−S11) S12S22 S13S23 + + pˆ h m1r 21 m2r 22 m3r 23

+ cot h

>

n

= 4 sin2 ~i sin2 ~j.

1 sin2 h 1− cos2 h cot

n

− i' 1 sin ~%3 − sin ~2(1−cos ~3) 2 m1r 21 s23

h ~%3 =



?



− i' 1 sin ~2 − sin ~%3(1−cos ~2) 2 m1r 21 s23

h ~2 =



S11 ~% S ~ cot 3 + 122 cot 1 2 m1r 1 2 m 2r 2 2



S13 ~% ~ cot 3 +cot 1 pˆ ~ 2 m3r 23 2 2

+ cot h +







J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

S11 ~2 ~2 S12 ~1 cot − cot + cot m1r 21 2 m2r 22 2 2



S13 ~1 cot pˆ ~ %3, 2 m3r 3 2

M. z = −

1 (pˆ ~ −pˆ ~ %3) + cot hM. x. m1r 21 sin2 h 2

The symmetry of the expressions obtained is apparent, so the partitioning of the Hamiltonian matrix into the symmetry blocks becomes easy. The only disadvantage of the proposed symmetrization method is the nonorthogonality of the angular basis set in the form of a triple product. However, the orthogonalization can be performed routinely.

11. Conclusions In this paper a new exact quantum mechanical rovibrational kinetic energy operator T. for ammonia-like molecules is presented. To derive T. we use the method exploiting the idea of a rotational invariance, analogous to the concept of Lukka [57]. However, we do not employ his infinitesimal rotational coordinates. Instead, the angular momentum operators are expressed directly by the elements of the rotation matrix. This method is practical in applications since all complications, which appear when the Euler angles are explicitly used, disappear. There is no necessity of inverting the large matrix K to avoid the work with the Euler angles. The general structure of the presented rovibrational operator T. is simple and its various terms have a clear physical meaning. The simplicity enables derivation of a useful construction rule formulated in Eq. (43). This rule allows us to build T. from simpler operators derived from the smaller parts of a molecule, the frame and the top. It can be viewed as a recursion principle which substantially simplifies the derivation of T. . This principle is valid irrespective of what vibrational coordinates and MSA are defined for the frame. The only necessary condition is the mutual

623

orthogonality of the frame and the top nuclear vectors. This result is a generalization of the previously derived recursion formula considered for particular sets of the orthogonal frame and top vectors [34]. Recently, the importance of the orthogonality condition has been stressed and a similar result has been obtained for a special case when the MSA of the frame is based only on two nuclear vectors [81]. The usefulness of the recursion principle has been demonstrated previously [34,81] so, here we only note that using this principle we can easily build T. for a molecule containing the amino group NH2 executing the inversion motion. For NH2 we can use the UL coordinates and the standard valence coordinates for the rest of the molecule. The most difficult stage in the derivation method consists in the calculation of the tensor d i. However, it can be determined only from the MSA basis orts defined as functions of the nuclear vectors. For the most frequently used MSA embedded to two nuclear vectors, the tensor d i can be calculated very easily. The expression of T. is given in the three-vectors embedded MSA. This is a new kind of the body reference frame, which is strictly connected to the vibrational coordinates. In the definition of our MSA, the collective inversion angle h is involved. The MSA z% axis makes the same angle h with each bond vector. The collective character of the MSA and of the angle h makes the determination of the tensor d i nontrivial. The vibrational UL coordinates, for which the pure vibrational operator T. V has been derived previously [22], are well adapted to the description of the inversion motion and satisfy most of the optimum criteria for the internal coordinates. Additionally, these coordinates together with the MSA reveal an important and remarkable property. They cancel exactly the total vibrational angular momentum. This simplifies T. and substantially reduces most of the Coriolis terms in a large volume of the configuration space. Only the Coriolis interaction terms involving the angles ~2 and ~3 are not small. The most separable form of T. is achieved for the orthogonal Radau vectors that enable intro-

624

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

duction of a similar MSA and the internal coordinates. Although the operator T. in the Radau or the bond vectors representation contains unfactorizable angular terms, its matrix elements can be factorized under appropriate choice of the weight factor in the variational basis functions. The introduced weight factor also cancels all singular terms in Hamiltonian matrix elements. We proposed a method of explicit symmetrization of T. treating all torsional coordinates ~i on equal footings. This enables a full exploitation of the symmetry of ammonia-like molecules in the variational calculations. According to the advised requirement of Ref. [58] ‘‘Now, because there is so much work in progress using the secular equation method for tetra-atomic molecules, it is important that the operators used for such systems are correct’’, we proved our operator by a fully numerical procedure.

1 mi|~Rki = mi|~mz€imhlkRhzRl€ 2 1 = (l|zl~€ − l|€l~z )mhlkRhzRl€ 2 = mhlkRh|Rl~.

Appendix B The linear coefficients vh defined in Eq. (42) contain two components v Jh and v ph. The first term in Eq. (42) is the contribution v Jh from the operator N. z = Jz − Lz. Since L. z 0 one has 1 1 j v Jh = % % d ikz(J. zd kh ). 2 i j mij

It is convenient to rewrite v Jh in another form. Because vhi, given in Eq. (38), are independent of the Euler angles, after operating by J. z on vhi one obtains

Appendix A In order to prove the relation (24) we start from the expressions for basis vectors e%z × e%€ =mz€le%l.

(A1)

After convoluting the both sides of this equation with the tensor mz€i we find mz€ie%z ×e%€ =mz€imz€le%l =2e%i. Substitution of e%l = Rhleh into the left-hand side of the above equation yields 1 e%i = mz€iRhzRl€eh × el. 2 Dot multiplication of the both sides of this equation on ek and application of Eq. (A1) to the LSA basis vectors, leads to 1 Rki = mz€imhlkRhzRl€. 2 If we again convolute both sides of this equation with the tensor mi|~, then we can write

(B1)

%% i

j

1 i 1 j j d kh(J. zd ki )= − % % (J. zd ikh)d ki . mij i j mij

As a result, v Jh can be written in the form v Jh = −

1 3 3 1 % % (J. d i )d j . 2 i = 1 j = 1 mij z iz ih

Using the coefficients d iiz as functions of the structural elements and of the rotation matrix R, and acting by J. z on this matrix (J. zRhk )= i'mzikRhi, we can transform (J. zd iiz) to the form (J. zd iiz)



= i' − a2i Rh2 + a2im1kiRhk cos h − a1im2kiRhk cos h

r%ii − a1i Rh1 ri

r%ii li1Rh1 − ri ri sin h

+ cot ha2im3kiRhk cos h



r%ii . ri

Substitution of this expression into v Jh allows us to obtain consequently

v J1 =



J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

i' a2j r%j 1 cot h % 2 m j=1 j 1r1rj 3



3

(a1i a2, j −a1j a2i )r%j 1 , mij rj i=1 j=1

+sin2 h cos h % % v J2 =



3

3

(a1,i a2, j −a1, j a2,i )r%i2 mij ri i=1 j=1 3

v J3 =



(a1,i a2, j − a1, j a2,i )r%i1r%j 2 , mij ri rj i=1 j=1

− cot h cos2 h % %



3 i' a sin h − % 2, j 2 m j=1 j 1r1 3



 

i' cot h 2

v J2 = −

 



(B2)



cos ~3 1 i' cos2 h cot h + S13 + , m0r1r3 m0r1r2 2 m0r2r3

+ v J3 = −

i' 2

 



−cos h

r%jk 3 1 % (pˆ ln a2, j ) rj i = 1 mij rj ik



=−

!

"

,

 

1 3 3 1 r% % % pˆ a l − cos h jk 2 i = 1 j = 1 mij ik 1j k3 rj

n

1 3 % (a r ) 2 j = 1 1j j

3

%

−cos h

1 1 + S m11r12 m0r2r3 11

  n    n

1 (pˆ i3 ln a2j ) i = 1 mij rj

1 (pˆ i3 ln a1, j ) m i=1 ij rj

cos ~2 1 + S m0r1r2 m0r1r3 12

+

+

  n

cos ~3 sin2 h + S , m0r1r3 m0r1r2 23

+

+



cos ~2 sin2 h + S m0r1r2 m0r1r3 22

+

3

%

v p2 = −

1 sin2 h + S21 m11r12 m0r2r3

i' cot h 2

!

r% 1 pˆ ik cos h jk m r rj i=1 ij j

Employing the relations (53) for the structural elements aii we find

n

1 3 % (a r ) 2 j = 1 2j j

−%

(a a −a1, j a2,i )r%j 1 +cos2 h % % 1,i 2, j . mij rj i=1 j=1

v J1 =

 

1 3 3 1 r%jk % % pˆ ika2j lk3 − cos h 2 i = 1 j = 1 mij rj

3



3

(B3)

Their components are the following

=

+sin2 h cos h % % 3

1 1 v ph = − % % R (pˆ d j ). 2 i j mij ik i,k ih

v p1 =

3 a1j r%j 1 i' cot h % 2 j = 1 mj 1r1rj



625

3

1 cos2 h − S m11r12 m0r2r3 21

cos ~2 cos2 h − S m0r1r2 m0r1r3 22

cos ~3 cos2 h − S . m0r1r3 m0r1r2 23

Now, let us consider the contributions v ph to vh from the vibrational momenta

r%jk 3 1 % (pˆ ln a1j ) rj i = 1 mij rj ik



3

1 r%jk pˆ ik cos h m r rj i=1 ij j

−%

v p3 = cot hv p1 + +



"

,

1 1 3 1 % pˆ 2 i = 1 mi1 ik r1 sin h







1 3 3 1 r% % % a2j lk3 − cos h jk (pˆ ik cot h). rj 2 i = 1 j = 1 mij

We will not report detailed complex forms of these coefficients since they are not used themselves, but they enter vh only together with v Jh. We have shown only the derivation method, which can be useful in an attempt to use this approach for other Hamiltonians. We report here the final results for vh : v1 = − cot hv3 +





i' sin ~2 sin ~3 cot h S12 + S , 2 m0r1r2 m0r1r3 13

(B4)



626

+

sin ~2 (cos2 h −sin2 h −cos ~3) + m r J% m0r1r2 1

2 33 3

n

+

2f(~3) 1 sin ~3 − S13 m0r1r2 m0r2r3 J%

+

(sin ~2 sin ~3 − 1) 1 1 − J% m0r1r2 m0r1r3



v2 = i' cot h 2



!  1

1+

2 11 1

m r

"

+

+



sin ~1 2f(~2) sin ~3 − J% m0r1r2 J%

n

2f(~3) sin ~2 2(cos2 h −sin2 h) + S11 + m0r1r3 J% m0r2r3



+ −

sin ~3 2f(~2) + m r J% m0r1r2 1



n

+

n 

+

1 (sin t2 sin t3 −1) 1 − J% m0r1r2 m0r1r3

+

2sin2 a− cos2 a m0r2r3

+i'



"

1 sin ~3 S23 m0r2r3 J%

+

i' 1 1 1 − . 2 J% m0r1r2 m0r1r3

bjh =

1 w 1 % % s (pˆ s v )+ % bjhs vjh, 2 i = 1 j = 1 mij ih w j,h j=1 1 1 % R (J. d i ). 2 i = 1 mij ih z iz

m11r 21

2f(~3) sin ~2 cos h + S m0r1r3 J% m0r2r3 21



+ − +

 n

sin ~1 2f(~2) sin ~3 1+ − J% m0r1r2 J% 2

+

sin ~3 2f(~2) cos2 h + − m r J% m0r1r2 m0r1r3 1

2 22 2

n

(C2)



i' 3 1 l l % − a2ilh2 − a1ilh1 − i1 h1 2 i = 1 mij ri sin h

=

1 sin ~2 S22 + m0r2r3 J%

! 

  

1 i' 1 sin ~1 1+ sin2 h 1− 2 sin h m1j r1 J%

+

1 sin ~3 sin ~3 f(~2)− cos2 h 1+ m2j r2 J% J%



sin ~2 sin ~2 f(~3)+ cos2 h 1− m3j r3 J% J%



i' 1 1 cos ~3 − cos ~2 2 sin h m1j r1 J%

+

1 1− cos ~3 f(~2) m2j r2 J%



1− cos ~2 f(~3) , m3j r3 J%

1

 

! 

bj 2 =

1



"



r%il . ri

n n n"

bj 1 = −

~ ~ ~2 cot 3 cot 1 2 2 2 1 1 − cos ~2 + − − − J% m0r1r2 m0r1r3 m0r2r3 m33r 23

! 

(C1)

+( a2im1hl − a1im2hl + cot ha2im3hl ) cos h

cot h 1 1 − cos ~1 1 1 −cos ~3 + 2 J% m11r 1 J% m22r 22 J%

1



In the explicit form, they are

cot

i' v3 = − 2



Substituting the already derived d iiz into Eq. (C2) and using the orthogonality of the rotation matrix, we obtain the following general formula for bjh :

fcos2 h −sin2 a −cos t3) 2f(~3) + S m0r1r2 m0r2r3 13





Here we determine the coefficients at the linear Hamiltonian terms containing the vibrational momenta

bjh =

1 sin ~3 1 sin t2 − m0r2r3 J% m33r 23 J%

n "

2 33 3

where

1 sin ~2 + S12 m0r2r3 J% +

sin ~2 cos2 h 2f(~3) + − m r J% m0r1r2 m0r1r2 1

Appendix C

hv=

2 22 2

(cos2 h −sin2 h −cos ~2) m0r1r3

+



J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

,

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

627

Table 1 The coefficients H vij of the linear terms h v ij

h

~2

11

cot h[sin2 ~1+2(1−cos ~1)]J%−2

1 ~2 ~3 ~3 cot +cot +cot sin ~1 J%−1 2 2 2 2



~3





+cot









22

cot h[sin2 ~3+2(1−cos ~3)]J%−2

~1 1 ~2 ~1 cot −cot +cot sin ~1 J%−1 2 2 2 2

33

cot h[sin2 ~2+2(1−cos ~2)]J%−2

1 ~1 ~3 cot −cot sin ~2 J%−1 2 2 2



12

−2 sin ~2+f(~2)(cot

13

2 sin ~3+f(~3)(cot

23





~3 2

n

~2 +sin ~3 sin ~1) J%−1 2

n



bj 3 = − +

n

~1 +sin ~2 sin ~3 2

−2 sin ~1+f(~1) cot

! 

cot

J%−1

 



.

After differentiation of the coefficients s vjh and summation in Eq. (C1), we can finally write h v in the form hv= −





1 1 i' % H vii + % H vij 2 2 i mii r i i " j m0ri rj ri ii

ri ij

(C3)

where H =2/ri and H =0. The coefficients H for v= h, ~2, ~3 are collected in Table 1.



~3 2

sin2 h−cos2 h+cot2

[3] [4] [5] [6]

"

1 sin ~3 1 sin ~2 1+ + 1− m2j r2 m3j r3 J% J%



~1 ~3 sin2 h−cos2 h+cot2 2 2

−cot

i' 1 sin ~1 cos h 1− 2 m1j r1 J%





v ij

[7] [8] [9] [10]

[11] [12]

[13]

References [1] C.G.J. Jacobi, Vorlesungen u¨ ber Dynamik, Reimer, Berlin, 1866. [2] G.D. Carney, L.L. Sprandel, C.W. Kern, Adv. Chem. Phys. 37 (1978) 305.

[14] [15] [16]

~3 sin ~1 J%−1 2









1 ~1 ~2 cot +cot sin ~3 J%−1 2 2 2 1 ~3 ~1 cot +cot 2 2 2

−4f(~2)J%−1

+sin ~2 sin ~1) J%−1



1 ~2 ~3 cot +cot 2 2 2

~1 2

+cot

~1 sin ~2 J%−1 2

−cot

~1 ~2 sin2 h−cos2 h+cot2 2 2





4f(~3)J%−1



−cot



~2 2

sin2 h−cos2 h+cot2



~1 2

J. Tennyson, Comput. Rep. 4 (1986) 1. S. Carter, N.C. Handy, Comput. Phys. Rep. 5 (1986) 115. P. Jensen, J. Mol. Spectrosc. 128 (1988) 478. Z. Bacˇ ic´ , J.C. Light, Ann. Rev. Phys. Chem. 40 (1989) 469. N.C. Handy, Int. Rev. Phys. Chem. 8 (1989) 275. D.W. Schwenke, Comput. Phys. Commun. 70 (1992) 1. X. Chapuisat, C. Iung, Phys. Rev. A 45 (1992) 6217. P. Jensen, Calculation of molecular rotation – vibration energiers directly from the potential energy function, in: S. Wilson, G.H.F. Diercksen (Eds.), Methods in Computational Molecular Physics, Plenum, New York, 1992. M.J. Bramley, J.W. Tromp, T. Carrington Jr., C.G. Corey, J. Chem. Phys. 100 (1994) 6175. J. Tennyson, Variational calculations of rotation – vibration spectra, in: P. Jensen, P.R. Bunker (Eds.), Computational Molecular Spectroscopy, Wiley, New York, 2000. M.J. Bramley, N.C. Handy, J. Chem. Phys. 98 (1993) 1378. M.J. Bramley, T. Carrington Jr., J. Chem. Phys. 99 (1993) 8519. M.J. Bramley, S. Carter, N.C. Handy, I.M. Mills, J. Mol. Spectrosc. 157 (1993) 301. S. Carter, N. Pinnavaia, N.C. Handy, Chem. Phys. Lett. 240 (1995) 400.

628

J. Makarewicz, A. Skalozub / Spectrochimica Acta Part A 58 (2002) 601–628

[17] D.W. Schwenke, J. Chem. Phys. 100 (1996) 2867. [18] S. Carter, N.C. Handy, J. Demaison, Mol. Phys. 90 (1997) 729. [19] J. Koput, S. Carter, N.C. Handy, J. Phys. Chem. 102 (1998) 6325. [20] F. Gatti, C. Iung, M. Menou, Y. Justum, A. Nauts, X. Chapuisat, J. Chem. Phys. 108 (1998) 8804. [21] F. Gatti, C. Iung, M. Menou, X. Chapuisat, J. Chem. Phys. 108 (1998) 8821. [22] N.C. Handy, S. Carter, S.M. Colwell, Mol. Phys. 96 (1999) 477. [23] M. Mladenovic´ , J. Chem. Phys. 112 (2000) 1070. [24] M. Mladenovic´ , J. Chem. Phys. 112 (2000) 1082. [25] A. Viel, C. Leforestier, J. Chem. Phys. 112 (2000) 1212. [26] C. Leforestier, A. Viel, F. Gatti, C. Munoz, C. Iung, J. Chem. Phys. 114 (2001) 2099. [27] C. Eckart, Phys. Rev. 47 (1935) 552. [28] A. Sayvetz, J. Chem. Phys. 7 (1939) 383. [29] C.C. Lin, J.D. Swalen, Rev. Mod. Phys. 31 (1959) 841. [30] J.T. Hougen, P.R. Bunker, J.W.C. Johns, J. Mol. Spectrosc. 34 (1970) 136. [31] H.M. Pickett, J. Chem. Phys. 36 (1972) 1715. [32] V. Szalay, J. Mol. Spectrosc. 102 (1984) 13. [33] Y. Guan, C.R. Quade, J. Chem. Phys. 84 (1986) 5624. [34] J. Makarewicz, A. Bauder, Mol. Phys. 84 (1995) 853. [35] J. Makarewicz, W. Łodyga, Mol. Phys. 64 (1988) 899. [36] J. Makarewicz, J. Phys. B 21 (1988) 1803. [37] A. Ernesti, J.M. Huston, J. Chem. Phys. 101 (1994) 5431. [38] H. Wei, T. Carrington Jr., J. Chem. Phys. 107 (1997) 3674. [39] H. Wei, T. Carrington Jr., J. Chem. Phys. 107 (1997) 9493. [40] E.B. Wilson, J.C. Decius, P.C. Cross, Molecular Vibrations, McGraw-Hill, New York, 1955. [41] J.K.G. Watson, Mol. Phys. 15 (1968) 479. [42] E.A. Hylleraas, Z. Phys. 48 (1928) 469. [43] L.G. Bonner, Phys. Rev. 46 (1934) 458. [44] R.D. Bardo, M. Wolfsberg, J. Chem. Phys. 67 (1977) 591. [45] S. Carter, N.C. Handy, Mol. Phys. 47 (1982) 1435. [46] H. Ruder, Z. Naturforsch. A 23 (1968) 579. [47] J. Tennyson, B.T. Sutcliffe, J. Chem. Phys. 77 (1982) 4061. [48] B.T. Sutcliffe, J. Tennyson, Int. J. Quantum Chem. 39 (1991) 183. [49] J. Zuniga, A. Bastida, A. Raquena, J. Chem. Soc. Faraday Trans. 93 (1997) 1681. [50] W. Bu¨ ttner, H. Ruder, Z. Naturforsch. A 23 (1969) 1163. [51] N.C. Handy, Mol. Phys. 61 (1987) 207. [52] M.J. Bramley, W.H. Green, N.C. Handy, Mol. Phys. 73 (1991) 1183. [53] B.T. Sutcliffe, Current Aspects of Quantum Chemistry, Elsevier, New York, 1982. [54] A.G. Csa´ sza´ r, N.C. Handy, Mol. Phys. 86 (1995) 959. [55] A.G. Csa´ sza´ r, N.C. Handy, J. Chem. Phys. 102 (1995) 3962.

[56] [57] [58] [59] [60] [61] [62] [63] [64]

[65] [66] [67]

[68]

[69]

[70] [71] [72] [73] [74] [75]

[76] [77] [78] [79] [80] [81]

.

S.B. Rempe, R.O. Watts, J. Chem. Phys. 108 (1998) 10084. T.J. Lukka, J. Chem. Phys. 102 (1995) 3945. S.M. Cowell, N.C. Handy, Mol. Phys. 92 (1997) 317. B. Podolsky, Phys. Rev. 32 (1928) 812. J.C. Decius, J. Chem. Phys. 16 (1948) 1025. J. Makarewicz, A. Skalozub, Chem. Phys. Lett. 36 (1999) 352. J.H. Frederick, C. Woywod, J. Chem. Phys. 111 (1999) 7255. J.D. Louck, H.W. Galbraith, Rev. Mod. Phys. 48 (1976) 69. L.C. Biedenharn, J.D. Louck, Encyclopedia of Mathematics: Angular Momentum in Quantum Physics, Addison – Wesley, Reading, 1981. A. Skalozub, A.Ya. Tsaune, Optics Spectrosc. 50 (1981) 247. J.P. Leroy, R. Wallace, Chem. Phys. 118 (1987) 379. J. Makarewicz, A. Skalozub, Exact Quantum Mechanical Kinetic Energy Operator for Internal Motions of a Molecule, The Proceedings of the XVIth International Conference on High Resolution Molecular Spectroscopy, Prague, 2000, p. 143. J. Makarewicz, The calculation of rotation – vibration energies for molecules with large amplitude vibrations, in: P. Jensen, P.R. Bunker (Eds.), Computation Molecular Spectroscopy, Wiley, New York, 2000. D. Papousˇek, V. S& pirko, A new theoretical look at the inversion problem in molecules, in: M.J.S. Dewar (Ed.), Topics in Current Chemistry, vol. 68, Springer, Berlin, 1976. V. S& pirko, J.M.R. Stone, D. Papousˇek, J. Mol. Spectrosc. 60 (1976) 159. D. Papousˇek, M.R. Aliev, Molecular Vibrational Rotational Spectra, Academia, Prague, 1982. V. S& pirko, J. Mol. Spectrosc. 101 (1983) 30. V. S& pirko, W.P. Kraemer, J. Mol. Spectrosc. 133 (1989) 331. P.R. Bunker, P. Jensen, Molecular Symmetry and Spectroscopy, NRC Research Press, Ottawa, 1998. G.O. Sørensen, A new approach to the hamiltonian of nonrigid molecules, in: J.M.S. Dewar (Ed.), Topics in Current Chemistry, vol. 82, Springer, Berlin, 1978. J.O. Hirschfelder, J.S. Dahler, J. Chem. Phys. 42 (1956) 363. F.T. Smith, J. Math. Phys. 3 (1962) 735. Yu.F. Smirnov, K.V. Shitikova, Sov. J. Part. Nucl. 8 (1977) 344. V. Aquilanti, S. Cavalli, J. Chem. Phys. 85 (1986) 1355. G. Littlejonh, M. Reinsh, Phys. Rev. 52 (1995) 2035. X.G. Wang, T. Carrington Jr., J. Chem. Phys. 113 (2000) 7097.