Pair-correlated speed distributions for the OH+CH<sub>4</sub>/CD<sub>4</sub> reactions: Further remarks on their classical trajectory calculations in a quantum spirit

Laurent Bonnet; José Carlos Corchado; Joaquin Espinosa-Garcia

doi:10.1016/j.crci.2016.03.010

Pair-correlated speed distributions for the OH+CH₄/CD₄ reactions: Further remarks on their classical trajectory calculations in a quantum spirit

Laurent Bonnet ^1,² ; José Carlos Corchado ³; Joaquin Espinosa-Garcia ³

¹ CNRS, Institut des Sciences Moléculaires, UMR 5255, 33405 Talence, France
² Univ. Bordeaux, Institut des Sciences Moléculaires, UMR 5255, 33405 Talence, France
³ Departamento de Química Física, Universidad de Extremadura, Avenida de Elvas S/N, 06071 Badajoz, Spain

Comptes Rendus. Chimie, Volume 19 (2016) no. 5, pp. 571-578.

Abstract

Ten years ago, Liu and co-workers measured pair-correlated speed distributions for OH+CH₄/CD₄ reactions by means of velocity map imaging (VMI) techniques at a collision energy of ∼10 kcal/mol [B. Zhang, W. Shiu, J. J. Lin and K. Liu, J. Chem. Phys 122, 131102 (2005); B. Zhang, W. Shiu and K. Liu, J. Phys. Chem. A 109, 8989 (2005)]. Recently, two of us could semi-quantitatively reproduce these measurements by performing full-dimensional quasi-classical trajectory calculations in a quantum spirit on an ab-initio potential energy surface of their own [J. Espinosa-Garcia and J. C. Corchado, Theor. Chem. Acc. (2015) 134: 6; J. Phys. Chem. B 120, 1446 (2016)]. The goal of the present work is to show that these results can be significantly improved by adding a few more constraints in order to better comply with the restrictions imposed by VMI. Overall, the level of agreement between theory and experiment is remarkable owing to the large dimensionality of the reactions under scrutiny. This is an encouraging result considering the computational challenges of quantum scattering calculations for such large processes.

Metadata

Received: 2015-10-22
Accepted: 2016-03-14
Published online: 2016-04-20

DOI: 10.1016/j.crci.2016.03.010

Keywords: Reaction dynamics, Velocity map imaging, QCT method, Bohr quantization

Author's affiliations:

Laurent Bonnet ^{1, 2}; José Carlos Corchado ³; Joaquin Espinosa-Garcia ³

¹ CNRS, Institut des Sciences Moléculaires, UMR 5255, 33405 Talence, France
² Univ. Bordeaux, Institut des Sciences Moléculaires, UMR 5255, 33405 Talence, France
³ Departamento de Química Física, Universidad de Extremadura, Avenida de Elvas S/N, 06071 Badajoz, Spain

@article{CRCHIM_2016__19_5_571_0,
     author = {Laurent Bonnet and Jos\'e Carlos Corchado and Joaquin Espinosa-Garcia},
     title = {Pair-correlated speed distributions for the {OH+CH\protect\textsubscript{4}/CD\protect\textsubscript{4}} reactions: {Further} remarks on their classical trajectory calculations in a quantum spirit},
     journal = {Comptes Rendus. Chimie},
     pages = {571--578},
     publisher = {Elsevier},
     volume = {19},
     number = {5},
     year = {2016},
     doi = {10.1016/j.crci.2016.03.010},
     language = {en},
}

TY  - JOUR
AU  - Laurent Bonnet
AU  - José Carlos Corchado
AU  - Joaquin Espinosa-Garcia
TI  - Pair-correlated speed distributions for the OH+CH4/CD4 reactions: Further remarks on their classical trajectory calculations in a quantum spirit
JO  - Comptes Rendus. Chimie
PY  - 2016
SP  - 571
EP  - 578
VL  - 19
IS  - 5
PB  - Elsevier
DO  - 10.1016/j.crci.2016.03.010
LA  - en
ID  - CRCHIM_2016__19_5_571_0
ER  -

%0 Journal Article
%A Laurent Bonnet
%A José Carlos Corchado
%A Joaquin Espinosa-Garcia
%T Pair-correlated speed distributions for the OH+CH4/CD4 reactions: Further remarks on their classical trajectory calculations in a quantum spirit
%J Comptes Rendus. Chimie
%D 2016
%P 571-578
%V 19
%N 5
%I Elsevier
%R 10.1016/j.crci.2016.03.010
%G en
%F CRCHIM_2016__19_5_571_0

Laurent Bonnet; José Carlos Corchado; Joaquin Espinosa-Garcia. Pair-correlated speed distributions for the OH+CH₄/CD₄ reactions: Further remarks on their classical trajectory calculations in a quantum spirit. Comptes Rendus. Chimie, Volume 19 (2016) no. 5, pp. 571-578. doi : 10.1016/j.crci.2016.03.010. https://comptes-rendus.academie-sciences.fr/chimie/articles/10.1016/j.crci.2016.03.010/

Version originale du texte intégral

1 Introduction

The fast-developing velocity map imaging (VMI) technique is increasingly used in molecular beam experiments on gas-phase polyatomic reactions. VMI is particularly interesting for reactions leading to two molecules [1–6], for it makes the accurate measurement of the angle-velocity distribution for a given quantum state of one of the two products possible. This pair-correlated angle-velocity distribution allows us to probe the dynamics of polyatomic processes at an amazing level of details, generally possible only for triatomic reactions before the invention of VMI. Consequently, this distribution is emerging as one of the most important experimental data on gas-phase chemical reactivity. In the present work, we mainly focus our developments on the theoretical description of the pair-correlated velocity distribution (PCSD), deduced from the previous density by integration over the angle. It should be noted that at the moment, this distribution appears to have been much more often measured than theoretically predicted [7–11]. There is thus an important gap to be filled. We also briefly approach the more complex issue of the calculation of the pair-correlated angular distribution.

We shall perform our study within the framework of the reactions OH+CH₄→H₂O+CH₃ and OH+CD₄→HOD+CD₃, studied ten years ago by Zhang et al. [12–14] using the previous techniques. These seven-atom reactions leading to two polyatomic products are among the largest bimolecular processes considered to date in this type of experiment. By varying both the collision energy and the probed quantum state of the CH₃/CD₃ radical, they could obtain accurate informations on the subtle way the vibrational state distribution of the co-fragment H₂O/HOD and the angular distribution may evolve.

If a qualitative interpretation of such results is possible from simple physical and chemical arguments, their quantitative analysis requires the development of theoretical models capable of accurately reproducing them. These models can then be used for detailed understanding of the dynamics. They can also represent (inexpensive) alternative approaches to predict data needed by atmospherical chemists and astrophysicists in order to model planetary atmospheres [15] and interstellar clouds [16], respectively.

The previous reactions take place in the electronic ground state. The first step of their accurate theoretical description requires performing ab-initio electronic structure calculations and using fitting procedures to construct the corresponding potential energy surface (PES) [17,18]. This work has recently been performed by two of us [19]. The surface obtained was called PES-2014. Note that the latter takes into account the full dimensionality of the system, i.e. 15 (reducing the dimensionality generally makes far easier the construction of PESs, but this may induce artificial dynamical behaviors ; treating a molecular collision in its full-dimensionality is often preferable). The second step of the simulation consists in moving nuclei in the previous PES. Ideally, one would like to do it quantum mechanically [20,21], and exactly. For the title reactions, however, the basis sizes necessary to converge quantum scattering calculations are huge, thus leading to exceedingly large computation times. The alternative is to classically move nuclei, a method traditionally called the quasi-classical trajectory method (QCTM) [22,23]. However, it is known that in the case where only a few vibrational states are available to the products, as is often the case when VMI is used, velocity distributions should preferably be estimated only from those classical trajectories ending in the products with integer vibrational actions, i.e. satisfying Bohr's quantization principle in the separated products [24]. During the last decade, the QCTM in this old quantum spirit proved to lead to very satisfying predictions as compared to exact quantum scattering calculations or experimental measurements [24]. This approach has been recently used by two of us [19,25] in order to reproduce some of the measurements of Liu and co-workers. In the present work, we focus on those obtained for a collision energy of ∼10 kcal/mol published in Refs. [12] and [14].

The corresponding QCTM results semi-quantitatively reproduced experimental evidence [19,25], which is already quite satisfying given the size of the processes under scrutiny. In particular, the theoretical vibrational state distributions of the co-fragment H₂O/HOD were found to be in good agreement with the experimental ones [25]. On the other hand, the predicted speed distributions associated with given vibrational states of the co-fragment tend to be broader than their measured counterparts and not always well located [25]. As a consequence, the PCSD obtained by summing the previous densities appears to partially or totally loose the structures seen in the experimental PCSD [12,14] (see Section 2.3). Without any further calculations, one might attribute these discrepancies to a possible lack of accuracy of the PES-2014, but we shall see that this is not the case. The goal of the present note is indeed to show that the previous disagreement can be significantly reduced by adding a few more constraints to the analysis of QCTM results so as to better comply with the restrictions imposed by VMI. In Section 2, these changes are detailed within the framework of the OH+CD₄ reaction, and comparison between the predictions of the new analysis and measurements is made for both title processes. We end up this section discussing the high numerical cost of the pair-correlated angular distribution (PCAD), as compared to the PCSD. We conclude in Section 3.

2 Calculations

2.1 Trajectories

We used the VENUS code [26] to run 10⁶ trajectories on the PES-2014 for each process. Full computational details on these calculations can be found in refs. [19,25] (note that in these previous studies, only 10⁵ trajectories were run for each reaction). As far as OH+CD₄→HOD+CD₃ is concerned, 90,646 trajectories turned out to be reactive. The product quantities involved in the next developments are the vibrational energy $E_{C D_{3}}^{v i b}$ of CD₃, its rotational energy $E_{C D_{3}}^{r o t}$ , the vibrational energy $E_{H O D}^{v i b}$ of HOD, its rotational energy $E_{H O D}^{r o t}$ , the vibrational actions a₁, a₂ and a₃ associated with the symmetric stretching, antisymmetric stretching and bending vibrational modes of HOD, respectively, the relative translational energy E_trans between HOD and CD₃, the speeds (with algebraic sign) v of CD₃ and v′ of HOD within the center-of-mass system, and the scattering angle θ. The VENUS code provides directly all these quantities except the actions and the velocities. The actions are obtained by means of the normal-mode analysis (NMA) method [27], which includes anharmonicity and Coriolis-coupling terms. The velocities are deduced from E_trans as follows; the total linear momentum within the center-of-mass system being zero, we have

m_{C D_{3}} v + m_{H O D} v^{'} = 0,

(1)

where

m_{C D_{3}}

is the total mass of CD₃ and m_HOD is the analogous mass for HOD. Moreover,

E_{t r a n s} = \frac{1}{2} m_{C D_{3}} v^{2} + \frac{1}{2} m_{H O D} v^{'^{2}} .

(2)

From the two previous equations, one obtains

v = {[\frac{2 E_{t r a n s}}{m_{C D_{3}} (1 + \frac{m_{C D_{3}}}{m_{H O D}})}]}^{1 / 2} .

(3)

This velocity is the one measured by Zhang et al. using the VMI techniques [12–14].

2.2 Vibrational state populations

The vibrational state of HOD is denoted (n₁, n₂, n₃). These quantum numbers are associated with the same vibrational modes as the previous actions (a₁, a₂, a₃). The corresponding frequencies are ω₁, ω₂, and ω₃ (their values are given in Ref. [25]). The population of state (n₁, n₂, n₃) is denoted as $p_{n_{1} n_{2} n_{3}}$ . In the experiments of Zhang et al. [12–14], frozen reagents meet with a collision energy of ∼10 kcal/mol (among others) and P(v) is measured for CD₃ in its vibrational ground state and its lowest rotational levels (N ≤ 5) [3]. The rotational energy $E_{C D_{3}}^{r o t}$ of the probed CD₃ is thus lower than ∼1 kcal/mol. Three vibrational states are found to be available to HOD (at a collision energy of 10 kcal/mol, to be exact [12]). These states are (1,0,0), (1,0,1) and (2,0,0), by increasing the order of energy.

We shall now estimate their populations from the reactive trajectory results. Since these results are in much larger amount than in Refs. [19,25], the present method will better cope with the constraints imposed by the VMI techniques than before. These constraints will be taken into account as follows:

(i) the zero point energy (ZPE) of CD₃ is $E_{C D_{3}}^{0} = 13.8$ kcal/mol. We assume that all the trajectories violating the ZPE, i.e. such that $E_{C D_{3}}^{v i b} < 13.8$ , do not contribute to the reactivity. We thus discard them from the statistics.
(ii) the energy of the first excited vibrational state of CD₃ is 14.9 kcal/mol [25], corresponding to one quantum in the lowest umbrella mode of CD₃ (396 cm⁻¹ or 1.1 kcal/mol) above the ZPE. Calling δ the previous quantum, the intermediate energy between the two previous energies is $E_{C D_{3}}^{0} + δ / 2 = 14.35 kcal/mol$ . We now assume, in the spirit of standard binning (SB) [24], that all the trajectories leading to $E_{C D_{3}}^{v i b} > 14.35 kcal/mol$ do not contribute to the vibrational ground state of CD₃ probed by Zhang et al. [12–14], and should thus be ignored.
(iii) finally, we also reject all the trajectories leading to $E_{C D_{3}}^{r o t} > 1 kcal/mol$ .

In Fig. 1, the cloud of points of coordinates $(E_{C D_{3}}^{v i b}, E_{C D_{3}}^{r o t})$ is represented for the 90,646 reactive trajectories. Gray/red points correspond to rejected/accepted trajectories. Only 722 paths turn out to be useful, thus representing only ∼0.07% of the whole set of calculated trajectories. This clearly illustrates the fact that simulating polyatomic reactions studied by VMI is numerically demanding, even with the QCTM.

Fig. 1
Cloud of points of coordinates $(E_{C D_{3}}^{v i b}, E_{C D_{3}}^{r o t})$ for the 90646 reactive trajectories relative to the OH+CD₄ reaction. Gray points correspond to the rejected trajectories (89,924) while red points correspond to the accepted ones (722) from which the vibrational state populations and the pair-correlated speed and angular distributions will be constructed.

A purely classical treatment of vibrational motions is realistic when many vibrational states (ideally, several tens) are available to the products. In the present case, however, only three vibrational states of HOD are found to be populated [12] and it is preferable, as stated in the introduction, to take into account Bohr's quantization principle in the analysis of the final trajectory results [24]. For polyatomic molecules, this is usually done by means of the energy-based Gaussian binning procedure, often called the 1GB procedure [24,25,28–31]. The latter puts strong emphasis on those trajectories leading to HOD with its classical vibrational energy close to one of the quantum energies of the available vibrational states. In practice, the 1GB estimation of $p_{n_{1} n_{2} n_{3}}$ is performed as follows: we consider the set of trajectories for which ${\bar{a}}_{1} = n_{1}, {\bar{a}}_{2} = n_{2} and {\bar{a}}_{3} = n_{3}$ , where ${\bar{a}}_{i}$ is the nearest integer of a_i, i = 1, 2 or 3. In other words, trajectories for which (a₁, a₂, a₃) does not belong to the unit cube centered at (n₁, n₂, n₃) are rejected. This selection is common to both the 1GB and SB procedures. However, contrary to SB which assigns the democratic unit weight to each trajectory, the 1GB weight, within the harmonic treatment of HOD vibrational motions, is given by

p_{n_{1} n_{2} n_{3}}^{1 G B} (a_{1}, a_{2}, a_{3}) = e x p (- {[\frac{\sum_{i = 1}^{3} ω_{i} (a_{i} - n_{i})}{ɛ \sum_{i = 1}^{3} ω_{i}}]}^{2})

(4)

[24,25,28–31]. The sum

\sum_{i = 1}^{3} ω_{i} (a_{i} - n_{i})

is proportional to the difference between the harmonic classical and quantum vibrational energies. When anharmonicities are not negligible, the previous difference can be replaced by the one between the exact classical and quantum vibrational energies [30,31]. This is just the case for HOD, and we shall thus use for states (1,0,0), (1,0,1) and (2,0,0) the experimental eigenvalues 2723.66, 4100.05 and 5363.59 cm⁻¹ above the ZPE of HOD, respectively [32]. Note that the set of trajectories used in the anharmonic case is still the one for which (a₁, a₂, a₃) belongs to the unit cube centered at (n₁, n₂, n₃).

p_{n_{1} n_{2} n_{3}}

is obtained by summing the 1GB weights over the previous paths and normalizing to unity the set of available populations. The full width at half maximum of the Gaussian in Eq. (4)

(in terms of the energy difference \sum_{i = 1}^{3} ω_{i} (a_{i} - n_{i}))

is proportional to ɛ. In principle, strictly applying Bohr's quantization would thus require taking ɛ as close to 0 as possible [24,29], so the trajectories contributing to

p_{n_{1} n_{2} n_{3}}

would be those leading to classical vibrational energies nearly equal to the quantum energy of state (n₁, n₂, n₃). In the following, we identify these paths as those for which

p_{n_{1} n_{2} n_{3}}^{1 G B} (a_{1}, a_{2}, a_{3})

is larger than 0.5 (see Eq. (4)) and call them 1GB paths. Their number, among the 722 paths complying with constraints (i–iii), is represented in Fig. 2 as a function of ɛ (blue curve). At ɛ = 0.041, the number of 1GB trajectories is equal to 100, which can be roughly considered as the minimum statistical ensemble for a reasonable estimation of the three 1GB populations. The latter are represented in terms of ɛ in Fig. 3. The vertical red dashed line is defined by ɛ = 0.041. On the right side of this line, we note that the populations vary almost linearly with ɛ. On the left side, the smaller ɛ, the smaller the statistics and the more irregular the variations of the populations. Since the experimental populations are assumed to be the best described by the 1GB populations corresponding to ɛ = 0, and the only satisfying 1GB populations are for ɛ ≥ 0.041, we deduce the former by linearly extrapolating the latter from ɛ = 0.041 down to 0. Approximating the 1GB populations by straight lines passing through the points of the curves corresponding to ɛ = 0.041 and 0.1 (see Fig. 3), we arrive at p₁₀₀ = 23, p₁₀₁ = 12 and p₂₀₀ = 65% (rounded to their nearest integer), to be compared with the experimental ones, i.e. p₁₀₀ = 19, p₁₀₁ = 11 and p₂₀₀ = 70% [12,33]. Note that at ɛ = 0.1, the number of 1GB trajectories is equal to 218.

Fig. 2
Variation in terms of ɛ of the number of 1GB trajectories effectively contributing to the 1GB populations. These paths have been identified as those for which the 1GB statistical weight is larger than 0.5 (see Eq. (4)). The blue and green curves correspond to the OH+CD₄ and OH+CH₄ reactions, respectively. The magenta dotted line shows that the values of ɛ for which 100 paths contribute to the previous processes are slightly larger than 0.04 and 0.05, respectively (their exact values are 0.041 and 0.053).

Fig. 3
1GB vibrational state populations (in percentage) of the product HOD of the OH+CD₄ reaction in terms of ɛ. The vertical red dashed line corresponds to the value of ɛ (0.041) such that 100 trajectories effectively contribute to the populations.

Constraints (i) and (ii) define an energy-based windowing, which is related somehow to SB and 1GB [24,28–31]. The main difference is that the former procedure takes into account paths for which the nearest integers of the vibrational actions of CD₃ are non 0, contrary to the latter. This is probably a bit less physical, but significantly increases the number of trajectories contributing to the populations, in addition to avoiding the determination of the six vibrational actions of CD₃.

In order to improve the statistics, we also made calculations with the maximum value of $E_{C D_{3}}^{r o t}$ kept at 2 kcal/mol (instead of 1 kcal/mol), but this affected the accuracy of the vibrational state populations. Moreover, this artificial shift, whenever used, should be as small as possible, since it allows for more energy channeled into the rotational motion of CD₃, at the expense of its translation motion. This is expected to shift the speed distribution towards small velocities.

2.3 Pair-correlated speed distribution

The contribution of state (n₁, n₂, n₃) to the PCSD P(v) is denoted as $P_{n_{1} n_{2} n_{3}} (v)$ . This quantity will be normalized to unity. From the previous developments, we have

P (v) = p_{100} P_{100} (v) + p_{101} P_{101} (v) + p_{200} P_{200} (v) .

(5)

We now concentrate on the calculation of P₁₀₀(v), P₁₀₁(v) and P₂₀₀(v). We use the same constraints on the rovibrational energy of CD₃ as previously used. Like $p_{n_{1} n_{2} n_{3}}$ , $P_{n_{1} n_{2} n_{3}} (v)$ should better be calculated in a quantum spirit. In order to illustrate this, we first focus on the most probable state (2,0,0) and calculate P₂₀₀(v) in the standard following way ignoring any quantization: we apply the SB procedure, i.e. we focus on trajectories for which (a₁, a₂, a₃) belongs to the unit cube centered at (2,0,0). We then divide in 80 equal intervals the range [0.5, 2.5] comprising the velocities measured in km/s [12] and count the number of trajectories per interval (the resulting curve is made of points separated by 0.025 km/s, thus leading to a satisfying resolution). The resulting SB distribution, normalized to unity, is represented by the blue curve in Fig. 4, to be compared with the magenta experimental curve [12]. The SB curve appears to be much broader than the magenta curve. We thus calculate P₂₀₀(v) by means of the 1GB procedure, following the same route as above, with the difference that trajectories are now assigned 1GB weights with ɛ = 0.041 (see Eq. (4)) instead of unit weights. The resulting distribution is given by the green curve in Fig. 4. The agreement with the experiment appears to be much more satisfying than for the SB procedure. Moreover, we repeat this calculation with ɛ = 0.022 (not shown) and observe that despite more irregularities in P₂₀₀(v) due to worse statistics, its width is nearly unchanged. The distribution at ɛ = 0.041 is thus expected to be close to its extrapolation at ɛ = 0. The analogous distributions for the states (1,0,1) and (1,0,0) are displayed in Figs. 5 and 6. As a matter of fact, the 1GB distributions are slightly shifted towards the large velocities, revealing some inaccuracies in the theoretical treatment. The three distributions, weighted by their corresponding 1GB vibrational state populations, are represented in Fig. 7 by the blue, green and orange curves for the (2,0,0), (1,0,1) and (1,0,0) states, respectively. Summing them up leads to P(v) (see Eq. (5)).

Fig. 4
Contribution of the state (2,0,0) to the pair-correlated speed distribution for the OH+CD₄ reaction. Blue curve: P₂₀₀(v) obtained from QCTM calculations and the SB procedure. Green curve: same as previously obtained with the 1GB procedure. Magenta curve: corresponding experimental distribution of Zhang et al. (see Fig. 2a in Ref. [12]).

Fig. 7
Pair-correlated speed distributions for the OH+CD₄ reaction. Blue, green and orange curves: weighted 1GB distributions p₂₀₀P₂₀₀(v), p₁₀₁P₁₀₁(v) and p₁₀₀P₁₀₀(v), respectively. Red curve: full 1GB distribution obtained from a convolution of the sum of the previous quantities (see Eq. (6)). Magenta curve: experimental distribution of Zhang et al. (see Fig. 2a in Ref. [12]).

At this point, a comment on the experimental measurement of P(v) is in order. In a hypothetical experiment where both reagents would be in the rovibrational ground state and the collision energy would be perfectly controlled, the total energy available to the products would take the same value for the whole set of reactions detected. Since E_trans is equal to the previous energy minus the quantized internal energy of the products, E_trans would also be quantized. Consequently, P(v) would be given by a set of Dirac peaks with given weights. In our calculations, however, we pseudo-quantize the vibration motion through the 1GB procedure, but not the rotation motion. P(v) is thus a rough approximation of what would be observed if the measurement of the final velocities were introducing a blurring of the peaks slightly exceeding the average rotational spacing between them. In real experiments, however, the blurring is much larger, for the collision energy and the reagent rotational states are not fully controlled and VMI itself significantly contributes to the blurring. We thus take this into account in the calculation of the PCSD through the Gaussian convolution

P_{c} (v) = \int d v^{'} \frac{1}{\sqrt{π} η} e x p [- {(\frac{v^{'} - v}{η})}^{2}] P (v^{'}) .

(6)

The value of η is chosen so as to reproduce as satisfyingly as possible the experimental speed distribution at the threshold and cut-off. With η kept at 0.08, one obtains the red curve in Fig. 7, to be compared with the experimental one (magenta curve), extracted from Fig. 2a of Ref. [12]. The convolution is welcome, as a million trajectories is not enough to get a smooth P(v). Without blurring, good agreement between experiment and theory would likely be beyond reach. Compared with the previous predictions, displayed in the lower left panel in Fig. 2 of Ref. [25], the PCSD is clearly improved. In particular, when summing up the three bells in the previous figure, only two bumps are found, instead of the three bumps observed both in the experimental and present theoretical results.

As far as the OH+CH₄ reaction is concerned, five states of H₂ O are found to be available (at a collision energy of 10.2 kcal/mol, to be more precise [14]). These states are (0,0,1), (0,1,0), (1,0,0), (0,1,1) and (1,0,1), by increasing the order of energy. We checked that the harmonic normal mode description of water vibrations is very satisfying for the energies of the previous states. For the PES-2014, the vibrational frequencies of the OH stretching, OH^∗ stretching and HOH bending are ω₁ = 3806, ω₂ = 3762 and ω₃ = 1651 cm⁻¹, respectively. Moreover, the ZPE of CH₃ is 18.8 kcal/mol, and the lowest umbrella frequency is 510 cm⁻¹ or 1.5 kcal/mol. The same type of calculation as before has been performed. Over 10⁶ trajectories run on the PES-2014 [19,25], 123567 are found to be reactive, and 1013 comply with the three constraints imposed by VMI. The dependence of the 1GB populations on ɛ is displayed in Fig. 8. The region of statistical confidence roughly lies on the right side of the red dashed line, defined by ɛ = 0.053. By linearly extrapolating the 1GB populations from ɛ = 0.053 down to 0, we arrive at p₀₀₁ = 6, p₀₁₀ = 53, p₁₀₀ = 19, p₀₁₁ = 10 and p₁₀₁ = 12%. The speed distributions P₀₀₁(v), P₀₁₀(v), P₁₀₀(v), P₀₁₁(v) and P₁₀₁(v) weighted by the previous populations are represented in Fig. 9. With η kept at 0.09, Gaussian convolution of the sum of these densities leads to the PCSD represented by the red curve in Fig. 9, to be compared with the experimental magenta curve extracted from Fig. 2a of Ref. [14]. Note that states (1,0,1) and (0,1,1) have similar energies, so the speed distributions associated with these states overlap (see the blue curves in Fig. 9). This comment holds for states (1,0,0) and (0,1,0) (see the green curves in Fig. 9). This is the reason why the vibrational populations deduced from the speed distributions were previously denoted $p_{01^{*} 1}$ for states (1,0,1) and (0,1,1), $p_{01^{*} 0}$ for states (1,0,0) and (0,1,0) and p₀₀₁ [25]. From the present calculations, $p_{01^{*} 1} = 22$ , $p_{01^{*} 0} = 72$ and p₀₀₁ = 6 (see above), to be compared with the experimental populations, i.e. $p_{01^{*} 1} = 19$ , $p_{01^{*} 0} = 75$ and p₀₀₁ = 6 [14]. Like previously observed, the PCSD is clearly improved as compared with its previous estimation, displayed in the upper left panel in Fig. 2 of Ref. [25]. In particular, when summing up the three bells in the previous figure, a single bell is found, instead of the three bumps observed both in the experimental and present theoretical results.

Fig. 8
1GB vibrational state populations (in percentage) of the product H₂O of the OH+CH₄ reaction in terms of ɛ. The vertical red dashed line corresponds to the value of ɛ (0.053) such that 100 trajectories effectively contribute to the populations.

Fig. 9
Pair-correlated speed distributions for the reaction OH+CH₄. Solid blue, dotted blue, solid green, dotted green and orange curves: weighted 1GB distributions p₁₀₁P₁₀₁(v), p₀₁₁P₀₁₁(v), p₁₀₀P₁₀₀(v), p₀₁₀P₀₁₀(v) and p₀₀₁P₀₀₁(v), respectively. Red curve: final theoretical result obtained from a convolution of the sum of the previous quantities (see Eq. (6)). Magenta curve: experimental distribution of Zhang et al. (see Fig. 2a in Ref. [14]).

Last but not least, we wish to emphasize that standard QCT-SB calculations do not even allow us to reveal the structures seen in the experimental PCSDs, as can be seen in Fig. 10. The QCT-SB distributions are indeed simple bell-shape curves broader than the experimental distributions. In comparison, QCT-1GB distributions are much more realistic (see Figs. 7 and 9).

Fig. 10
Comparison between the experimental pair-correlated speed distributions and those obtained by means of the QCTM within the standard binning procedure for the two reactions of interest.

2.4 High numerical cost of the pair-correlated angular distribution

Following Ref. [19], we estimate the PCAD by limiting the restrictions imposed by VMI to the second constraint in Section 2.2 and using the 1GB procedure with ɛ = 0.041. This distribution is called P_a(θ). For the OH+CD₄ reaction, the number of 1GB trajectories is then equal to 4417. P_a(θ) is found to nearly reproduce the previous prediction represented by the solid black curve in Fig. 7 of Ref. [19]. This distribution turns out to be in qualitative agreement with the experimental PCAD of Zhang et al. represented by the dashed black curve in the same figure (see also Fig. 3a in Ref. [12]). The main difference is that the former involves a forward contribution not seen in the latter. One may thus wonder whether calculating the PCAD in a way consistent with the estimation of the PCSD might make this forward contribution disappear. Moreover, it is logical to use the same approach for estimating different observations measured under the same conditions. This calculation would require adding constraints (i) and (iii) (see Section 2.2) to the method leading to P_a(θ) and deduce the PCAD at ɛ = 0, called P_b(θ), from its values at ɛ = 0.1 and 0.041 by linear extrapolation. However, the PCAD is determined by means of 218 and 100 trajectories at ɛ = 0.1 and 0.041, respectively. Now, these numbers are too small for the PCAD to be reliable, especially in the neighborhood of the poles, as is clearly shown in the study by Truhlar and Blais [34]. In particular, it is shown in Figs. 7 and 8 of the previous work that at least one or two thousand trajectories are necessary to get significant results around the poles. The accurate calculation of P_b(θ) for the title processes would thus imply running at least 10 million trajectories.

3 Conclusion

The development of the velocity map imaging (VMI) technique [1] over the last twenty years [2] has made possible the study of polyatomic reactions at an unprecedent level of details through the measurement of pair-correlated speed distributions. For the title processes, for instance, one may extract from the previous observations accurate data on the subtle way H₂O/HOD vibrates and rotates for a given quantum state (or set of closely spaced quantum states) of CH₃/CD₃ [12–14]. If one is able to accurately reproduce these data by an adequate theoretical treatment of nuclear motions on a high-quality ab-initio potential energy surface, one may then use this treatment to get a deep understanding of the dynamics of the processes under scrutiny and improve the synergy between experiment and theory. We have confirmed in this work that the QCT method can be such an approach provided that Bohr's quantization of product vibration motions is taken into account in the final analysis of trajectory results [24,25,28–30]. This is an encouraging result considering (a) the scarcity of theoretical predictions on pair-correlated speed distributions, and (b) the computational challenges of quantum scattering calculations for processes involving large systems. Moreover, we have confirmed that the PES-2014 correctly describes the dynamics of the title processes. We plan to apply this QCT method in a quantum spirit to several polyatomic bimolecular reactions and photodissociations for which pair-correlated speed distributions have never been theoretically reproduced. It should be noted that pair-correlated angular distributions can also be predicted by this method, though the number of trajectories required should be much larger than for the previous distributions.

Acknowledgments

LB is grateful to Prof. Kopin Liu for clarifying explanations on the cross-beam experiments performed in his laboratory.

References

[1] D.W. Chandler; P.L. Houston J. Chem. Phys., 87 (1987), p. 1445

[2] A.T.J.B. Eppink; D.H. Parker Rev. Sci. Instrum., 68 (1997), p. 3477

[3] J.J. Lin; J. Zhou; W. Shiu; K. Liu Rev. Sci. Instrum., 74 (2003), p. 2495

[4] D. Townsend; W. Li; S.K. Lee; R.L. Gross; A.G. Suits J. Phys. Chem. A, 109 (2005), p. 8661

[5] A.I. Chichinin; S. Kauczok; K.-H. Gericke; C. Maul Int. Rev. Phys. Chem., 28 (2009), p. 607

[6] S.K. Lee; Y.F. Lin; S. Lingenfelter; L. Fan; A.H. Winney; Wen Li J. Chem. Phys., 141 (2014), p. 221101

[7] D. Townsend; A. Lahankar; S.K. Lee; S.D. Chambreau; A.G. Suits; X. Zhang; J. Rheinecker; L.B. Harding; J.M. Bowman Science, 306 (2004), p. 1158

[8] M.L. González-Martínez; L. Bonnet; P. Larrégaray; J.-C. Rayez Phys. Chem. Chem. Phys., 12 (2010), p. 115

[9] M. Monge-Palacios; E. González-Lavado; J. Espinosa-Garcia J. Chem. Phys., 141 (2014), p. 094307

[10] R. Martínez; P.A. Enríquez; M.P. Puyuelo; M. González Chem. Phys., 461 (2015), p. 98

[11] B. Zhang; K. Liu; G. Czakó J. Phys. Chem. A, 119 (2015), p. 7190

[12] B. Zhang; W. Shiu; J.J. Lin; K. Liu J. Chem. Phys., 122 (2005), p. 131102

[13] B. Zhang; W. Shiu; J.J. Lin; K. Liu J. Phys. Chem. A, 109 (2005), p. 8983

[14] B. Zhang; W. Shiu; J.J. Lin; K. Liu J. Phys. Chem. A, 109 (2005), p. 8989

[15] Y.L. Yung; W.B. De More Photo-chemistry of planetary atmospheres, Oxford University Press, Oxford, 1999

[16] P. Pilling; D.P.P. Andrade; R. Neves; A.M. Ferreira-Rodrigues; A.C.F. Santos; H.M. Boechat-Roberty Mon. Not. R. Astron. Soc., 375 (2007), p. 1488

[17] J.M. Bowman; G. Czakó; B. Fu Phys. Chem. Chem. Phys., 13 (2011), p. 8094

[18] J. Espinosa-Garcia, M. Monge-Palacios, J. C. Corchado, Adv. Phys. Chem. Volume 2012, Article ID 164752, 19 pages, Hindawi Publishing Corporation.

[19] J. Espinosa-Garcia; J.C. Corchado Theor. Chem. Acc., 134 (2015), p. 6

[20] G. Nyman; H.-G. Yu Int. Rev. Phys. Chem., 32 (2013), p. 39

[21] H. Song; J. Li; B. Jiang; M. Yang; Y. Lu; Hua Guo J. Chem. Phys., 140 (2014), p. 084307

[22] R.N. Porter; L.M. Raff Dynamics of molecular collisions (W.H. Miller, ed.), Plenum, New York, 1976

[23] T.D. Sewell; D.L. Thomson Int. J. Mod. Phys. B, 11 (1997), p. 1067

[24] L. Bonnet Int. Rev. Phys. Chem., 32 (2013), p. 171

[25] J. Espinosa-Garcia; J.C. Corchado J. Phys. Chem. B, 120 (2016), p. 1446

[26] W.L. Hase; R.J. Duchovic; X. Hu; A. Komornicki; K.F. Lim; D.-H. Lu; G.H. Peslherbe; K.N. Swamy; S.R. Van de Linde; A.J.C. Varandas; H. Wang; R.J. Wolf QCPE Bull., 16 (1996), p. 43

[27] J.C. Corchado; J. Espinosa-Garcia Phys. Chem. Chem. Phys., 11 (2009), p. 10157

[28] G. Czakó; J.L. Bowman J. Chem. Phys., 131 (2009), p. 244302

[29] L. Bonnet; J. Espinosa-Garcia J. Chem. Phys., 133 (2010), p. 164108

[30] G. Czakó J. Phys. Chem. A, 116 (2012), p. 7467

[31] R. Conte; B. Fu; E. Kamarchik; J.M. Bowman J. Chem. Phys., 139 (2013), p. 044104

[32] W.S. Benedict; N. Gailar; E.K. Plyler J. Chem. Phys., 24 (1956), p. 1139

[33]

Note that a linear extrapolation from ɛ = 0.041 down to 0 using the tangent to the 1GB vibrational state populations at 0.041 leads to the QCTM populations p₁₀₀ = 20, p₁₀₁ = 12 and p₂₀₀ = 68%. These are in slightly better agreement with the experiment than those deduced from the previous linear approximation. However, this extrapolation is based on a local approximation of the populations around ɛ = 0.041, while the previous extrapolation is based on a more global approximation between ɛ = 0.041 and 0.1.

[34] D.G. Truhlar; N.C. Blais J. Chem. Phys., 67 (1977), p. 1532

Comments - Policy