Theoretical Investigation on H2O2-Ng (He, Ne, Ar, Kr, Xe, and Rn) Complexes Suitable for Stereodynamics: Interactions and Thermal Chiral Rate Consequences

Although molecular collisions of noble gases (Ng) can be theoretically used to distinguish between the enantiomers of hydrogen peroxide - H2O2 (HP), little is known about the effects of HP-Ng interactions on the chiral rate. In this work, the chiral rate as a function of temperature (CRT) between enantiomeric conformations of HP and Ng (Ng=He, Ne, Ar, Kr, Xe, and Rn) are presented at MP2(full)/aug-cc-pVTZ level of theory through a fully basis set superposition error (BSSE) corrected potential energy surface. The results show that: (a) the CRT is highly affected even at a small decrease in the height of trans-barrier; (b) its smallest values occur with Ne for all temperatures between 100 and 4,000 K; (c) that the decrease of CRT shows an inverse correlation with respect to the average valence electron energy of the Ng and (d) Ne and He may be the noble gases more suitable for study the oriented collision dynamics of HP. In addition to binding energies, the electron density ρ and its Laplacian ∇2ρ topological analyses were also performed within the atoms in molecules (AIM) theory in order to determine the nature of the HP-Ng interactions. The results of this work provide a more complete foundation on experiments to study HP's chirality using Ng in crossed molecular beams without a light source.


INTRODUCTION
Hydrogen peroxide -H 2 O 2 (HP) is a molecule of interest in a large and diverse number of fields in addition to its industrial uses. For example, it has emerged as a major metabolite in redox signaling and regulation (Antunes and Brito, 2017;Sies, 2017), and its presence was observed in Martian atmosphere (Encrenaz et al., 2004) and also on the surface of Jupiter's moon Europa (Carlson et al., 1999). The HP is interesting since it is simplest molecule that exhibits internal (torsional) rotation and chirality. Furthermore, this molecule can form dimers (Dobado and Molina, 1993;González et al., 1997), clusters (Yu and Yang, 2011), complexes with water (Mo et al., 1994;González et al., 1997) and with biologically important molecules such as adenine (Dobado and Molina, 1999), DNA (Piatnytskyi et al., 2016), glycine (Shi and Zhou, 2004) or nitrosamines (Roohi et al., 2010).
On the other hand, the hydrogen peroxide seems to be a prototypical model to be used into experiments to observe chirality in crossed molecular beam without a light source , a frontier in research of stereodynamics which is still at early stages (Su et al., 2013;Lombardi and Palazzetti, 2018). In these kind of experiments, the molecular orientation control on the intense continuous beam is mandatory to the phenomena of chiral selectivity to be demonstrated (Aquilanti et al., 2005). For this reason, the interaction between HP and atoms, molecules and ions is so relevant to sterodynamics studies (Barreto et al., 2007(Barreto et al., , 2010Lombardi et al., 2011;Roncaratti et al., 2014;Leal et al., 2016).
In the present paper we investigated the dynamics of the chiral molecule HP interacting by van der Waals forces with noble gases Ng (Ng=He, Ne, Ar, Kr, Xe and Rn) in order to obtain the chiral rate as a function of temperature (CRT) for these complexes, as well as its consequences for the chiral kinetic interconversion when tunneling effect is included. The aim of this work was to understand how the kinetic interconversion of the two HP's isomers is affected along collisional events through a basis set superposition error (Boys and Bernardi, 1970) (BSSE) corrected potential energy surface (PES). In addition, the electron density ρ(r) and its Laplacian ∇ 2 ρ(r) topological analyses were also performed within the atoms in molecules (AIM) theory in order to determine the nature of the intermolecular interactions.

Computational Details
All calculations were performed using the Gaussian09 package (Frisch et al., 2009). The structures were optimized without constraints at MP2(full) method in conjunction with aug-cc-pVTZ-PP for Xe and Rn (Peterson et al., 2003) and aug-cc-pVTZ for the remaining atoms (Dunning, 1989;Woon and Dunning Jr, 1993;Wilson et al., 1999). Vibrational frequencies at the same level of theory were also performed in order to ensure that each minimum has only positive frequencies and that each transition state has only a single imaginary frequency, as well as to obtain the zero point vibrational energy (ZPE). The counterpoise method of Boys and Bernardi (1970) was used to correct the BSSE for binding energy.

Overview of the Transition State Theory
The transition state theory (TST) (Truhlar et al., 1996) was developed primarily by Henry Eyring (Eyring, 1935) and Michael Polanyi (Polanyi and Wigner, 1928) between 1928and 1935. The TST is an improvement over the so-called theory of collisions (Lewis, 1967), and it is widely used to calculate the rate constants of chemical reactions.
The start point of TST is the existence of a transition state (TS) between the reagents and products. Located at the top of the potential energy barrier and it assumes a quasi-equilibrium between reactants and activated transition state complexes. For a bimolecular reaction given by (1) The TS is characterized by a single imaginary frequency along the reaction coordinate of the molecular system which is represented here byν 1 . In its turn, the reaction coordinate can be represented by angular changes in bond distances during the chemical reaction (Henkelman et al., 2002). The equation that determines the reaction rate is known as the Eyring equation, given by where 0 < κ(T) ≤ 1 is the so-called transmission coefficient, k B is the Boltzmann constant, h is the Planck constant, q • m is the standard molar partition function, N A is the Avogadro constant, R is the gas constant and E 0 b is the barrier energy with zero-point energy correction. In addition, the TS, R 1 and R 2 subscripts stand for the transition state and reagents, respectively. Thus, the rate constant is determined by the parameters that characterize both reagents and the TS.
The general partition function is formed by the product of translational q trans. , rotational q rot. , vibrational q vib. and electronic q ele. partition functions. The translational partition function for a free particle with mass m moving along the length dimension l x can be evaluated by considering that the separation of energy levels is small and that a large number of states are accessible at room temperatures. Therefore, the energy levels should be continuous and the sum contribution of the translational partition function becomes an integral. Which the solution for the three-dimensional case is (Atkins et al., 2013) Although the system can be excited at normal modes, the energy levels are discrete for the rotational mode. The three degrees of freedom of spatial rotation and the three moments of inertia I A , I B and I C must be taken into account for a non-linear molecule (Atkins et al., 2013), thus where σ is the so-called number of symmetry. The vibrational mode has reasonably spaced energy which must be taken into account since they are partially occupied. As a consequence, the vibrational partition function is strictly calculated as a sum over the occupied states. In the case of n vibrational degrees of freedom, the vibrational partition function is given by the product of n partition functions, where ν i is each of the fundamental vibrational frequencies. In most cases, only the lowest energy state is occupied and the electronic energies should not contribute considerably to the total partition function (Atkins et al., 2013). A good approximation is to disregard the contributions of the nuclear and electronic spins and to vanish the fundamental energy level for the electronic partition function. Under these considerations the electronic partition function should be equal to unity (Atkins et al., 2013) q ele. = 1.
On the other hand, the coefficient κ(T) represents the tunneling effect of the reaction coordinate of the chemical system and it is usually important for light atoms or molecules at low temperatures. Thus, tunneling estimates were made using both Wigner (Polanyi and Wigner, 1928) and Eckart (Eckart, 1930) methods. The Wigner tunneling correction proposes a parabolic potential, where E b corresponds the energy potential barrier of MEP,ν 1 is the imaginary frequency of transition state and s is the coordinate reaction. This implies in a transmission coefficient given by Bell (1959) κ Wigner (T) = 1 − 1 24 For very low temperatures, the Wigner tunneling effect is not very effective, and for this reason, it was also employed Eckart tunneling correction (Truhlar et al., 1985). The Eckart tunneling correction uses a potential of the type where α is a parameter described by and µ is the reduced mass of the system. These parameters determine the barrier width. Here it is important to note that the A and B can be positive, negative or zero. The A parameter corresponds to the energy difference V(x → −∞) and V(x → +∞), and B is a parameter that measures the height of the barrier given by So the most usual form for the Eckart's potential in the study of reaction rates is (Truhlar et al., 1985) where s is the coordinate of the reaction and s 0 is the reaction coordinate corresponding to the maximum of the barrier, which is given by Finally, the transmission probability (Bell, 1980), obtained through the solution of the Schrödinger equation with Eckart's potential, is expressed by the following equation where k, β and δ depend onν 1 , A, B and energy (E). The quantum tunneling correction κ(T) can thus be calculated from the ratio between the quantum rate k quan. (T) and the classical rate k class. (T) in which the particles cross the barrier. Thus, the Eckart tunneling correction with transmission coefficient is given by (15) where integration is performed over all possible energies.

Geometric Parameters, Interactions and AIM Analysis
The details about the generation of the potential energy surface are described in another work of our group (Roncaratti et al., 2014), so it will be commented briefly here. First, all HP geometry  Topological studies performed on this adjusted potential energy surface (PES) showed that the HP and HP-Ng complexes have two overall minimum configurations, termed cis (labeled as θ − ) and trans (labeled as θ + ), separated by two potential barriers, denoted here as cis-barrier and trans-barrier. The potential energy curves (PEC) obtained from the PES are then presented in Figure 2.
The PECs obtained for HP-Ng complexes are similar in shape and depth. The cis-barriers for the HP-Ng complexes are all smaller than the respective value for the free HP. The free HP has a cis-barrier of 7.5594 kcal/mol whereas the values for the complexes increase monotonically from 6.9828 kcal/mol for HP-Rn up to 7.5107 kcal/mol for HP-He. In addition the trans-barrier values are also lower than the respective value for the free HP, which is 1.0427 kcal/mol, and their values are 1.0928 kcal/mol for HP-He, 1.0817 kcal/mol for HP-Ne, 1.0651 kcal/mol for HP-Ar, 1.0676 kcal/mol for HP-Kr, 1.0736 kcal/mol for HP-Rn and 1.0749 kcal/mol for HP-Xe complexes. For the free HP, the cis-barrier and trans-barrier experimental energies (Hunt et al., 1965) are 7.0334 ± 0.0715 kcal/mol and 1.1036 ± 0.0114 kcal/mol, respectively. These values are in a good agreement with our results. However, for TABLE 3 | Bond critical point (BCP) data for charge density ρ (in ×10 −3 e/a 3 0 ), Laplacian of the charge density ∇ 2 ρ (in ×10 −2 e/a 5 0 ), electronic energy density H(r) and ellipticity ε for configurations 1(cis), 2(cis-barrier), 3(trans) and 4(trans-barrier) of the HP-Ng complexes. the HP-Ng complexes we did not find experimental data for comparison. The results concerning geometric parameters, interactions and their characterization are summarized in Tables 1-3. The geometrical parameters obtained at MP2(full)/aug-cc-pVTZ, optimized without any constraints for two minimum structures and transition states (presented as cis and trans barriers) are given in Table 1 together with the graphical representation in Figure 3. Table 2 lists the binding energies corrected and uncorrected for the BSSE. Table 3 shows the numerical results for AIM analysis and Figure 4 depicts the ∇ 2 ρ(r) contour plots for cis and trans barrier configurations.

Complexes Configuration Description
The PESs yield θ − and θ + as true minima, i.e. without any imaginary frequencies in accordance with results from literature (Maciel et al., 2006;Roncaratti et al., 2014). In addition, all transition state structures displayed a well characterized imaginary frequency around 600 cm −1 for cis and 400 cm −1 for trans barriers (see Supplementary Information for further details). Figure 5 describes a schematic representation of vibrational modes of the isolated HP in the transition state with the actual frequencies and an imaginary FIGURE 3 | Geometric representation of the cis-barrier, θ − , trans-barrier and θ + structures. frequency, which represents the frequency along the reaction coordinate.
For the HP-Ng complexes, the geometrical parameters are almost the same when compared with isolated HP in agreement with the weak interaction of these systems. The HP-Ng distances increase from He up Rn. On average, they are close to 2.55Å(He), 2.65Å(Ne), 2.75Å(Ar), 2.85Å(Kr), 3.00Å(Xe), and 3.06Å(Rn).
Regarding the binding energies, HP-He and HP-Ne are all repulsive, being less repulsive for the cis barrier configuration. This can be understood as a consequence of the fact that the noble gases turn out to be the hardest elements (Furtado et al., 2015) and this hardness decreases when the Ng atomic number is increased (the hardness in this context is a resistance to changes in its electronic population Furtado et al., 2015 coupled to Ng's high electronegativity Allen and Huheey, 1980). Although the BSSE increases monotonically from He to Rn, the binding FIGURE 4 | Laplacian of the electron density, ∇ 2 ρ, contours map in the molecular plane for cis-barrier (top) and trans-barrier (bottom) conformations of HP-Ne at MP2(full)/aug-cc-pVTZ level. The thick solid lines represent the molecular graph that joints the nuclei, the bond critical point and ring critical point, and also represent the zero flux surface. energies also become more attractive as the atomic number increases.
For the four structures of each HP-Ng, the higher binding values are always observed for the cis-barrier configuration. As it will pointed latter, the decrease of the rate through the two barriers are not correlated with the binding energy, suggesting the hyperconjugation effects on HP may be important for the decrease of the interaction rate.
Regarding the AIM analysis, the existence of (3,−1) bond critical point (BCP) and its associated atomic interaction line indicates that electronic charge density is accumulated between the linked nuclei (Bader, 1991). In its turn, the values of the charge density ρ(r) in BCP are small while their corresponding ∇ 2 ρ(r) are positive in accordance with a closed shell type of interaction. As a consequence, all configurations of all complexes show an interaction of a van der Waals type. Since higher ellipticity suggests conjugation and hyperconjugation effects of electron delocalization, these effects seem more pronounced in the HP-He and HP-Ne complexes. Another interesting feature is that all cis-barrier configurations of all complexes show a (3,+1) BCP indicating a cyclic nature.

Thermal Chiral Rate Analysis
The temperature dependence of the rate constant for cis to trans (i.e., through trans-barrier) and trans to cis (i.e., through cis-barrier) conformations for HP and HP-Ng complexes are presented in Figure 6. These results, in addition to conventional rate, are also exhibited with Eckart's or Wigner's tunneling corrections.
It was found that for the entire 100 K up to 4,000 K range the HP-Ne has the lowest rates for both barriers among all noble gas complexes, followed by HP-He. This result suggests that Ne and He are the noble gases more suitable for study the oriented FIGURE 5 | Schematic representation of vibrational modes of the H 2 O 2 isolated in the transition state. There are 6 vibrational modes, 5 of them correspond to the actual frequencies (ν 2 , ν 3 , ν 4 , ν 5 , and ν 6 ) and one of them corresponds to the imaginary frequencyν 1 , which represents the frequency along the reaction coordinate. collision dynamics with HP. In fact, the decrease of CRT shows an inverse correlation with respect the average valence electron energy (Allen, 1989), which follows the sequence (from higher to lower values): Ne, He, Ar, Kr; with Xe and Rn having very close values.
Nevertheless, there is a trend of rate increase as are move from Ar up to Rn. It is interesting to note that although this behavior is very similar regarding the cis-barrier for all rates (conventional, Wigner e Eckart), it seems that the tunneling is more important to describe the trans barrier's rate, where there is a significant difference for Eckart's values specially in the 100-200 K range when compared to respective Wigner and conventional results.
The final thermal rate constant can be expressed in the two familiar Arrhenius forms. In this work, the first is the Arrhenius modified form given by where A is the pre-exponential factor, T a temperature, n is a real number, R is the universal gas constant and E a is the activation energy. The second is the d-Arrhenius form Silva et al., 2013;Carvalho-Silva et al., 2017) expressed by where d is a parameter that yield the degree of deformation of the exponential function. The curve obtained by the reaction rate constant vs. the temperature can be fitted (Ramalho et al., 2011) to obtain the parameters A, n and E a for the Arrhenius modified form, as presented in Table 4, and the parameters A, d, and E a for d-Arrhenius form, as presented in Table 5. This feature confirms the trend of lower k(T) observed for HP-Ng complexes for both barriers when compared to isolated HP.
It can be also observed in Figure 6 that the trans to cis conformation rate of HP is lower (in the range 100-200 K) than the corresponding ones for HP-Ar, HP-Kr, HP-Xe, and HP-Rn. In the case of the chiral transition from cis to trans, the rates of all HP-Ng complexes are lower than that of the isolated HP. These results showed that the transition rate from cis to trans is greater than the corresponding trans to cis for both the isolated HP molecule and for all HP-Ng complexes. This suggests that the most important barrier that separates the chiral configurations of the isolated HP and the HP-Ng complexes is the trans-barrier, since it is the smallest. The energy of the HP's trans-barrier is relatively small (1.0427 kcal/mol) compared to its cis-barrier (7.5595 kcal/mol) as already seen in Figure 2.
An interesting result is presented in the Table 6. Although the increase in the trans-barrier of the HP-Ng complexes relative to HP is considerably small (see Figure 7), the change in the transition rate from cis to trans is relatively high. This is verified for high (4,000 K), room (298.15 K) and also for low temperatures (100 K). The most pronounced decrease in the rate corresponds to the HP-Ne complex, in which the decrease of the trans-barrier of just 0.0389 kcal/mol (see Figure 6) corresponds to a decrease     of over 60% for the cis-trans transition rate, followed by HP-He. It is also interesting to note that this small change in energy barrier but with a substantial change in rate was also observed for other HP-Ng complexes. For example, HP-Ar complex showed a decrease of just 0.0224 kcal/mol but a 28.68% decrease of rate at 100 K. Finally, at a temperature close to 300K, the Boltzmann distribution shows that about 16% of HP's population has higher energy than the trans-barrier with thermal fluctuations of approximately 1.7686 kcal/mol (Ball and Brindley, 2016). It has also been found that at low temperatures the chiral interconversion quantum encapsulation time of HP is very small. At a temperature of 100 K this time is <1 pico-second (Bitencourt et al., 2008), and at temperatures close to 0 K which can reach 3 pico-seconds.

CONCLUSIONS
The obtained results indicate that the chiral transition rate of trans to cis configuration of hydrogen peroxide in the presence of the noble gases He and Ne were the lowest over the entire temperature range of 4,000-100 K.
The AIM analysis shows that the interaction between H 2 O 2 and the noble gases should be a van der Waals type. Although the H 2 O 2 acts as an acid in the context of this investigation, the high hardness and high electronegativity of the nobles gases hold their electrons very tight to permit a covalence bond between H 2 O 2 and Ng. On the other hand, it seems that both He and Ne are better able to affect the hyperconjugation effect and destabilizing repulsion among the lone pairs that are responsible for rotational barriers (Song et al., 2005). This may explain why the chiral transition rate decreases more for the complexes composed by Ne and He atoms, the hardest and more electronegative noble gases (Furtado et al., 2015).
Finally, the trans-barrier plays an important role because it is much smaller than the cis-barrier. The results showed that a small increase in the trans-barrier height in the complexes is responsible for a significant decrease in the rate of transition from cis to trans. Thus, these effects may contribute to the feasibility of separating one or the other enantiomer of the H 2 O 2 molecule.

AUTHOR CONTRIBUTIONS
RG conceived and supervised the study. RG also helped write the paper. YS performed the H 2 O 2 -Ng electronic and thermal chiral rate calculations. PN determined the H 2 O 2 -Ng minimum and transition state configurations and LdM used the AIM theory to perform the H 2 O 2 -Ng topological analyses and wrote the manuscript, which was reviewed by all authors.