Six-Dimensional Potential Energy Surface and Rovibrational Energies

Accurate calculation of anharmonic vibrational frequencies of medium sized molecules using local coupled cluster methods. Tomica Hrenar , Hans-Joachim...
1 downloads 0 Views 295KB Size
Download PDF
5520

J. Phys. Chem. A 2006, 110, 5520-5529

Six-Dimensional Potential Energy Surface and Rovibrational Energies of the HCCN Radical in the Ground Electronic State† Mirjana Mladenovic´ * and Peter Botschwina‡ Institut fu¨r Physikalische Chemie, UniVersita¨t Go¨ttingen, Tammannstrasse 6, D-37077 Go¨ttingen, Germany

Cristina Puzzarini| Dipartimento di Chimica “G. Ciamician”, UniVersita` di Bologna, Via Selmi 2, 40126 Bologna, Italy ReceiVed: NoVember 21, 2005; In Final Form: January 12, 2006

We report large-scale quantum mechanical calculations for the HCCN radical in its ground electronic state. A six-dimensional potential energy surface based on MR-ACPF/cc-pVQZ ab initio energy points is developed and adjusted to reproduce experimental findings for ν15 and ν1 of HCCN. Rovibrational energy levels of HCCN and DCCN are computed for total rotational angular momentum J ) 0-4 by making use of combined (functional + point wise) coordinate representations together with contraction schemes resulting from several diagonalization/truncation steps. The classical barrier to linearity is determined to be 287 cm-1. Spectroscopic parameters are calculated for low lying states and compared with available experimental data. Energy patterns attributed to the ν4 bending mode and to the quasilinear ν5 bending mode are identified. It has been also found that ν2 and ν3 + (ν14,ν15)0,0 are coupled in HCCN, while the mixing between ν3 and (2ν04, 2ν05)0,0 is seen in DCCN.

1. Introduction Cyanocarbene, HCCN, has been subjected to extensive experimental1-14 and theoretical analysis15-22 in the last 40 years. The ground electronic state of the radical was found to have triplet multiplicity.1-4 In the early work on HCCN, great effort was put into answering the question whether the radical is linear or bent.1-,5,15-19 A low barrier to linearity was first estimated by theoretical means.15-17 From a detailed microwave study of isotopically substituted HCCN, Brown et al. deduced in 1990 an unusually short CH bond length, what led to the conclusion that HCCN may possess a quasilinear structure.5 The ν5 fundamental transition was observed at 128.907968(40) cm-1 for HCCN11 and at 2243.72708(21) GHz (74.843 cm-1) for DCCN13 by means of the far-infrared laser magnetic resonance technique. Malmquist et al.18 and Koput21 determined onedimensional optimum energy profiles along the HCC angle and employed them in conjunction with the semirigid bender approach.23 However, this model was found to be insufficient for quantitatively describing the nuclear motion in HCCN.18,21 We refer to refs 11 and 21 for an overview of the previous experimental and theoretical studies on HCCN. In the present paper, we have undertaken a detailed fulldimensional quantum mechanical study of HCCN and DCCN in the ground electronic state. The motivation for this work is provided by the recent high-resolution spectra6-13 of HCCN and DCCN, which have established the frequencies of the ν5(HCC) bending and ν1(HC) stretching vibrations to high accuracy. The †

Part of the special issue “John C. Light Festschrift”. * Corresponding author. Present address: Papendiek 17, D-37073 Go¨ttingen, Germany. E-mail: [email protected]. ‡ E-mail: [email protected]. | E-mail: [email protected].

purpose of our study is to construct a potential energy surface and to compute rovibrational energy levels. Another aim is to study the quasilinear ν5 mode under realistic conditions of the vibrating-rotating molecule and to clarify characteristic energy patterns. Our numerical results for H/DCCN are of interest for studies on tetratomic molecules involving a hydrogen atom attached to a nearly linear chain of heavy atoms, such as, e.g., HCNO, HNCO, and HNCS. In section 2, we address the potential energy surface (PES) constructed in the present work. The computational method used to obtain rovibrational energies for J ) 0-4 is briefly explained in section 3. Our results for HCCN and DCCN are discussed in section 4. Section 5 contains our concluding remarks. 2. Potential Energy Surface A. Construction of the PES. The six-dimensional potential energy surface for the ground (3A′′) electronic state of HCCN was scanned by the multireference averaged coupled-pairfunctional (MR-ACPF) method,24-27 using the correlation consistent valence quadruple-ζ (cc-pVQZ) basis set. The Complete Active Space SCF (CASSCF) method was employed to optimize the orbitals. All calculations were performed with the MOLPRO suite of programs.28 The nuclear coordinate space was parametrized in terms of the bond-distance-bond-angle coordinates r1, r2, r3, R, β, and τ, shown in Figure 1a. Here, r1 ) r(H-C), r2 ) r(C-C), and r3 ) r(C-N). The two inplane bending angles, ∠HCC and ∠CCN, are denoted respectively by R and β, whereas τ stands for the dihedral angle measured from the cis side. The PES was calculated for the region specified by 1.4a0 < r1 < 4.1a0, 1.9a0 < r2 < 4.6a0, 1.7a0 < r3 < 4.3a0, 80° e R e 180°, 110° e β e 180°, and 0° e τ e 180°. In total, 1311 ab initio points were computed.

10.1021/jp056743u CCC: $33.50 © 2006 American Chemical Society Published on Web 02/18/2006

HCCN Radical in the Ground Electronic State

J. Phys. Chem. A, Vol. 110, No. 16, 2006 5521 TABLE 1: Expansion Coefficients Cijklmn (in Atomic Units) of Eq 1 for the Adjusted Six-Dimensional MR-ACPF PES Derived in This Work for the 3A′′ Electronic State of HCCNa

Figure 1. Description of the internal geometry of HCCN (a) in terms of the bond-distance-bond-angle internal coordinates r1, r2, r3, R, β,τ and (b) in terms of the orthogonal (diatom + diatom) internal coordinates R, d1, d2, θ1, θ2, and χ.

The calculated potential energies were in the range from -131.21429Eh to -130.98376Eh. The ab initio points were fitted to the six-dimensional analytical expression

V)

Cijklmnxi1xj2xk3Pnl (cR)Pnm(cβ) cos(nτ) ∑ ijklmn

(1)

where cR ) cosR, cβ ) cosβ, and l, m g n. For the stretching part, we use a polynomial expansion in the modified Morse coordinates xi (i ) 1, 2, 3) ref 1 xi ) [1 - e-ai(ri/ri -1)] ai

(2)

recommended by Meyer et al.29 The angular part of V is given by a multipole-like expansion, involving associated Legendre functions Pnl (cR), Pnm(cβ) and Chebyshev polynomials cos(nτ), as suggested previously.30 The functional form of eq 1 is invariant under the spatial inversion (τ f 2π - τ) and independent from the torsion angle τ for linear arrangements R,β ) 0,π. The analytical expansion of eq 1 has six nonlinear parameters, rref i and ai, which define the Morse coordinates xi of eq 2. In the preliminary stages of the fitting, we employed a nonlinear least-squares technique to obtain optimum values for both rref i and ai. In the final fitting, rref i and ai were held constant. For given rref i and ai, the expansion parameters Cijklmn were determined by means of a linear least-squares procedure. Several sets of {i, j, k, l, m, n} were tested in our fitting of eq 1 to the ab initio points. The complete quintic expansion, augmented by the one-dimensional radial contributions up to the seventh degree and by the three-dimensional angular part of the sixth degree, was found to be sufficient for achieving a reasonable balance between the number of fitting parameters, number of input ab initio points, and associated standard deviation σ. The chosen ansatz of {i, j, k, l, m, n} defines in total 336 fitting parameters Cijklmn. This number of parameters was reduced by an extensive search for all significantly determinable parameters. Eliminating stepwise all parameters smaller than twice their standard deviation and most of the parameters smaller than three times their standard deviation, we arrived at a functional six-dimensional expansion consisting of 125 terms.

ijklmn

Cijklmn

ijklmn

Cijklmn

ijklmn

Cijklmn

000000 200000 300000 400000 500000 600000 700000 020000 030000 040000 050000 060000 002000 003000 004000 005000 006000 007000 011000 012000 013000 014000 021000 022000 023000 031000 032000 101000 110000 120000 201000 210000 220000 000100 000200 000300 000400 000500 000600 000010 000110 000210

-0.0000836 0.7423250 -1.2529641 1.2940035 -1.1970301 0.9368890 -0.4416142 1.2419591 -2.9420632 3.9165475 -2.4605068 -0.6459143 1.6724940 -5.8652217 7.8379300 -5.7335005 2.0925619 -0.5721966 0.8697992 2.0180161 -2.1399322 -2.1759310 -1.1684469 3.2219910 3.6368392 -2.6150563 -3.3335227 0.0217615 0.0021936 -0.0357270 0.0144122 0.0086634 -0.1064960 -0.0868363 0.1745690 -0.0139063 -0.0101806 -0.0043744 -0.0008691 -0.3129661 -0.4723798 0.2439514

000310 000020 000120 000220 000030 000130 000230 000040 000140 000111 000211 000311 000411 000511 000121 000221 000321 000421 000131 000231 000331 000322 000422 000232 001010 001020 001030 001100 001110 001200 001210 001300 001400 001111 001211 001311 001222 002010 002020 002100 002111 003200

-0.0016181 -0.4014165 -0.4477221 0.1329813 -0.1959100 -0.2138155 0.0274065 -0.0387080 -0.0431464 0.0110927 0.0085826 0.0059613 -0.0010206 0.0003250 0.0078587 0.0084313 0.0033287 -0.0004627 0.0015849 0.0019083 0.0009184 0.0003446 0.0000690 -0.0000945 0.5458700 0.5965830 0.1533578 -0.5312928 -0.2864925 -0.3406292 -0.1293005 -0.0681403 -0.0189303 -0.0670670 -0.0731704 -0.0235973 0.0016467 -0.7734884 -0.2828110 0.0868854 0.0366603 0.1573389

010010 010020 010100 010110 010120 010200 010210 010220 010300 010400 010111 010121 010211 010311 011010 011100 011111 012100 020010 020100 020200 020300 020111 021010 021200 030010 030100 030200 100010 100020 100100 100200 100300 100111 100211 101100 110100 200100 200200 200300 300100

0.1333863 0.1001659 0.8068241 1.1190861 0.3713975 0.5152400 0.4998266 0.1619246 0.0558572 0.0233452 0.0762339 -0.0115625 0.0817841 0.0277756 0.1585303 -0.2973279 0.0491293 1.4686769 0.0403994 -0.6353160 -0.5613150 -0.1753884 -0.0991483 0.2756852 0.3383835 0.1357801 -0.2013531 -0.1253704 0.1117913 0.0411506 -0.1076529 -0.0385778 -0.0095890 -0.0076520 -0.0051029 0.0284819 0.0455224 -0.1073827 -0.0629688 -0.0198132 0.0346831

ref The reference distances rref 1 ) 2.01964a0, r2 ) 2.49804a0, and ) 2.25743a0 and the exponential parameters a1 ) 0.4434, a2 ) 0.8009, and a3 ) 0.3948 are used for the Morse coordinates xi defined by eq 2. a

rref 3

For the final fitting, we set a standard deviation of a ith date point of energy Ei to be Ei + 2000 cm-1, as previously done, e.g., by Meyer et al.29 The standard deviation of the weighted 125-term expansion was 4.5 cm-1. This σ is only 0.9 cm-1 larger than the σ of 3.6 cm-1 obtained for the initial 336-term representation. Preliminary rovibrational calculations using the 125-term PES expansion in combination with the DVR(+R)+FBR method described below (section 3) gave the ν1 and ν5 fundamentals for HCCN, which deviate by 14.5 and 5.4 cm-1 from experiment. Our analytical representation was hence refined by adjusting the parameters C200000 and C000200 to available experimental data for the ν1 stretching6 and ν5 bending11 frequencies. The expansion coefficients Cijklmn of the refined 125-term PES are listed in Table 5. The original values of C200000 and C000200 were 0.7468463 and 0.1744197, respectively. For the equilibrium configuration, the MR-ACPF PES gives the harmonic frequencies of 3355, 1803, 1173, 553, 441, and 400 cm-1 for HCCN and of 2488, 1794, 1134, 506, 431, and 331 cm-1 for DCCN.

5522 J. Phys. Chem. A, Vol. 110, No. 16, 2006

Mladenovic´ et al.

TABLE 2: Fundamental Transitions of HCCN and DCCN (in cm-1)a HCCN this work ν1 ν2

3246.66 [25.94] 1733.71

ref

3246.6573(5) [26.80] ∼1750 1735

ν14

[49.90] 1178.57 [23.03] 430.02

ν15

[-20.39, -47.50] 128.96

ν3

DCCN expt

[17.88, -21.98]

this work

6 6 8 3

1178.5

3

365(15) 383(20) 458 [-21.05, -48.95] 128.907968(40) 145(15) [18.16, -22.40]

7 6 3 7 11 7 7

2436.64 [24.44] 1730.27 [67.84] 1150.65 [-9.86] 407.48 [-18.36, -53.68] 75.65 [8.67, -41.12]

expt

ref

2436.3723(7) [25.42]

9 9

1729.5

3

1127

3

367(15)

7

405 [-19.07, -53.86] 74.845(2) 90(15) [8.70, -41.93]

3 7 9 7 7

The number in brackets gives the difference B0 - BV (in MHz) between the ground-state rotational constants B0 and effective rotational constant BV for the vibrational state V. For the MR-ACPF PES, B0 is calculated to be 0.364872 and 0.329029 cm-1 for HCCN and DCCN, respectively. The corresponding experimental values due to McCarthy et al.7 are 0.366465 and 0.330441 cm-1 for HCCN and DCCN, respectively. a

TABLE 3: Transitions Involving the Bending ν5 and ν4 Vibrationsa HCCN this work 2ν25 2ν05 (ν14,ν15)0, 1 (ν14,ν15)0, 0 (ν14,ν15)2, p 3ν35 3ν15 4ν45 4ν25 4ν05 2ν04 2ν24 3ν14 3ν34 4ν04 4ν24 4ν44 a

343.44 [-27.53] 460.97 [-12.86] 564.77 [-34.12] 575.81 [-26.59] 560.10 618.23 [-44.79] 725.75 [-9.43, -38.50] 939.71 1007.17 [-45.00] 1088.61 [-37.83]

DCCN expt

ref

341.73 340(15) [-1.01] 435(20) [-16.60] 525(20) [-35.84] 540(20) [-25.55]

12 7 7 7 7 7 7 7 7

625(20) [-2.53]

7 7

851.09 860.41 1273.46 1290.79 1686.79 1695.14 1721.89

this work 210.89 [-32.03] 317.92 [-15.91] 502.35 [-40.35] 516.77 [-107.8] 495.97 386.84 [-43.69] 514.90 [56.56, -53.20] 593.89 697.41 [-96.43] 744.02 [-58.15]

expt

ref

207.95 215(15) [-30.71] 311(15) [-18.40]

12 7 7 7 7

400(20) [-43.12]

7 7

837.76 821.85 1261.92 1240.88 1689.92 1681.57 1663.58

For more details, see Table 2.

B. Topography of the PES. The two-dimensional maps of the adjusted six-dimensional MR-ACPF PES are displayed in Figure 2. These are made by freezing the remaining four coordinates at their equilibrium values. The contour plots involving {R,τ} and {β,τ} (the second and third plot in the second row of Figure 2) clearly show that HCCN is a planar molecule possessing only a trans minimum. The geometrical parameters of the PES minimum are given by re1 ) 2.020a0, re2 ) 2.499a0, re3 ) 2.257 a0, Re ) 146.6°, βe ) 174.7°, and τe ) 180°. As seen in Figure 2, the coupling between the HCC bending angle R and the CC bond length r2 is especially prominent. Figure 2 also shows the variation of r1, r2, r3, β, τ along the minimum energy path (MEP) in the direction of R, as well the variation of τ along the MEP for β. The MEP along a chosen coordinate is computed by minimizing V with respect

to the remaining five coordinates. Along the MEP for R, the optimum triple (r1, r2, r3) assumes the value of (2.042a0, 2.565a0, 2.236a0) for R ) 120° and of (2.009a0, 2.447a0, 2.281a0) for R ) 180°, yielding changes of -0.033a0, -0.118a0, and +0.045a0 for r1, r2, and r3, respectively. Malmquist et al.18 and Koput21 also found that the straightening of the HCC bending angle decreases r2(CC) and increases r3(CN) in HCCN. The variation of the potential energy along the MEP for R, β, and τ is displayed in Figure 3. The barrier to linearity for HCC and CCN is calculated to be 287 and 84 cm-1 along the respective MEPs. The torsion barrier, i.e., the energy difference between the optimum cis and trans HCCN conformations, is 84 cm-1. In Figure 3, we also show the one-dimensional cut, Vcut, of the PES along the angle R, which is obtained by freezing the other five coordinates at their equilibrium values. On Vcut,

HCCN Radical in the Ground Electronic State

J. Phys. Chem. A, Vol. 110, No. 16, 2006 5523

TABLE 4: Combination Bands Involving the Stretching ν1 Vibrationa HCCN this work ν1 +

ν15

3355.504(2) [45.33, 2.76] 3226.597(2)

ν1 + ν14 - ν14

3356.40 [43.73, 1.70] 3227.44 3690.06 [7.48] 3229.10 3552.69 [-2.68] 3423.73 3809.80 3464.17 3209.25 3673.97 [5.41, -22.63] 3243.95

ν1 + ν15 - ν1 ν1 + 2ν25 - ν1 ν1 + 3ν35 - ν1 ν1 + 2ν05 - ν1

109.01 303.84 558.77 443.40

ν1 + ν15 - ν15 ν1 + 2ν05 ν1 + 2ν05 - 2ν05 ν1 + 2ν25 ν1 + 2ν25 - ν15 ν1 + 3ν35 ν1 + 3ν35 - 2ν25 ν1 + 2ν25 - 2ν25 ν1 + ν14

a

DCCN expt

ref

this work

expt

10 10 6

2500.26 [37.18, -16.53] 2424.60 2750.61 [11.85] 2432.68 2623.67 [-8.74] 2546.71 2788.33 2575.46 2412.78 2836.22 [9.03, -36.38] 2428.74

2499.572(1) [32.72, -19.11] 2424.727(1)

9 9 9

[12.29] 2432.1316(11)

9 9

[-5.07] 2544.7432(9)

9 12

2573.4795(7) 2411.636(3)

12 9

63.20 183.22 345.06

12 12 12

[24.68] 3420.6663(6)

12 12

3460.9116(8) 3207.8482(6)

12 12

[8.76, -23.32] 3244.084(1)

6 6

108.85 302.92 555.98

12 12 12

62.96 185.07 347.76 313.97

ref

For more details, see Table 2.

-1 TABLE 5: Bending Frequencies νeff 5 (in cm ) of HCCN and DCCN for the Ground-State (gs), ν1, ν2, and ν3a

HCCN gs ν2 ν3 ν1

DCCN

νeff 5

Vlin eff

θmin 1

νeff 5

Vlin eff

θmin 1

129 137 116 109

257 179 323 332

148 150 147 146

76 92 59 63

274 193 333 346

150 152 149 148

a -1 The height of the barrier to linearity Vlin eff (in cm ) and the (in deg) are additionally shown. location of the minimum θmin 1

cm-1

the arrangement with R ) 180° is lying 528 above the absolute minimum and is thus 241 cm-1 higher than the corresponding MEP value. A similar one-dimensional cut of the PES in the direction of β shows no noticeable deviation from the corresponding MEP.

expression of the kinetic energy operator Tˆ from ref 31 can be rewritten as a sum of several contributions

Tˆ ) Tˆ vib + Tˆ rot + Tˆ cv + Tˆ cr

(3)

where Tˆ vib and Tˆ rot stand for the vibrational and rotational parts. The term Tˆ cv is Coriolis-like kinetic coupling in the vibrational angular subspace, whereas Tˆ cr is Coriolis-like kinetic coupling between the rotation and vibrational angular motion. The four contributions are explicitly given by

Tˆ vib ) Tˆ str + Tˆ bend

(4)

2Tˆ bend ) f (R, d2)lˆ22 -f (R, d1)[p2(∂2θ1 + cotθ1∂θ1) - (Jˆ z - ˆl 2z)2/sin2θ1] (5) 2Tˆ rot ) fR[Jˆ 2 - 2(Jˆ z - ˆl 2z)2]

3. Method for Rovibrational Calculations

Tˆ cv ) fR(- iplˆ2y∂θ1 + cotθ1ˆl 2xˆl 2z)

In the bound state variational calculations for the rovibrational energies of H/DCCN, we have employed the orthogonal (diatom + diatom) internal vectors R, d1, d2 in the body-fixed formulation.31 The vectors d1 and d2 are chosen to be the HC and CN bond distance vectors, whereas R is the vector joining the centers of mass of the HC and CN subunits. The internal vectors are parametrized in terms of the three distances R, d1, and d2, the two bending angles θ1 and θ2, and the dihedral angle χ, schematically shown in Figure 1b. The reduced masses µR, µd1, and µd2 associated with R, d1, and d2 are given for HCCN explicitly by µd1 ) mHmC/(mH + mC), µd2 ) mCmN/(mC + mN), and µR ) mHCmCN/M, where mA stands for the mass of the unit A and M for the overall mass of HCCN. The two-vector embedded body-fixed (BF) reference frame is chosen such that the z-axis of the BF is aligned with the vector R and the z ∧ x BF plane is the plane spanned by R and d1. The MR-ACPF PES minimum in terms of the (diatom + diatom) coordinates is given by Re ) 3.845a0, θe1 ) 149.6°, θe2 ) 2.31°, and χe ) 180°, while de1 ) re1 and de2 ) re3. For an entirely orthogonal set of the internal vectors, such as the (diatom + diatom) vectors of Figure 1b, the body-fixed

Tˆ cr ) fR[ipJˆ y∂θ1 - Jˆ xˆl 2x - Jˆ yˆl 2y - Jˆ zˆl 2z

(6) (7)

+ cotθ1(Jˆ xJˆ z - Jˆ xˆl 2z - ˆl 2xJˆ z)] (8) where Jˆ is the angular momentum of the nuclear frame and Tˆ str ) Tˆ str(R) + Tˆ str(d1) + Tˆ str(d2) stands for the kinetic energy of the three stretching vibrations. The term Tˆ cv affects only vibrational energies and may lead to vibrational l-type doubling/ resonance effects. Rotational l-type doubling/resonance effects may arise from Tˆ cr. In eqs 4-7, ∂θ1 stands for ∂/∂θ1, ∂2θ1 for ∂2/∂θ12, and fR for 1/µRR2. The function f (R, di), defined as 1/µRR2 + 1/µdidi2, is the inverse of the reduced mass associated with the bending coordinate θi. The vibrational angular momentum operator associated with d2 is ˆl2 ) ˆl2(θ2,χ). The body-fixed projections of Jˆ are Jˆ x, Jˆ y, Jˆ z. The standard rotation-angular basis in the body-fixed formulation is given for tetratomic molecules by the following product31

P h l1|k-K|(cosθ1)Yl2k(cosθ2,χ)|JKM〉

(9)

5524 J. Phys. Chem. A, Vol. 110, No. 16, 2006

Mladenovic´ et al.

Figure 2. Two-dimensional contour plots of the six-dimensional MR-ACPF PES, obtained by freezing the remaining four coordinates at their equilibrium values. Contour lines are drawn in intervals of 250 cm-1 with the first contour placed at 50 cm-1. The dashed lines show the variation of the internal coordinate shown on the y-axis along the minimum energy path for the coordinate displayed on the x-axis.

Figure 3. Minimum energy path (solid line) along the HCC angle R, the CCN angle β, and the torsion angle τ. The one-dimensional cut Vcut along R is given by the dashed line.

where P h kl (cosθ) are normalized associated Legendre functions k and Yl2(cosθ2,χ) are spherical harmonics with the CondonShortley phase convention.32 The quantum number of the projection of ˆl2 onto the z-axis of the BF is k. The symmetric top eigenfunctions33 are denoted by |JKM〉, where M and K are the quantum numbers for the space-fixed Z-projection and the body-fixed z-projection of Jˆ , respectively. The action of the kinetic energy operator Tˆ of eq 3 on the basis functions of eq 9 can be solved analytically.31 The total rotational angular momentum J, its space-fixed projection M, and the parity p are strictly conserved quantum numbers for the eigenstates of a tetratomic molecule. For a given J, a parity-adapted angular basis can thus facilitate the computation. A detailed description of the parity-adapted rotation-angular basis functions employed in the present work can be found in ref 31. In actual computations, the discrete variable representation (DVR) for the bending angles θ1and θ2 is used together with the finite basis representation (FBR) of the dihedral angle χ, the eigenfunction basis in d1, d2, and the discretized Jacobi coordinate R. The DVR+FBR is combined with a sequential diagonalization/truncation scheme. This approach is termed DVR(+R)+FBR34 and involves no dynamical approximation. In the DVR(+R)+FBR method, contributions to the rovibra-

tional Hamiltonian are taken sequentially into account. After the inclusion of new term(s), the higher-dimensional Hamiltonian matrix is calculated and diagonalized, while the basis set is truncated by retaining only a preselected number of eigenvectors for subsequent computations. For chosen K, k, Rz and parity p, where Rz stands for a discrete point in R, the terms Tˆ bend + Tˆ str(d1) + Tˆ str(d2), and Tˆ rot are first considered together with the potential energy contribution diagonal in k. In the next step, the term Tˆ cv is added along with the potential energy contribution nondiagonal in k, resulting in five-dimensional eigenvalue problems for each of the pairs (K, Rz). The inclusion of the kinetic energy in R leads to the formulation of six-mode eigenvalue problems for a fixed K. The final nine-dimensional rovibrational Hamiltonian matrix is constructed by including the term Tˆ cr and is of relatively modest size owing to a very efficient contraction scheme of the DVR(+R)+FBR approach. At different stages of calculation, the Hamiltonian matrices are diagonalized by conventional dense matrix algorithms to give both eigenvalues and eigenfunctions. A. Adiabatic Projection Scheme. The computational strategy of the DVR(+R)+FBR approach readily permits the construction of the adiabatic (zero-order) basis, consisting of eigenvectors computed in the adiabatic R-stretch approximation.

HCCN Radical in the Ground Electronic State

J. Phys. Chem. A, Vol. 110, No. 16, 2006 5525

Figure 4. Low-lying states of HCCN calculated for total rotational angular momentum J and parity p. The quantum number assignments are made according to the linear molecular model. The quantum labels of the states denoted by a, b, and c are 3ν15, (ν14, 2ν25)1,p, and (ν14, 2ν05)1,p, respectively. For the states shown, the quantum label K for the projection of Jˆ on the body-fixed z-axis is essentially a good quantum number.

The five-mode eigenstates calculated for a grid of the discrete points Rz form adiabatic potentials Viadi for R. The stretching states |R〉, adiabatically separated from the remaining five vibrational modes, can be obtained for the ith five-mode state by employing the adiabatic Hamiltonian H ˆ iadi ) H ˆ iadi(R), consisting of the kinetic energy Tˆ str ) Tˆ str(R) and the adiabatic potential Viadi. To monitor the evaluation of the adiabatic representation within the DVR(+R)+FBR computational scheme, the quantum state character correlation scheme is used. The latter scheme was originally developed for the adiabatic torsion approximation and the DVR(6) approach.35 The adiabatic vectors |R, i) ) |R〉| i〉 are used in the DVR(+R)+FBR for the calculation of adiabatic expansions of numerically exact wave functions | n〉

|n〉 )

CR,in|R, i) ∑ R, i

(10)

The latter expansion supplies a zero-order origin of |n〉 by locating the dominant adiabatic contribution among all |R, i), providing, thus, an important foundation for quantum number assignments. The adiabatic projection of eq 10 is easily extended to include also the J * 0 situation. B. Numerical Details. The numerical parameters used to compute the rovibrational energy levels of H/DCCN for J ) 0-4 by means of the DVR(+R)+FBR are summarized here. The pointwise representation of the Jacobi distance R is obtained in the spirit of the potential optimized DVR method of Echave and Clary.36 For HCCN, we used 12 discrete Rz points distributed nonevenly between 3.54a0 and 4.43a0. Three eigenfunctions for d1 and five eigenfunctions for d2 were constructed for the reference potential determined by the equilibrium H/DCCN geometry. The maximum value kmax for the projection quantum number k was 5 in J ) 0 and 9 in J ) 4 calculations. The integrals over the dihedral angle χ were solved by means of Gauss-Chebyshev quadrature of order 11. For the bending angles, Gauss-Legendre DVR points were used, determined from Plkimax of k ) 0, 1,‚‚‚, kmax with lmax of 35 for θ1 and lmax of 1 2 95 for θ2. The Gauss-Legendre DVR set was truncated, keeping only the points distributed between 90 and 180° for θ1 and between 0 and 40° for θ2. The primary basis included 347400 functions for J ) 0 and 3495960 for J ) 4. The truncation of the basis set was made with the help of the cutoff parameters

Nkmax, Nzmax, and NKmax, which give the number of eigenvectors to be saved after the diagonalization of the (K, Rz, k)-blocks, (K, Rz)-blocks, and K-blocks, respectively. The parameters Nkmax, NRmax, NKmax were chosen to be 150, 800, and 500, respectively. The size of the six-dimensional matrices of the K-blocks was thus 9600. The size of the final rovibrational Hamiltonian matrix was 1000, 1500, 2000, and 2500 for J ) 1, 2, 3, 4, respectively. 4. Results Selected results obtained for J ) 0-4 by means of the DVR(+R)+FBR method and the MR-ACPF PES are shown in Tables 2-5 and Figures 4-7. In these tables, the level energies are given relative to the ground vibrational state, calculated to be 3787.357 cm-1 for HCCN and 3289.913 cm-1 for DCCN. A. Quantum State Labels. Following the conventional linear molecule notation, the state assignments can be made in terms of ν1, ν2, ν3, νl44, and νl55. The mode labels ν1, ν2, ν3, ν4, and ν5 refer respectively to CH stretch, CN stretch, CC stretch, CCN bend, and HCC bend, whereas l4, l5 denote the vibrational angular momentum quantum numbers. The permitted values of li are Vi, Vi - 2, Vi - 4, ..., - Vi for i ) 4, 5. Since K, k g 0 holds in the parity-adapted formulation,31 we employ in the present work the quantum label (V4ν|l4 4|, V5ν|l5 5|)K,p, where p stands for parity (p ) 0 for even parity and p ) 1 for odd parity). In the body-fixed formulation used here, l5 correlates to the projection quantum number k. Angular momentum conservation additionally restricts l4, l5 to l4 + l5 ) l ) K e J, as seen in eq 9 and ref 31. For quantum states given by V4 ) 0, we use the label V5νK5 , and similar for the V5 ) 0 situation. The low-lying rovibrational states of HCCN calculated for J ) 0, 1, and 2 and both parities are found in Figure 4. One may note that the levels ν14 and ν15 occur only for J g 1, whereas the states 2ν25, 2ν24, and (ν14,ν15)2,p become accessible for J g 2. The two K ) 0 levels shown by the dashed lines in Figure 4 are assigned (ν14,ν15)0,p in the notation used here or -1 ν(1 4 + ν5 . The alternative assignment for the two K ) 2 levels 1 1 2,p (1 (ν4,ν5) of Figure 4 is ν(1 4 + ν5 . In the bent molecule notation, the levels 2ν05, (ν14,ν15),0,0 and 1 1 0,1 (ν4,ν5) are respectively νb5, νb4, and νb6, whereas ν15 and 2ν25 correlate with the K ) 1 and K ) 2 (a-type) transitions from

5526 J. Phys. Chem. A, Vol. 110, No. 16, 2006

Mladenovic´ et al.

Figure 5. Wave function probability amplitudes for the ground vibrational state (gs), 2ν05, (ν14,ν15),0,0 and 2ν04 of HCCN for J ) 0. The remaining four coordinates are kept constant at their equilibrium values. The contours seen in the (θ1,θ2) plane show where the probability amplitude takes 10%, 20%, 30%, 40%, 50%, 60%, 70%, 80%, 90% of its maximum value.

the ground state. The levels assigned ν14 in Figure 4 are, however, not easy to interpret by means of the bent molecular model. The two-dimensional (θ1,θ2) cuts of the wave function probability amplitudes for the ground vibrational state, 2ν05, (ν14,ν15),0,0 and 2ν04 of HCCN are shown in Figure 5. The welldefined nodal patterns seen for 2ν05 and 2ν04 in Figure 5 clearly imply the excitation of the θ1 bending motion in 2ν05 and the excitation of the θ2 bending motion in 2ν04, confirming the linear molecule convention for HCCN. A striking feature in Figure 5 is a large wave function amplitude seen at the linear geometry, θ1 ) π, even in the ground vibrational state. This property is also visible in Figure 6, where the wave function probability amplitude integrated over R, d1, d3, θ2, χ and three Euler angles are shown for the ground vibrational state, ν15, 2ν05, 3ν15, and 4ν05. Although the ground vibrational state at 3787 cm-1 lies ≈90 cm-1 below the barrier to linearity, the ground-state wave function amplitude differs significantly from zero at θ1 ) π. For the levels with l5 ) 1, the wave function amplitude at θ1 ) π is strictly zero as a consequence of the rotational symmetry, Figure 6 and eq 9. The states with l5 * 0 exceed the barrier to linearity. The vibrationally averaged bending amplitude 〈θi - θie〉 is calculated as arccos[〈cos(θi - θie)〉]. For the ground vibrational state, 2ν05, (ν14,ν15),0,0 and 2ν04 of Figure 5, the amplitude 〈θ1 θe1〉 describing the quasilinear mode is calculated to be 11.4,

16.0, 10.6, and 12.0°, respectively. For 〈θ2 - θe2〉, we obtained 4.2, 4.2, 6.3, and 8.0°. B. Fundamental Transitions. The fundamental transitions of HCCN and DCCN are summarized in Table 2. The MRACPF PES, adjusted to reproduce experimental wavenumbers for ν15 and ν1 of HCCN (section 2), provides a good description of ν5 and ν1 also for DCCN. As seen in Table 2, the ν1 and ν5 frequencies calculated for DCCN differ by +0.27 and +0.81 cm-1 from the high-resolution spectroscopic results of Sun et al.9 So far, no high-resolution spectroscopic study has been reported for the ν3 and ν4 vibrations of either HCCN or DCCN. Dendramis and Leroi3 derived from matrix isolation infrared spectra five fundamental transitions for various isotopically substituted forms. McCarthy et al.7 and Morter et al.6 estimated the transition energies of ν4 and ν5 from relative intensity measurements in the microwave and infrared region, respectively. Miller et al.8 recorded the high-resolution infrared spectrum of HCCN in the region of the ν2 stretching fundamental and found six separate vibrational bands. However, no vibrational identification could be made. Fermi resonance between ν2 and ν3 + (ν14,ν15)0,0 was suggested in ref 8 as a probable intensity borrowing mechanism in the region of ν2. In the present work, the energy spectra calculated for HCCN and DCCN are studied in great detail with the help of the adiabatic expansions of eq 10. This analysis led to the

HCCN Radical in the Ground Electronic State

J. Phys. Chem. A, Vol. 110, No. 16, 2006 5527

Figure 6. Wave function probability amplitudes integrated over the five internal coordinates R, d1, d2,θ2,χ and three Euler angles for the ground vibrational state (gs), ν15, 2ν05, 3ν15, and 4ν05 of HCCN. The levels shown are calculated for J ) 1 and p ) 1. The effective groundstate potential energy profile along θ1 is denoted by Veff0. The horizontal lines indicate the level energy.

identification of several resonances involving the fundamental stretching ν2 and ν3 vibrations. For HCCN, we have

|ν2〉 ≈ 0.799| 0, 9) - 0.439| 1, 2)

(11)

|ν3 + (ν14,ν15)0, 0〉 ≈ 0.534| 0, 9) + 0.725| 1, 2) and

|ν3〉 ≈ -0.428| 0, 5) + 0.818| 1, 0)

(12)

|(ν14, 3ν15)0, 0〉 ≈ -0.803| 0, 5) - 0.447| 1, 0) The HCCN levels ν2, ν3 + (ν14,ν15),0,0 ν3, and (ν14, 3ν15)0,0 occur at 1733.7, 1759.2, 1178.6, and 1185.8 cm-1, respectively. The ν2 state of DCCN displays no resonance mixing. For ν3 of DCCN, we obtain

|(2ν04, 2ν05)0, 0〉 ≈ 0.715| 0, 7) - 0.567| 1, 0)

(13)

|ν3〉 ≈ 0.560| 0, 7) + 0.720| 1, 0) where (2ν04, 2ν05)0,0 and ν3 lie respectively at 1131.2 and 1150.6 cm-1. The adiabatic projections of eqs 11-13 are given in terms of the two dominant components, providing at least 81% of the full-dimensional wave function. For the zero-order vectors, we employ the notation used by the code, where |R, i) stands for the adiabatic level with R quanta in the stretching (R) vibration and i quanta in the five-mode vibration, as described in section 3.1.

Adiabatic projection methods are useful for quantitative studies on separability of internal molecular motions and for making quantum number assignments to the calculated fulldimensional rovibrational states. The adiabatic (zero-order) basis is, however, coordinate dependent. To verify eqs 11-13, we studied several other orthogonal (Jacobi type) descriptions and found that the levels ν2 and ν3 occur as zero-order mixtures also for Jacobi type coordinates. These test calculations led additionally to the conclusion that internal motion of HCCN is most separable in the (diatom + diatom) description of Figure 1b. C. Transitions Involving the Bending ν5 and ν4 and Stretching ν1 Vibrations. The transitions involving the bending ν5 and ν4 vibrations are summarized in Table 3, whereas Table 4 displays the transitions involving the stretching ν1 vibration. The level energies in Tables 3 and 4 are in excellent agreement with high-resolution spectroscopic data of refs 6, 9, 10, and 12, showing a maximum deviation of +3.3 cm-1 for the combination band ν1 + 3ν35 - 2ν25 calculated at 3464 cm-1 for HCCN. One may note that in Table 4 we also list the energies ν1 + nν(n - ν1 (n ) 1, 2, 3), which are obtained in the 5 semirigid bender calculations of Hung et al.12 The energy levels of the pure ν4 and pure ν5 bending vibrations of HCCN are graphically presented in Figure 7. The states of ν4 follow in particular closely the energy pattern of the two-dimensional (anharmonic) oscillator. To quantify the degree of quasilinearity of molecular vibrations, Yamada and Winnewisser37 introduced the parameter γ0, given by

γ0 ) 1 - 4[E(ν1i ) - E0]/[E(2νi0) - E0]

(14)

Our results from Tables 2 and 3 give for the quasilinear ν5 bending mode γ0 of -0.12 for HCCN and of +0.05 for DCCN. On the other hand, we obtain γ0 of -1.02 and -0.95 for the ν4 of HCCN and DCCN, respectively, as expected for a degenerate bending mode in well-behaved linear molecules.37 In Figure 7, the states with li * 0 are doubly degenerate. Pronounced l-type splitting is, however, found only for the K ) 1 states. For ν15 and ν14 of HCCN, the l-type doubling constant qv is determined to be 40.1 and 27.6 MHz, which compare well to 40.6 and 27.9 MHz predicted by McCarthy et al.7 For the ν15 and ν14 states of DCCN, the calculated qv values of 49.9 and 35.6 MHz agree within 0.8 MHz with the experimental data7 of 50.6 and 34.8 MHz. For the states ν1 + ν5, we obtained qv of 42 MHz for HCCN and of 54 MHz for DCCN, which should be compared to qv of 44 and 52 MHz estimated by Morter et al.6 and Sun et al.9 As clearly seen in Figure 7, the level separation between V4νl40 and V4νV44 is much smaller than the separation between V5νl50 and V5νv55, where l0 ) 0 for Vi even and l0 ) 1 for Vi odd. For HCCN and DCCN, the level 2ν05, e.g., lies 118 and 108 cm-1 above 2ν25, respectively. An interesting difference between HCCN and DCCN is found in the ordering of the 2ν04 and 2ν24 states, Table 3. Whereas the levels 2νl44 of HCCN show the expected ordering with the level 2ν04 being 9 cm-1 below 2ν24, the opposite is found for DCCN. In the latter case, the level 2ν04 is 16 cm-1 above 2ν24, The bending HCC frequency obtained for four different stretching excitations are found in Table 5. The effective frequencies νeff 5 are calculated as Eνi+ν5 - Eνi using the (J ) 1, p ) 1) values. As seen in Table 5, the excitation of ν1 and ν3 -1 and lowers νeff 5 of HCCN approximately by 20 and 13 cm

5528 J. Phys. Chem. A, Vol. 110, No. 16, 2006

Mladenovic´ et al.

Figure 7. Energy levels of the pure ν4 (left) and pure ν5 (right) bending vibrations of HCCN, obtained for J ) 4. One may note that the levels associated with li * 0 are doubly degenerate.

rotational constant Bv for the vibrational state V. The Bv values were determined by a least-squares fit of the calculated rovibrational energies to the following approximate expression

Ev ) Tv + BvJ(J + 1)

Figure 8. Effective potential energy profiles and the minimum energy path along the angle θ1. For the purpose of comparison, each of the curves shown is measured relative to its minimum. The minimum energy path is denoted by MEP. -1 νeff 5 of DCCN by 13 and 17 cm , yielding changes up to 22%. increases by 12 and 16 cm-1 upon On the other hand, νeff 5 excitation of ν2 in HCCN and DCCN, respectively. To gain a more general insight into the dynamics of the quasilinear ν5 mode, we have also carried out detailed investigations of effective (adiabatic) potential energy profiles along the angle θ1. For various excitations of other vibrational modes, the effective bending profiles were calculated in the spirit of the DVR(+R)+FBR strategy (section 3) by solving fivedimensional eigenvalue problems for a chosen θ1 grid. The effective θ1 profiles obtained for the five-mode ground state, ν1, ν2, and ν3 are depicted in Figure 8, whereas the height of the barrier to linearity, Vlin eff, and the location of the , are presented in Table 5. As seen there, the minimum, θmin 1 height of the barrier to linearity on the five-mode ground-state profile is 30 cm-1 lower than the corresponding MEP value of 287 cm-1. Upon excitation of ν1 and ν3, the barrier to linearity increases by 75 and 66 cm-1 for HCCN and by 72 and 59 cm-1 for DCCN relative to the ground-state value. At the same time, -1 for HCCN ν2 excitation introduces lowering of Vlin eff by 78 cm and by 81 cm-1 for DCCN. According to Table 5, the effective bending frequency νeff 5 nicely correlates with the effective , i.e., the higher the barrier to linearity, the barrier height Vlin eff lower the effective bending frequency. D. Rotational constants. In addition to the transition energies, Tables 2-4 also give the difference ∆Bv ) B0 - Bv between the ground-state rotational constant B0 and the effective

(15)

where Tv stands for the term energy. The ground-state rotational constant B0 is calculated to be 0.364872 cm-1 for HCCN and 0.329029 cm-1 for DCCN. They are 0.0016 and 0.0014 cm-1 smaller than the experimental values of 0.366465 and 0.330441 cm-1 derived from the millimeter spectra.7 On the other hand, the theoretical results of Koput,21 which read 0.367468 and 0.331676 cm-1, overestimate B0 by 0.0010 and 0.0012 cm-1. One may also note that the equilibrium parameters (re1, r22, re3,Re,βe) obtained here (section 2.2) differ by (0.0a0, 0.011a0, -0.016a0, -1.7°, 0.7°) from the best estimate of Koput (2.020a0, 2.510a0, 2.241a0, 144.9°, 175.4°) determined at the CCSD(T) level of theory. The ∆Bv calculated for the fundamental transitions in Table 2 agree remarkably well (within 1.5 MHz) with the experimental values3,6-9,11 for both HCCN and DCCN. Equally good agreement is also found for ν15, 2ν25, and 3ν35 of DCCN in Tables 2 and 3. In the case of HCCN, the experimental and theoretical rotational constants for ν15, 2ν25, and 3ν35 exhibit, however, different trends. While the experimental rotational constants7 are nearly constant, our results in Table 3 for 2ν25 and 3ν35 of HCCN show rather prominent changes of ∆Bv with the l5 excitation. Similar behavior is consequently seen for the state ν1 + 2ν25 in Table 4. McCarthy et al.7 also reported three l ) 0 states at 435, 525, and 540 cm-1 for HCCN and one l ) 0 state at 311 cm-1 for DCCN, each possessing a rotational constant somewhat larger than the ground-state B value, Table 3. Making use of the ∆Bv computed in the present work, the observed l ) 0 states of HCCN are readily assigned 2ν05, (ν14,ν15),0,1 and (ν14,ν15),0,0 while the l ) 0 state of DCCN is labeled 2ν05. For these four states, the calculated and experimental ∆Bv results agree to better than 4 MHz in Table 3. 5. Summary In this paper, we have presented a detailed study of the rovibrational energy spectra of HCCN and its isotopomer DCCN. For the ground electronic state we first developed a six-dimensional MR-ACPF potential energy surface, which is subsequently adjusted to give best agreement with experimental data for ν1 and ν5 of HCCN. Quantum mechanical calculations were carried out for total rotational angular momentum J ) 0-4 using the DVR(+R)+FBR method. In addition to the

HCCN Radical in the Ground Electronic State transition energies, our study also includes the determination of spectroscopic parameters, such as effective rotational constants and l-type doubling constants. The rovibrational energy levels of both HCCN and DCCN are studied in great detail by means of vibrationally averaged geometries and adiabatic projection methods. Resonance interaction of the ν2 stretching mode with ν3 + (ν14,ν15)0,0 is found in HCCN, whereas the ν3 stretching mode of DCCN is coupled with (2ν04, 2ν05).0,0 To gain a more general insight into the dynamics of the quasilinear ν5 bending mode, the effective (adiabatic) energy profiles along θ1 have been calculated for various five-mode excitations. For both HCCN and DCCN, we have found a clear increase in the height of the effective barrier to linearity for the excited ν1 and ν3, resulting in the lower effective ν5 bending frequency. On the other hand, ν2 excitation lowers the barrier to linearity and increases the effective ν5 frequency. Acknowledgment. This work has been supported by the Deutsche Forschungsgemeinschaft through SFB 357 (Molekulare Mechanismen unimolekularer Prozesse). Computer time was made available to M.M. through a grant from the John-vonNeumann Institut fu¨r Computing (NIC) at Forschungszentrum Ju¨lich. References and Notes (1) Bernheim, R. A.; Kempf, R. J.; Humer, P. W.; Skell, P. S. J. Chem. Phys. 1964, 41, 1156. (2) Wasserman, E.; Yager, W. A.; Kuck, V. J. Chem. Phys. Lett. 1970, 7, 409. (3) Dendramis, A.; Leroi, G. E. J. Chem. Phys. 1977, 66, 4334. (4) Saito, S.; Endo, Y.; Hirota, E. J. Chem. Phys. 1984, 80, 1427. (5) Brown, F. X.; Saito, S.; Yamamoto, S. J. Mol. Spectrosc. 1990, 143, 203. (6) Morter, C. L.; Farhat, S. K.; Curl, R. F. Chem. Phys. Lett. 1993, 207, 153. (7) McCarthy, M. C.; Gottlieb, C. A.; Cooksy, A. L.; Thaddeus, P. J. Chem. Phys. 1995, 103, 7779. (8) Miller, C. E.; Eckhoff, W. C.; Curl, R. F. J. Mol. Struct. 1995, 352/353, 435. (9) Sun, F.; Kosterev, A.; Scott, G.; Litosh, V.; Curl, R. F. J. Chem. Phys. 1998, 109, 8851.

J. Phys. Chem. A, Vol. 110, No. 16, 2006 5529 (10) Han, J.-X.; Hung, P. Y.; DeSain, J.; Jones, W. E.; Curl, R. F. J. Mol. Spectrosc. 1999, 198, 421. (11) Allen, M. D.; Evenson, K. M.; Brown, J. M. J. Mol. Spectrosc. 2001, 209, 143. (12) Hung, P. Y.; Sun, F.; Hunt, N. T.; Burns, L. A.; Curl, R. F. J. Chem. Phys. 2001, 115, 9331. (13) Jones, W. E.; Sun, F.; Curl, R. F.; Evenson, K. M.; Brown, J. M. Can. J. Phys. 2001, 79, 389. (14) Nimlos, M. R.; Davico, G.; Geise, C. M.; Wenthold, P. G.; Lineberger, W. C.; Blanksby, S. J.; Hadad, C. M.; Petersson, G. A.; Ellison, G. B. J. Chem. Phys. 2002, 117, 4323. (15) Zandler, M. A.; Goddard, J. D.; Schaefer, H. F., III. J. Am. Chem. Soc. 1979, 101, 1072. (16) Kim, K. S.; Schaefer, H. F., III; Radom, L.; Pople, J. A.; Binkley, J. S. J. Am. Chem. Soc. 1983, 105, 4148. (17) Rice, J. E.; Schaefer, H. F., III. J. Chem. Phys. 1987, 86, 7051. (18) Malmquist, P. A.; Lindh, R.; Roos, B. O.; Roos, S. Theo. Chim. Acta 1988, 73, 155. (19) Seidl, E. T.; Schaefer, H. F., III. J. Chem. Phys. 1992, 96, 4449. (20) Francisco, J. S. Chem. Phys. Lett. 1994, 230, 372. (21) Koput, J. J. Phys. Chem. A 2002, 106, 6183. (22) Koput, J. J. Phys. Chem. A 2003, 107, 4717. (23) Bunker, P. R.; Landsberg, B. M. J. Mol. Spectrosc. 1977, 67, 374. (24) Gdanitz, R. J.; Ahlrichs, R. Chem. Phys. Lett. 1988, 143, 413. (25) Cave, R. J.; Davidson, E. R. J. Chem. Phys. 1988, 89, 6798. (26) Szalay, P. G.; Bartlett, R. J. Chem. Phys. Lett. 1993, 214, 481. (27) Werner, H.-J.; Knowles, P. J. Theo. Chim. Acta 1990, 78, 175. (28) MOLPRO is a package of ab initio programs designed by H.-J. Werner and P. J. Knowles. Amos, R. D.; Bernhardsson, A.; Berning, A.; Celani, P.; Cooper, D. L.; Deegan, M. J. O.; Dobbyn, A. J.; Eckert, F.; Hampel, C.; Hetzer, G.; Knowles, P. J.; Korona, T.; Lindh, R.; Lloyd, A. W.; McNicholas, S. J.; Manby, F. R.; Meyer, W.; Mura, M. E.; Nicklass, A.; Palmieri, P.; Pitzer, R.; Rauhut, G.; Schu¨tz, M.; Schumann, U.; Stoll, H.; Stone, A. J.; Tarroni, R.; Thorsteinsson, T.; Werner, H.-J. (29) Meyer, W.; Botschwina, P.; Burton, P. J. Chem. Phys. 1986, 84, 891. (30) Mladenovic´, M. J. Chem. Phys. 2003, 119, 11513. (31) Mladenovic´, M. J. Chem. Phys. 2000, 112, 1070. (32) Arfken, G. Mathematical Methods for Physicists, 3rd ed.; Academic Press: New York, 1985. (33) Zare, R. N. Angular Momentum; J. Wiley & Sons: New York, 1988. (34) Mladenovic´, M. In NIC Symposium 2001; Rollnik, H., Wolf, D., Eds.; NIC Series 9; John von Neumann Institite for Computing: Ju¨lich, Germany, 2002; pp 85-94. (35) Mladenovic´, M. Spectrochim. Acta, Part A 2002, 58, 795. (36) Echave, J.; Clary, D. C. Chem. Phys. Lett. 1992, 190, 225. (37) Yamada, K.; Winnewisser, M. Z. Naturforsch. 1976, 31A, 139.