THE TORSIONAL SPECTRUM OF THE HYDROGEN TRIOXIDE MOLECULE G.A.Pitsevich, A.E.Malevich, U.U.Sapeshka

The energies of the stationary torsional levels of the hydrogen trioxide molecule were calculated at the B3LYP, MP2, and CCSD(T) levels of theory using augmented correlation consistent acc-pVTZ basis set. The molecular symmetry group whose elements are inherent to both equilibrium conformer's symmetry elements (C2 and CS) of the HT molecule was found. Different methods of the molecular parameters calculations were suggested and analyzed. The torsional, spin and total wave functions were classified by irreducible representations of the C2V(M) molecular symmetry group. The 2D dipole moment surface was calculated too. The energies of the stationary torsional states were found using DVR and Fourier methods. With this and acceptable combinations of spatial and spin wave functions, the IR torsional spectrum was calculated at different temperatures. The tunneling frequencies in the ground and some excited torsional states were estimated too.


INTRODUCTION
The hydrogen trioxide (HT) molecule HООOH in many ways is a benchmark object.
Firstly, this molecule is the second after hydrogen hydroxide representative of the series of hydrogen polyoxides specified by the formula НОnH (n=1,2,3,4,5….), the researchers' interest to which is undoubtedly increasing in recent times [1][2][3][4]. Secondly, the НТ molecule is the simplest representative of a whole series of molecules having two equivalent non-coaxial internal tops.
Thirdly, the molecule is capable of forming the cluster structures due to the formation of hydrogen bonds. Besides, the НТ molecule has two equilibrium configurations, which differ in the mutual orientation of hydroxyl groups. Finally, the НТ is potentially non-rigid molecule since each of the two conformers can exist in two configuration equivalent versions, even though, the matter of the real rigidity or non-rigidity has yet to be analyzed. Many researchers also point out the significant role of the НТ molecule in oxidizing processes that are taking place in atmospheric phenomena [3,5,6] and biological systems [3,7,8], as well as in other significant chemical reactions [9][10][11][12][13][14][15]. In addition, hydrogen and oxygen atoms are included in a set of 11 2 elements of the "astronomical periodic table of elements" [16]. Consequently, an HT molecule can be of interest in astronomical applications.
The low stability of the molecule being studied makes the performance of the spectral investigations considerably more difficult. The first success in registering the IR and Raman spectra of the HT molecule belongs to a group of scientists under the leadership of Giguėre in the early 1970s [17][18][19]. However, the presence of a number of other polyoxides in the samples made it possible to make only some preliminary assignments in the experimental spectra. The detailed interpretation of the IR spectrum of HT registered in argon matrix was proposed by Engdahl and Nelander in 2002 [20], and in 2005, Suma registered a purely rotational spectrum of the HT molecule [21].
The theoretical investigations of the potential energy surface (PES) and spectra of the HT molecule were started simultaneously with the beginning of the spectral investigations [22][23][24].
In 1978, Cremer examined in detail how the basic set choice and the way of accounting of the electron correlation affects the structural and energetic parameters of the molecule [25]. On the PES, the global and local minima, transition states were revealed as well as the characteristics of НТ were compared with those of hydroxide peroxide [26]. The investigations of the PES were continued in the papers [27][28][29], while the mechanisms [30] and enthalpies [31,32] of the formation of the HT molecule and methods of its synthesis [33,34] were determined by quantumchemical methods in a number of papers.
As we know, the calculation of the normal vibrations of the HT molecule was firstly performed in paper [19]. Later the calculations of the IR spectra of the molecule using the quantum-chemical methods were continued in the papers [29,35,36]. In the first of these papers, the structures and IR spectra were calculated using the CASSCF method and augmented by diffuse and second sets of polarization functions basis set. In the paper [36], the respective calculations were performed on the CCSD(T)/cc-pVXT (X=T,Q,5,6) levels of the theory although the anharmonic calculations of the IR spectrum were performed using the DFT methods. In the paper [29], the calculations of the structure and IR spectra of the ТН were performed using the explicitly correlated techniques (CCSD(T)-F12/cc-pVTZ-F12) taking into account the anharmonicity effects. Many interesting aspects concerning the structure, spectral as well as physical and chemical characteristics of the HT molecule can be extracted from recent reviews [37,38].
As follows from the review presented, the researchers' interest in the HT molecule does not wane. Many spectral and structural characteristics of the molecule are becoming more accurate as new improved methods in quantum chemistry appear. However, one of the most interesting and important aspects, which is the internal rotation of the two hydroxyl groups, has not yet been 3 studied. The question of whether the HT molecule is rigid or nonrigid object has yet to be answered. These problems will be the main objective of this paper. To give the answers to these questions, first of all, in section 2 we analyze the symmetry properties of the molecule. Then in Section 3, we will discuss the mathematical backgrounds of calculations and different methods of data preparation. In Section 4 we compare the results obtained using DVR and Fourier methods for the numerical solution of the Schrödinger equation. In Section 5 we discuss the results of torsional levels energy calculations and simulate the torsional IR spectra of both conformers of the HT molecule. In Section 6 some conclusions will drow.

MOLECULAR SYMMETRY GROUP
The rigidity or nonrigidity of the HT molecule is determined by its capability of tunneling from one configuration to another, which depends on the height and shape of the potential barrier separating equivalent configurations. If there is no tunneling, the order of a group of the molecular symmetry [39,40] HT is equal to that of its point symmetry group (С2) [41]. If the tunneling takes place, the order of HT molecular symmetry group can be determined using the results of the paper [42], according to which the order of the molecular symmetry group is equal to the product of the order of the point symmetry group by the number of the configuration versions of the nonrigid molecule. Since the group order С2 is equal to two and there are two HT versions, the order of the molecular symmetry group is four. One can find that the molecular symmetry group of the HT molecule contains four symmetry elements: E, (12)(34), E * , (12)(34) * and, in accordance with [40], this molecular symmetry group may be called C2V(M). Note that all these symmetry elements are inherent to both equilibrium conformers symmetry elements (C2 on the left in Fig.1 and CS in the middle in Fig.1) of the HT molecule.   P1  P2  P3  P4  A1  1  1  1  1  A2  1  1  -1  -1  B1  1  -1  1  -1  B2  1  -1  -1  1 The symmetry properties of the 2D surfaces of the potential energy, dipole momentum and kinematic coefficients of the HT molecule are very close to those of the methanediol molecule (MD), for which the internal rotation was recently analyzed [43]. In [43] the symmetry elements of the 2D surfaces were connected with symmetry elements of the C2V(M) molecular symmetry group. According to [43] let's call the plane which is perpendicular to the coordinate plane ,  and intersect it along diagonal (0 0 ,0 0 -360 0 ,360 0 ) plane 1. The plane which is perpendicular to the coordinate plane ,  and intersects it along diagonal (0 0 ,360 0 -360 0 ,0 0 ) will be plane 2. The intersection of these two planes will be the C2 symmetry axis. Then these elements will be connected with symmetry elements of the C2V(M) group according to data from  Table 2 too.
The results of the transformations of some molecular characteristics under C2V(M) symmetry operations are collected in Table 3. Table 3 The symmetry species and results of the transformations of some HT molecular characteristics (see (6) and (15)) under C2V(M) molecular symmetry group.

Symmetry operators
Irreducible

CALCULATION DETAILS
The equilibrium configurations HT are presented in Fig.1 at the global (trans-conformer) and local (cis-conformer) energy minima as well as in the configuration being the initial one for performing the calculation of the 2D potential energy surface (PES). The hydroxyl groups rotate around the O-O bonds. Further, we will be denoting the equivalent torsional coordinates as  and  . The zero values  and  were chosen in the configuration presented on the right of Fig.1 where the hydroxyl atoms of hydrogen (Н1 and Н2) lie in the plane, which is formed by the О3-O5-О4 atoms.
The potential energy satisfies the condition (see Table 2 and 3):  Using symmetry operation P2 and then P4 it is easy to obtain values of all molecular parameters in full square. Let's it will be the 1-st method of data preparation. For the reasons stated below sometimes we also calculated the data in the full square shown in Fig.2. In this 6 case, we used this data as they are (the 2-nd method, when symmetry conditions are not satisfied), or averaged them by four equivalent points (the 3-d method). Preparing the data by the 1-st method the torsional coordinates satisfied the conditions: Preparing the data by the 2 and 3 methods the torsional coordinates satisfied the conditions: The values of the energy of the torsional states were calculated using the DVR and Fourier . Let's it will be the 4-th method of data preparation (when interpolation is used). To check interpolated values quality the less dense equidistant 2D K*K grid (K=15) was extracted from 2D 45*45 grid calculated data for symmetrized methods 1 and 3, where torsional coordinates were varied from 12 0 to 348 0 in step 24 0 . This data was used to obtain again equidistant 2D K*K grid (K=45) by interpolation using the same method [44].
In the nodes of the equidistant 2D K*K grid (K=45), keeping fixed values of the torsional coordinates, HT geometry was optimized on all other coordinates using MP2 [45,46] and CCSD(T) [47,48] levels of theory and acc-pVTZ [49][50][51] basis sets in frame of quantumchemical package [52]. Configurations optimized in each node were used to calculate the kinematic coefficients using Wilson vectors ( s -vectors) [53]. The vibrational Schrödinger equation with reduced dimensionality for the torsional vibrations of two hydroxyl groups can be written as follows [43,[54][55][56]: The kinematic parameters can be determined from the following relations [43,[54][55][56]: To solve the Schrödinger equation (6) by the DVR method [57][58][59][60][61][62] the Hamiltonian matrix was calculated by the formula: Here: To solve the Schrödinger equation (6) where , , coefficients were found by the fitting coefficients in the equation (6) by 2D complex Fourier series using [44]: The wave function is sought in the form: Since the Cartesian coordinate system varied its orientation in relation to the molecule as the hydroxyl groups rotate, components of dipole moment vector at each node of the 2D grid were recalculated for the Cartesian coordinate system related to the molecular skeleton, in which the Z-axis was directed along О3-O5-О4 angle bisector, the X-axis lied in the plane of the О3-O5-О4 atoms and the Y-axis supplemented the X and Y to the right-hand triple and was, therefore, normal to the plane of three heavy atoms. 8 Squares of matrix elements of the dipole moments operator were found using [44] and next formulae: (15) Here ,, x y z p p p are the dipole moment vector components in the molecular fixed Cartesian are the torsional wave functions in the initial (i) and final (f) states.

AND FOURIER METHODS
The calculated using DVR method values of the torsional stationary states energies for the data sets, obtained by 1-4 methods on the grids with different densities are represented in Table 4  Before proceeding to the discussion of the data in Table 4 it is important to stress that due to energy optimization in every node of the grid the calculated values of the potential energies and kinematic coefficients in four symmetrical points, presented in Fig.2, are not equal to each other. This is due to the fact that the finish of the calculations occurs when the conditions that satisfy the convergence criteria are reached. Moreover, at different equivalent points, these 9 conditions can be achieved for several different configurations of molecules with slightly different energies and geometric parameters. In spite of these differences are very small they can affect the accuracy of some results of calculations. On another hand, it is clear that using method 3 one can get by averaging the data at four equivalent points more accurate values of the energy and kinematic parameters that in addition satisfy the symmetry conditions.
Turning to the analysis of the data given in Table 4 Table 4) and for different M and N (see formulas (11) and (12)) are represented in Table 5. As one can see from Table 5 all values but tunneling frequency in the ground (and exited) states are nearly insensitive to M and N. The differences in the values of the latter are five orders of magnitude. Comparing the data in Table 5 and data in the 4-th column of  Table 4. This data collected in Table 6. Comparing data in Table 4 and Table 6 it is possible to state that results of calculations of the tunneling frequency in the ground state using DVR and Fourier methods (in the last case   /2 MK  ) become more similar. It should be emphasized that when we perform calculations on the 2D equidistant K * K grid (K = 135), the results of calculation of the tunneling frequency in the ground state using the Fourier method are weakly depend on M and N values. In addition, since results represented in columns 6 and 7 of Table 6 are very similar, we performed calculations of the energies and kinematic coefficients at CCSD(T)/acc-pVTZ level of theory only in triangle represented in Fig. 2 and then spread over the whole square using the 1-st method.

DISCUSSION OF THE RESULTS.
In Table 7 we collected data on the parameters of equilibrium trans-and cis-configurations of the HT molecule calculated at MP2/acc-pVTZ and CCSD(T)/acc-pVTZ levels of theory. In addition, based on literature data, we included in Table 7 the calculated geometrical parameters for the trans-conformer obtained at more advanced levels of theory.  [29] According to represented in Table 7    The most accurate from our point of view values of energies of the stationary torsional states measured from the ground torsional state of the HT molecule are summarized in Table 8.
As the results, represented in Tabl Table 8).  Table 8.
We were also able to determine torsional quantum numbers for several lowest torsional states.
The as n and s n in columns 5 and 6 of Table 8 are antisymmetric and symmetric torsional quantum numbers in trans-and cis-conformers (see the Configuration column in Table 8).
Symmetry species of the torsional states analyzed with the help of the data from Tables 1 and 2, are represented in columns 5 and 10. The tunneling splittings of the torsional levels are indicated in column 11 of Table 8. As expected the splitting of the torsional levels increases with the increasing their energy. To calculate the intensities of the torsional IR spectrum, it is also necessary to consider Thereby, the torsional wave functions must be combined with the spin wave functions as it is shown in Table 9. As one can see from Table 9, the torsional wave functions that belong to A1 and A2 symmetry species are nondegenerate, while the torsional wave functions that belong to B1 and B2 symmetry 16 species are triple degenerate. In addition, the symmetry of the spin wave functions must be saved during the torsional transitions. Relative intensities of the IR bands for torsional transitions from the initial state i to the final state f were calculated using the next formulas [43,40,67- (17) where i g is the multiplicity of degeneracy of the initial torsional state, the value of which is coordinated with the type of symmetry of the torsional wave function in accordance with the data of Table 8 and 9, i s and f s are the specifiers of the initial and final states, the value of which is 1 for the torsional states of the А1 and А2 type and -1 for the torsional states of the В1 and В2 as shown in Table 9, () QT -partition function. With the use of these data, the relative intensities of the IR absorption bands at different temperatures were calculated from the formula (17). Fig. 9 presents the torsional IR spectrum of HT molecule at the temperatures of 300 and 30 K. It is clear that the HT cis-conformer is less stable than the trans-conformer. However, the potential barrier for cis  trans transition is higher than 800 cm -1 and the matrix insolation and low temperatures can stabilize cis-configuration. That is why the torsional spectrum of the cisconformer of the HT molecule was also calculated at the temperature of 300 and 30 K. In this case, in the formula (17) when calculating the partition function, the summation was performed from i = 13 to i = 50 with the step of 1. It is presented in Fig. 10  In addition as one can see from Table 8, the fifth torsional state is metastable due to all four lower torsional states have different symmetry of the spin functions. With this, the torsional spectrum of the trans-conformer of HT was calculated for which the 5th torsional level is ground at a temperature of 300 and 30 K. In this case, in the formula (17) when calculating the partition function, the summation was performed from i = 5 to i = 50 with the step of 1. It is presented in Fig. 11. In the calculated torsional spectrum of the trans-conformer which is stuck in the fifth excited torsional state at 30 K, there are only two absorption bands. The intensive IR band has a maximum near 330.7 cm -1 while weak one has a maximum near 746.6 cm -1 .
Let us compare calculated and experimental (in matrix isolation [20]) IR spectra of the HT molecule. It is obvious that the calculated values of frequencies of the most intense absorption bands in the torsional IR spectrum of the HT at low temperatures turn out to be below 18 500 cm -1 . Unfortunately, the paper [20] presents the torsional spectrum of the HT molecule in the matrix isolation within the interval 740-800 cm -1 only. According to [20] 346.4 and 387.0 cm -1 IR bands were assigned to symmetric and antisymmetric torsional vibrations of the hydroxyl groups, while 776.1 cm -1 band was assigned to antisymmetric stretching O-O vibration. In addition, in the IR spectrum of the HT molecule represented in [20], one can see the broad band with a maximum near 752 cm -1 which was not assigned by authors. According to the above, in the torsional IR spectrum of the ground state trans-conformer at 30 K, there are no absorption bands in the spectral region higher 400 cm -1 . Thus the 752 cm -1 IR band can indicate the presence in sample cis-conformer or trans-conformer in fifth excited torsional state and may be assigned to the one of three mentioned above (740.5, 746.6 and 786.8 cm -1 ) torsional vibrations in these conformers. Furthermore, since in the first four torsional states the HT molecule has to be the para-spin isomer while in fifth and six torsional states the HT molecule has to be the ortho-spin isomer the excitation of the antisymmetric torsional vibration from the two lowest states is forbidden. That means that assignments of the 387.0 cm -1 band to the antisymmetric torsional vibration is questionable. According to our calculations, this band may be assigned to the symmetric torsional vibration (400.7 cm -1 ) in the cis-conformer, however, it is obvious that additional investigation is required.

CONCLUSIONS
The calculation of the torsional IR spectrum of the HT molecule has been performed at several levels of theory. It is quite expected that the calculated energy of the torsion levels depends on the level of theory used. However, as can be seen from Table 7 (8), when calculating at various levels of the theory, the relative location of different torsion energy levels also changes. Performing calculations in frame of the DVR method, it is important correctly to choose the position of points of a 2D equidistant grid, taking into account the symmetry of the potential energy surfaces and kinematic coefficients. Accounting for the C2V(M) symmetry of the HT molecule made it possible to propose a minimum set of points on a 2D equidistant grid in which the molecular parameters should be calculated (see Fig. 2). An important parameter in calculations using the DVR method is the density of 2D grid points. In particular, it is shown that the use of 2D 15*15 mesh leads to unsatisfactory results (see column 1 in Table 4). At the same time, the Fourier method allows even in this case to obtain reasonable values of the energies of the torsion states (see column 1 in Table 5). Improving the accuracy of calculating the energy values of torsion states can be achieved by calculating the molecular parameters on a denser 2D equidistant grid (45 * 45 and 135 * 135). At the same time, as shown in the work, the correct use of interpolation methods that take into account the periodicity of molecular parameters and their 19 symmetry properties makes it possible to increase the accuracy of calculations of the energy values of torsion levels without increasing computational costs.
As one can see the splitting due to tunneling of the torsional levels, especially in the case of the ground state, is very sensitive to the quality of prepared data on 2D PES and 2D surfaces of the kinematic parameters. They have to satisfy the symmetry conditions and in addition, must be as smooth as possible. It should be noted that during data symmetrization and during the solving the Schrödinger equation using the Fourier method all 2D surfaces are smoothed. Due to the last circumstance, under equal conditions, the Fourier method shows better results in comparison with the DVR method.
The classification of the torsional and spin states of the HT molecule by the irreducible representation of the С2V(M) molecular symmetry group has been performed. The torsional IR spectra of the HT molecule at the temperatures of 300 and 30 K have been calculated. The preliminary comparative analysis of the calculated and experimental IR spectra of the НТ at the temperature of 30 K has been performed. This analysis allows us to assume that the argon matrix stabilizes the cis-conformer of the molecule and the latter can be presented in the sample.
Based on the results of the calculations presented in Section 4, the calculated values of the tunneling frequencies should be treated with caution if their values are in the order of 10 -9 -10 -10 cm -1 . It is obvious that the accuracy of calculations of the tunneling frequencies is directly related to the accuracy of calculations of the energy values of torsion states. In this case, it is also important to take into account the accuracy of calculations of energy and geometrical parameters at the nodes of uniform 2D grids using quantum chemical packages. In particular, tightening the convergence criteria when optimizing the geometry of a molecule should lead to an increase in the accuracy of calculations of the values of the tunneling frequencies. The results of this paper show that taking into account molecular symmetry in determining the values of molecular parameters is important in improving the accuracy of calculating the tunneling frequencies. The calculation of molecular parameters on the full 2D grid followed by averaging the values over four equivalent points on the 2D coordinate plane (see Fig. 2) should also lead to an increase in the accuracy of tunneling frequency calculations. Achieving agreement on the calculated values of the tunneling frequencies obtained using DVR and Fourier methods can serve as a good criterion for the reliability of the results obtained.
Based on the data presented in Table 4, and based on the analysis of the results of Section 4, it can be assumed that the tunneling frequency in the ground state is in the range of 10 -10 -10 -9 cm -1 or 3-30 Hz. According to the data of calculations at the CCSD (T) / acc-pVTZ level of the theory, the tunneling frequency in the cis-conformer is several orders of magnitude greater than in the trans-conformer and turns out to be 1.047 * 10 -6 cm -1 or 31388 Hz. Finally, the frequency 20 of tunneling in the metastable excited antisymmetric torsion vibrational state of the transconformer (see Table 8) is comparable in magnitude with the frequency of tunneling in the ground state and presumably must fall in the frequency range 5.5 -55 Hz. All tunneling frequencies as well as few allowed for trans-and cis-conformers torsional vibrations may serve as a reference value when searching the HT in the cosmic space and in comets where the cisconformer can theoretically exist at low temperatures too.