VOOZH about

URL: https://arxiv.org/html/2506.06088v1

⇱ High-Energy Reaction Dynamics of O3

High-Energy Reaction Dynamics of O3

JingChun Wang Department of Chemistry, University of Basel, Klingelbergstrasse 80, CH-4056 Basel, Switzerland Juan Carlos San Vicente Veliz Department of Chemistry, University of Basel, Klingelbergstrasse 80, CH-4056 Basel, Switzerland Meenu Upadhyay Department of Chemistry, University of Basel, Klingelbergstrasse 80, CH-4056 Basel, Switzerland Markus Meuwly m.meuwly@unibas.ch Department of Chemistry, University of Basel, Klingelbergstrasse 80, CH-4056 Basel, Switzerland
(June 6, 2025)
Abstract

The high-temperature atom exchange and dissociation reaction dynamics of the O(3P) + O system are investigated based on a new reproducing kernel-based representation of high-level multi-reference configuration interaction energies. Quasi-classical trajectory (QCT) simulations find the experimentally measured negative tempe-rature-dependence of the rate for the exchange reaction and describe the experiments within error bars. Similarly, QCT simulations for a recent potential energy surface (PES) at a comparable level of quantum chemical theory reproduce the negative dependence. Interestingly, both PESs feature a “reef” structure near dissociation which has been implicated to be responsible for a positive dependence of the rate inconsistent with experiments. For the dissociation reaction the dependence correctly captures that known from experiments but underestimates the absolute rates by two orders of magnitude. Accounting for an increased number of accessible electronic states reduces this to one order of magnitude. A neural network-based state-to-distribution model is constructed for both PESs and shows good performance in predicting final translational, vibrational, and rotational product state distributions. Such models are valuable for future and more coarse-grained simulations of reactive hypersonic gas flow.

\externaldocument

si

1 Introduction

Molecular-level characterization of high-energy collisions in shock-heated systems as they occur in atmospheric re-entry or combustion is a challenging undertaking. Important processes under such conditions include non-equilibrium excitations of internal degrees of freedom following inelastic or reactive collisions and full dissociation.1 The gas flow in hypersonics is often in a state of chemical and thermal non-equilibrium because the high translational energy is rapidly converted to internal degrees of freedom. The most prevalent chemical species in reactive air-flow include N2, O2, NO and atomic nitrogen and oxygen. Hence, the chemistry of oxygen-containing species is of particular relevance.1

For high-energy collisions a comprehensive characterization of the state-to-state cross sections is an essential ingredient for microscopic modeling hypersonic, reactive and rarefied gas flow at high temperatures.2, 3, 4, 5 Using high-level electronic structure calculations for full-dimensional, reactive and global potential energy surfaces (PESs) together with their representation and dynamics simulations provides - in principle - the necessary information for their computation. However, characterizing all state-to-state reaction cross sections and rates even for triatomics is a daunting task. This is due to the large number of ro-vibrational reactant states of the diatomic molecule that can combine any diatomic product state (also ), which yields state-to-state cross sections. For converging one such cross section from QCT simulations, typically a minimum number of to simulations are required, which leads to to QCT simulations that would have to be run at a given collision energy. Depending on the range of collision energies ( eV in hypersonics) to be considered the total number of QCT simulations is of the order of to which is neither possible nor desirable.

Over the past few years several machine learning-based treatments to predict reaction outcomes have been reported.6, 7, 8, 9, 10, 11, 12 One of the promising approaches is a state-to-distribution (STD) model. For the N+O2 O+NO reaction an STD model performed rather well (average ) to predict final state distributions for arbitrary initial conditions when compared with explicit QCT simulations.9 More importantly, thermal rates from using the STD model were also in excellent agreement with rates determined from explicit QCT simulations over a wide range of temperatures (1500 K to 20000 K).9

Dynamics studies of reactive processes require global and accurate PESs. A particularly contested feature of the PES for the O(3P) + O reaction concerns a so-called “reef” along the minimum energy path for the atomic oxygen approaching the O2 collision partner at long range.13, 14, 15, 16, 17, 18 Although electronic structure methods at different levels of theory find the reef19, 20, 21, removing it yielded improved agreement between computed and measured thermal rates.13 Also, the absence of a “reef” has been directly linked to a negative temperature dependence of the thermal rate from wavepacket calculations which is consistent with observations.17 Several full-dimensional and reactive PESs have been presented and used in dynamics studies in the past.22, 23, 13, 14, 15, 24, 25, 26, 16, 27, 28 In addition, more local PESs for spectroscopic studies have also been presented.29, 30

For the O(3P) + O collision system multiple reaction pathways are operative: 1) elastic and inelastic collisions, 2) the exchange of oxygen atoms (OA + OBOC OB + OAOC or OC + OAOB) and 3) dissociation of the O2 in the entrance channel (OA + OBOC OA + OB + OC). All these processes are important in atmospheric chemistry and in reactive hypersonic flow.31, 32 In addition, the atom exchange reaction (process 2) is known to exhibit particular isotope effects.33, 34 Both experimental and computational studies were undertaken over the past decades to elucidate the reaction dynamics and underlying intermolecular interactions governing the dynamics.22, 13, 35, 14, 36, 26, 15, 16, 30, 27, 28, 37, 18

The present work reports the thermal rates for atom exchange and O2 decomposition reactions using a reproducing kernel Hilbert space (RKHS)-represented38 PES based on multireference configuration interaction with Davidson correction (MRCI+Q) calculations for validation. In addition, such rates were also determined for an earlier27 permutationally invariant polynomial (PIP)-based PES for comparison and for better characterizing the relationship between the overall shape of the PES and the computed rates. Finally, for both PESs a neural network-based (NN-based) STD model is generated which is useful for more coarse-grained simulations of high-energy reactive airflow.

This work is structured as follows. First the methods are presented. Next, the PESs generated and used are characterized and the results of QCT simulations for the atom exchange and atomization reactions are described. Then, the STD models for the RKHS- and PIP-based PESs are presented and discussed, followed by a discussion of the results and conclusions.

2 Methods

2.1 The Reactive Potential Energy Surface

The ground state PES for O3 (1A′) was constructed at the MRCI+Q level39, 40 with the augmented Dunning-type correlation consistent polarized triple zeta (aug-cc-pVTZ, AVTZ)41 basis set. This level of theory has been found to adequately describe the electronic structure for global and reactive PESs of triatomic C-, N-, and O-containing neutral species.42, 43, 44 Furthermore, such a treatment is consistent with previous work on thermal rates and final state distributions for the [NNO], [OON], and [OOC] collision systems45, 46, 47, 48, 43 which will allow consistent modeling and incorporation of the relevant microscopic information - such as state-to-state or thermal rates - in reaction networks encompassing all these species.3

All electronic structure calculations were carried out in Jacobi coordinates using the Molpro suite of codes and in symmetry.49 The coordinates , , and are the separation between one of the oxygen atoms and the center of mass of O2, the O2 bond length, and the angle between the vectors and , respectively. Reference electronic structure calculations were carried out on a grid comprising a0, a0, and with the other half of the angles defined by symmetry. This leads to a total of 2269 geometrically valid and feasible ground-state geometries.

For the MRCI calculations, multistate CASSCF(12,9) calculations were carried out to correctly characterize the wave function which was then used as the starting point for the MRCI+Q calculations. Preference of MRCI+Q over conventional MRCI calculations was given because higher level of accuracy can be achieved due to the mitigation of size inconsistency of MRCI by introducing the Davidson correction.50, 51 It should be noted, however, that including the Davidson correction in MRCI calculations overestimates the dissociation energy for O3 O + O2 and underestimates the dissociation barrier height13 at the 20 meV (0.5 kcal/mol) level. A total of 8 states, two per spin state (singlet and doublet), and symmetry group were included in the state-averaged calculations.

As has been found previously for such calculations, either the CASSCF or the MRCI+Q calculations may not always converge to the correct electronic states or do not converge at all. Such energies were removed and the grid was reconstructed and completed (“cleaned”) using 1D RKHS and/or 2D RKHS interpolations to evaluate the missing points. Using this cleaned grid the 3-dimensional RKHS representation was generated using the kernel-toolkit.38 The permutationally invariant, reactive PES , where is the vector of three interatomic distances (, and ), was constructed by mixing the three possible channels OAOB+OC, OAOC+OB, and OBOC+OA using an exponential switching function

(1)

for a given structure of O3. Using the “mixing dataset” (see below), the mixing parameters were determined by a grid-based search to yield a0. In total, the RKHS-PES was constructed from 6300 reference energies. This compares with 1686 energies determined at the XMS-CASPT2/maug-cc-pVTZ level of theory from another, more recent global PES for O3.27

The “mixing dataset” was employed for optimizing the mixing parameters in the switching function, see Eq 1. This grid was defined by , a0 and a0 to cover the regions where the three channels overlap.

Finally, an “off-grid” dataset was constructed to validate the overall performance of the mixed RKHS-PES. The geometries were defined by , and and in the range of a0. None of the off-grid points was used either in constructing the single-channel PES or for optimizing the mixing parameters in the crossing regions. All these reference calculations were again carried out at the MRCI+Q/AVTZ level. For the “mixing dataset” energies larger than 300 kcal/mol relative to full dissociation energy of O3 were excluded.

2.2 QCT Simulations

For thermal rates (atom insertion and full dissociation), independent QCT simulations were carried out for a given temperature. Because 16O has nuclear spin only odd initial values for O2 are allowed.17 The maximum simulation time was 75 ps or until the interatomic distance between the initial diatomic O2 was larger than 20 a0 or the distance between the O-atom and either of the atoms in the diatomic O2 exceeded 24 a0. Semiclassical initial conditions were sampled from Boltzmann distributions of the total angular momentum , impact parameter , collision energy , and rovibrational states O. Thermal rates at given temperature were then obtained from

(2)

Here, is the electronic degeneracy factor, is the reduced mass of the collision system, is the Boltzmann constant, and is the number of reactive trajectories (either O-atom exchange or Oatomization). For the exchange reaction,52, 13, 26, 16 whereas for the full dissociation reaction,13, 53, 31 was adopted for all temperatures. The sampling methodology was discussed in detail in Ref. 43. Statistical errors were quantified through bootstrapping. For this, 10 batches of samples with 10 random shuffles of the data (100 times of resamples in total) were used to yield the expected thermal rates along with their standard deviations.

Construction of the STD model for the atom exchange reaction was based on QCT simulations that were carried out on a grid of initial conditions. These included , , and eV. State-specific QCT simulations were carried out for a total of trajectories per initial condition and stratified sampling was used for the impact parameter .

2.3 STD training and evaluation

A neural network model was employed to predict collision outcomes between an atom and a diatom.54 The model was trained using 11 input features that characterized the reactant state, and it produced 178 output nodes corresponding to the amplitudes of the three product states. The NN architecture consisted of seven residual blocks, each comprising two hidden layers. To enhance training efficiency, the input features were standardized (, ) and the outputs were normalized.55 The training process minimized the root-mean-squared deviation (RMSD) between the predicted and reference quasi-classical trajectory (QCT) final state distributions, namely . Since the NN outputs represent probabilities and must remain non-negative even after normalization, a softplus activation function was applied to the output layer.

For training, the weights and biases of the NN were initialized using the Glorot scheme,56 and optimized with the Adam algorithm,57 employing an exponentially decaying learning rate. Training was performed using TensorFlow 1.0,58 and the model parameters that yielded the lowest loss on the validation set were selected for final prediction. Overall, final state distributions from 2414 initial conditions on a grid defined by eV with eV, eV with eV, and eV with eV; with step size 1; and (odd values only) with step size were generated.

The results of the QCT simulations starting from 2414 initial conditions were collected for generating the dataset. Simulations from 34 initial conditions yield insufficiently low reaction probabilities for convergence and were thus excluded from the training set. Initial conditions for these low-probability final states are characterized by low and extremely high values. This is because for low initial values the atom exchange reaction is improbable to occur whereas for the largest translational and/or internal energies, atomization to 3 O(3P) is the dominant final state.

The final dataset (final state distributions from 2380 initial conditions) will be referred to as “on-grid” in the following. A 80:10:10 split of the dataset was randomly drawn for training, validation and testing.54 The “on-grid” dataset is to be distinguished from “off-grid” final state distributions that originate from initial conditions which differ in any of the quantum numbers used for training. For additional technical details, see Ref.54.

3 Results

3.1 Construction and Validation of the PES

First, the quality of the RKHS representation of the MRCI+Q/AVTZ reference calculations is assessed. For the single-channel PES (Figure 1A) the RKHS representation across 9 eV is of exceptional quality ( eV; ). Mixing the three PESs and cyclic permutations increases the RMSD to 0.047 eV with for “on-grid” points, see Figure 1B. The low-energy part of the energy range (between and eV) is still very accurately represented whereas for energies 7 eV above the minimum energy a few outliers appear.

For the mixing dataset, Figure 1C, the RMSD is 0.11 eV with . No obvious outliers occur but the region between –6 and –3 eV is somewhat widened, see inset. Finally, validation of the reactive 3d PES on the offgrid data yields with an RMSD of 0.13 eV which compares with an average RMSD of kcal/mol for a recent PIP-PES across a comparable energy range (up to 200 kcal/mol).27

In constructing the global 3-dimensional reactive PES it was found that certain geometries required particular attention. As a concrete example, the point with is considered. Because for this grid point the electronic structure calculations did not converge, the energy for this point was inferred from converged energies at nearby geometries. The extrapolated MRCI-predicted values is 0.4 kcal/mol lower than the 1D RKHS predicted value, as in Figure S1. The dashed line at stands for the RKHS predictions without the extra data of MRCI extrapolation at this point, while the solid grey line for the RKHS predictions with the extra one point of training data.

👁 Refer to caption
Figure 1: Validation of the RKHS-represented reactive, 3-dimensional PES for O+O2. Panel A: Correlation between single-channel RKHS representation and reference energies for on-grid training data. Panel B: Correlation between 3d mixed RKHS representation and reference energies for the on-grid training data. Panel C: Performance of the 3d mixed RKHS representation on the mixing dataset (see text). Panel D: Performance of the 3d mixed RKHS representation for 582 offgrid points. The RMSD between RKHS representation and reference data and the corresponding are given in each panel.

Next, the global reactive PES was characterized in terms of critical points, shape and comparison with previous PESs which include

  • •

    PES1 based on multi-reference CI with singles and doubles excitation using the AVTZ basis set and represented as three-dimensional cubic splines,59

  • •

    PES2 (referred to as SSB), with energies obtained at the MRCI+Q/cc-pVQZ level of theory and represented as a cubic spline22, 13 based on reference energies,

  • •

    PES3 (referred to as DLLJG) computed at the MRCI-F12/VQZ-F12 level of theory, and using interpolating moving least squares (IMSL)26, 16 to represent energies, and

  • •

    PIP-PES constructed from XMS-CASPT2 energies using the maug-cc-pVTZ basis set and represented as a permutationally invariant polynomial (PIP) surface.27

In what follows, additional details are provided for each of these PESs and the critical points from them are summarized in Table 1. PES1 used electronic structure calculations for the five lowest 1A′ states using the internally contracted multireference configuration interaction method with single and double excitations (MRD-CI) with the AVTZ basis set throughout. The CI wave functions were based on state-averaged (SA) CASSCF orbitals with 18 electrons in 12 orbitals (full-valence active space) and three fully optimized closed-shell inner orbitals. The averaging includes the five lowest 1A′ states with equal weights.

PES2 was constructed at the level of internally contracted MRCI+Q/cc-pVQZ based on CASSCF(12,9) reference wave functions.22, 13 The subsequent MRCI calculations included Davidson correction and PES2 was represented as a three-dimensional splines which reproduce the input energies exactly but no performance measures on off-grid points were provided.

For PES3 the MRCI-F12/VQZ-F12 level of theory was used.26, 16 To improve convergence at large separation, states correlating with the third molecular state of the diatomic oxygen molecule were included in the dynamical weighting SA-CASSCF(18,12) calculations, and 20 states (11 of symmetry 1A′ and 9 of symmetry 1A′′) were included. The MRCI-F12 calculations included 7 1A′ reference states and Davidson correction was applied. The data used for fitting PES3 included energies within 2.6 eV (21000 cm-1) of the global minimum and the reported RMS fitting error was 2 cm-1.16 A cap at kcal/mol with respect to their reference energy was applied for strongly repulsive arrangements.

PES427, referred to as PES-PIP in the remainder of the present work, used multi-state complete active space second-order perturbation theory (XMS-CASPT2) with minimally augmented correlation-consistent polarized valence triple-zeta(maug-cc-pVTZ) basis set based on reference states from SA-CASSCF(12,9) calculations. The 1s and 2s orbitals were fully optimized but were kept doubly occupied in all configurations. To improve convergence to the desired active space, restrictions were adopted on doubly occupied orbitals. In all SA-CASSCF calculations, states were averaged with dynamical weighting. A level shift of 0.3 Eh was applied to mitigate intruder state errors, and the extended multi-state approach was also used with fully invariant treatment of level shifts. The overall accuracy of this PES: RMSE 2.9, 4.5, 8.4, 14.5, 6.1, and 26.2 kcal/mol for energies 100, 100-200, 200-500, 500-1000, 0-1000, and kcal/mol.

Finally, there is also a very accurate representation of ic-MRCI/AV5Z reference data primarily geared towards spectroscopic applications. This analytical fit reproduces the ca. 5000 reference energies to within a few cm-1 over an energy range of 1.24 eV (10000 cm-1).30

Table 1: Minima (MIN1, MIN2) for the RKHS and different PESs from the literature and from experiment. The transition state (TS1) was found using the Nudged Elastic Band (NEB)60, 61 method. MIN1 is the global minimum (open structure) of O3, and MIN2 is the ring minimum eV above the O2()+O(3P) dissociation threshold21. The angle defined as (∘) has atom OA at its apex.
1A′ (a0) (a0) (∘)
MIN1
RKHS-PES 2.43 2.43
PES159 2.42 2.42
PES222, 13 2.41 2.41
PES326, 16 2.40 2.40
PIP-PES27 2.39 2.39
Expt.62 2.40 2.40
MIN2
RKHS-PES 2.74 2.74 29.2
PES2 2.72 2.72 28.7
PES3 2.72 2.72 30.8
PIP-PES 2.69 2.69 27.4
Theo.263 2.73 2.73 29.1
Theo.364 2.72 2.72 30.8
TS
RKHS-PES 2.66 2.66 48.0
PIP-PES 2.63 2.63 55.8
Theo.364 2.76 2.76 57.5
Theo.465 2.66 2.66 55.8

Table 1 summarizes and compares important characteristics of the present PES compared with surfaces PES1 to PIP-PES from the literature.59, 22, 13, 26, 16, 27 The global minimum energy geometry (MIN1) on the RKHS-PES is an isosceles triangle with a0 and an apex angle of at atom OB. There is a second, local minimum (MIN2) with an equilateral triangle structure and separations 2.74 a0 between all oxygen atoms. The energy difference between these two structures on the RKHS PES is 29.2 kcal/mol which compares with 30.7 kcal/mol from minimization of both structures at the MRCI+Q/AVTZ level of theory, and 29.0 kcal/mol on the single-channel RKHS-PES. The difference between the optimizations on the RKHS-PES and the electronic structure calculations are due to slightly different geometries. Optimizations using quantum chemistry compress the bond separations from 2.43 a0 to 2.42 a0 and increase the angle from to . Using the larger AVQZ basis set this energy difference increases by 0.5 kcal/mol.

For MIN1 the optimized structure from the RKHS-PES is consistent with previous work although differences of a0 in bond lengths and up to can be found, see top block in Table 1. The experimentally reported equilibrium structure has a0 and .62 The structure of the higher-lying minimum, MIN2, was also determined in previous work and yields bond lengths within a0 of the present calculations. MIN2 is 29.2 kcal/mol from the RKHS-PES which compares with a range of 27.4 kcal/mol to 30.8 kcal/mol on earlier PESs, calculated at comparable levels of quantum chemical theory but not represented as a RKHS. Finally, the TS between MIN1 and MIN2 features bond lengths of a0 and a bond angle of with an energy of 48.0 kcal/mol above MIN1. This compares with results from the literature that report bond lengths around a0 with an equilibrium angle between and and energies kcal/mol above MIN1 but with considerable variations for the PESs that reported this TS.

👁 Refer to caption
Figure 2: Contours for the O3 RKHS-PES and PIP-PES27. Energies are in kcal/mol with the global minimum (MIN1) energy as the zero of energy. Panels A/B and C/D for the RKHS-PES and PIP-PES, respectively. In panels A and C the two bond lengths forming the angle (OB-OA-OC) are identical. In panels B and D the bond length a0 while scanning the other bond length. Red circles, triangles and crosses designate MIN1 (global minimum), MIN2 (tight minimum), and TS1.

Another property of the PES that impacts the high-energy reaction dynamics of O is the dissociation energy to form O + O. There is no direct measurement of but analysis of thermochemical data yields an estimate of kcal/mol that is typically used for comparison.66, 67, 29, 22, 30, 27 The RKHS-PES yields 21.7 kcal/mol which compares with 19.8 kcal/mol from minimized structures for O3 and O2 at the MRCI+Q/AVTZ level of theory. Using the AVQZ basis set this increases to kcal/mol. This is consistent with early work that reported a pronounced dependence of on the basis set size. Compared with the complete basis set (CBS) limit at the CCSD(T) level, using an AVTZ basis set underestimates by kcal/mol and the CBS limit is lower by 1.3 kcal/mol compared with the predicted .67 Hence, CCSD(T)/AVTZ calculations underestimate the assumed by 4.2 kcal/mol.67 Subsequent icMRCI+Q calculations using the cc-pVQZ basis set reported29 a value of 24.2 kcal/mol, which is larger than the value of 22.6 kcal/mol from CCSD(T)/cc-pVQZ calculations.67 Selected more recent calculations reported values of 23.7 kcal/mol (icMRCI+Q/cc-pVQZ)22, 23.0 to 25.8 kcal/mol (icMRCI+Q with AVTZ to AV6Z basis set)30, and 26.7 kcal/mol ( included in fit; the XMA-CASPT2 value is 36.0 kcal/mol)27. Hence, all these calculations confirm that using the AVTZ basis set underestimates by up to 4 kcal/mol.

The shapes of the present RKHS-PES with those from the literature are compared in Figures 2 and S2. The RKHS-PES and PES-PIP shown in Figures 2A/C and B/D are deceptively similar from a visual comparison. The only major difference concerns the region around and a0 for which the RKHS-PES is tighter and does not extend to longer separations compared with PES-PIP. Given the rather different underlying grid and representation strategy (RKHS vs. PIP) such close agreement is notable. It is anticipated that simulations using these two PESs yield comparable results as will be discussed in the next subsection. Figure S2 compares for the RKHS-PES, PES-PIP and PES3. Again, the three PESs share common topographies in particular in the bound state region. For high-energy regions (yellow isocontours) the RKHS-PES and PES-PIP behave in a comparable fashion whereas PES3 features rather sharp edges.

👁 Refer to caption
Figure 3: Comparison of the dependence of the RKHS-PES (solid lines), the PIP-PES27 (dashed lines), PES220 (blue circles), and PES326 (magenta circles). PES2 features the “reef” whereas PES3 was designed to have the “reef” removed. The Jacobi angles for the RKHS and PIP-PES are (green), (red), and (orange), respectively, see Figure 1A. Depending on the angle considered, the RKHS-PES features a “reef” or not whereas the PIP-PES has a “reef” for all three values of .

Finally, it is of interest to compare 1-dimensional cuts in the near-dissociation region with a particular emphasis on the “reef” structure that has spurred intense discussions.13, 14, 15, 16, 17, 18, 19, 20, 21 Figure 3 compares 1-dimensional scans along the O+O2 dissociation coordinate within kcal/mol of the dissociation energy. Both, the RKHS-PES (solid lines) and PES-PIP (dashed lines) feature more or less pronounced (submerged) barriers depending on the angle of approach . On the other hand, the green solid line do not display this feature. The particular relevance of the “reef” is the fact that a positive temperature dependence for has been implicated to arise from PESs featuring such a reef whereas reef-free PES were found to lead to a negative dependence which is consistent with experiments.13, 14, 15, 16, 17, 18, 19, 20, 21 As will be shown below, from both, RKHS-PES and PES-PIP, lead to a negative dependence in agreement with experiment despite the reefs that are present in Figure 3.

3.2 The Exchange and Full Dissociation Reactions

Next, thermal rates and for the atom exchange (OA + OBOC OB + OAOC or OC + OAOB) and dissociation O(3P) + O 3O(3P) reactions, respectively, are discussed, see Figures 4 and 5.

Exchange Reaction: Thermal rates for the atom exchange reaction (OA + OBOC OB + OAOC or OC + OAOB) from simulations at 100, 200, 300, 400 and 500 K together with associated error bars from bootstrapping using the two reactive PESs, the RKHS PES and PIP-PES, are given in Figure 4. The thermal rates including the electronic degeneracy factor52, 13, 26, 16 from using the RKHS and PIP-PESs are reported as red and blue symbols in Figure 4, respectively. Both PESs lead to a negative dependence for , which is consistent with the experiments (black line with shaded area indicating measurement errors).68 In addition, earlier results from quantum wavepacket simulations using the SSB (green)22 and DLLJG (lila)16 PESs are shown, respectively.

👁 Refer to caption
Figure 4: Comparison of for the O2 + O O2 + O exchange reaction from the present QCT simulations using the RKHS (orange circles) and PIP27 (blue circles) PESs and earlier Jacobi coordinate-based quantum wave packet method simulations using the DLLJG16, 17 (violet) and SSB22, 14 PESs (green). The experimentally measured rates (black line) are shown together with the reported uncertainties (grey background).35, 68

The results indicate that the QCT simulations using the RKHS-PES yield consistent with experiments within the error bars but underestimate the reported values from experiment somewhat. This is also the case for the wavepacket simulations using the DLLJG PES (PES3). QCT simulations using the PIP-PES yield the correct dependence but overestimate somewhat. Finally, using the SSB PES (PES2) leads to a positive dependence which is not what the measurements report. It is of interest to mention that PES3 had the “reef” removed (see Figure 3) whereas the RKHS and PIP-PESs both feature “reefs” as does PES2. Nevertheless, the dependence of the RKHS and PIP-PESs follows that from the measurements which indicates that the presence or absence of the “reef” is not directly related to capturing the correct dependence known from the experiments.

It is of interest to note that the fraction of reactive trajectories ranges from to in the temperature range considered. The majority of the trajectories feature inelastic scattering with changes of () and . At K, 25% of the trajectories are elastic, and 65% of the trajectories are inelastic. The number of trajectories that remain in the O3 state for longer than 75 ps is 2766, 866, 486, 339, 252 out of for 100, 200, 300, 400, 500 K, respectively.

Dissociation Reaction: For the dissociation reaction O + O 3 O the thermal rates are reported in Figure 5. The calculated dissociation rates between 1000 K and 20000 K for both the RKHS-PES and PIP-PES from QCT simulations follow the negative dependence reported from measurements for K.69, 70 Fitting the measured dissociation rates to a linear regression (dashed black line) and shifting to best overlap with the computed rates (dashed red line) demonstrates that the simulations recover the correct dependence (slope) of . It should be noted that for both PESs the number of reactive trajectories decreases significantly as decreases, see Tables S1 and Table S2. The number of trajectories that remain as O3 for longer than 75 ps are 78, 23, 10, 2, 2, 3, 1, 2, 0, 1 out of for , 2000, 3000, 4000, 5000, 6000, 8000, 10000, 15000, and 20000 K, respectively.

👁 Refer to caption
Figure 5: Rates for the dissociation reaction from QCT simulations. Simulation results using the RKHS- and PIP-PESs are the red and blue symbols including error bars from bootstrapping. The measurements are in black69, 70 together with a linear regression (dashed black line). To visually underscore the correct dependence of the simulations, the dashed black line was shifted (dashed red line) to best match the computed data.

Although the dependence for is correctly described compared with the measurements, the dissociation rates are too small by orders of magnitude. There are several factors that potentially influence this finding. First, the dissociation reaction probes, inter alia, the strength of the O2 bond. Geometry optimization of O2 at the MRCI+Q level of theory using the AVTZ and AVQZ basis sets yield 113.5 and 116.4 kcal/mol, respectively, which compare with 115.0 kcal/mol from the RKHS-PES a value of kcal/mol from combining results from measurements.71, 72, 73 On the other hand, when conceiving the PIP-PES the O2 dissociation energy was included in the fit and the value is kcal/mol. Given that the dissociation rates from QCT simulations using the RKHS and PIP-PESs are within a factor of but the value of differs by kcal/mol suggests that is not particularly sensitive to the value of .

Secondly, the electronic degeneracy factor was throughout, as was suggested in previous work13, 53, 2 on high-temperature simulations for O3. Increasing the degeneracy factor - which amounts to including higher lying electronic states that become populated at the collision energies considered in the present work - will increase the dissociation rates and improve the quantitative agreement with experiments. Reconsidering the electronic degeneracy has already been proposed, in particular for high collision energies.31 Based on earlier work, in the high-temperature limit the electronic degeneracy should increase to at least which still includes population of O(1D) states or non-Born-Oppenheimer effects.74 Hence, an increase of the computed rates by one order of magnitude only be adopting a more likely value for is conceivable.

Finally, increasing the size of the basis set used, e.g. AVTZ to AVQZ, will further improve the quality of the PESs. Within transition state theory, a difference of orders of magnitude in the rate points to differences in energies of 2 to 3 kcal/mol which is reminiscent of the increase by 2.9 kcal/mol in the O2 dissociation energy when going from MRCI+Q/AVTZ to MRCI+Q/AVQZ, see above.

3.3 Final State Distribution and NN-based model

For hypersonic modeling, explicit computation of the state-to-state rates from QCT simulations is rather time consuming. More coarse-grained simulations benefit from machine-learned models based on a small subset of the reactant states. In the following, a state-to-distribution (STD) model and its performance for the atom exchange reaction is described. First, representative final state distributions , , , and for the O(3P) + O the reactions OA + OBOC OB + OAOC or OC + OAOB are considered, see Figure 6. Depending on the initial condition (see figure caption) the final state distributions differ appreciably. This is due to the nonequilibrium nature of the conditions under which the reaction occurs. For example, with increasing initial the final shifts to higher energies and broadens. On the other hand, for low but high initial , the final decays to 0 for primarily due to full atomization as a competition channel.

👁 Refer to caption
Figure 6: Product QCT distribution for the 1A′ RKHS-PES from 3 different initial conditions: [ eV] (blue), [ eV] (red) and [ eV] (green). The final , , , and are plotted as a function of the reaction probability. The probability is computed using histogram binning.

Finally, STD models were trained following the procedures described in the methods section for the RKHS and PIP-PESs, see Figures 7 and 8. Both figures report the performance of the trained model on off-grid initial conditions, which were not used for training the NN. The top, middle, and bottom rows report the best, average and worst predictions by the NN. Symbols labelled “grid” refer to the reference amplitudes obtained from running averages over the QCT simulations whereas the solid lines represent the prediction from the NN. For the top two lines the agreement between the NN-trained model and the true QCT simulations is excellent whereas for the worst case (bottom row) is still acceptable but for and the overall shape is captured but details differ.

👁 Refer to caption
Figure 7: NN-STD Model for final state distributions from QCT simulations using the RKHS-PES. Reference amplitudes (Grid) obtained from taking moving averages of the raw QCT data in comparison to STD predictions (NN). The final , , and are plotted as a function of the reaction probability. Three cases for the quality of the NN-trained models are distinguished: initial condition for which the prediction is best (A to C, highest ), is closest to the average (D to F, ), and worst (G to I, ). The corresponding initial conditions are [ eV, , ]; [ eV, , ]; [ eV, , ] for the best, mean and worst cases, respectively. Note that 34 samples were excluded from the dataset based on .

Hence, the present work provides two statistical models for predicting the entire state-to-state distributions based on rigorous QCT simulations. Such NN-based models can be used in more coarse-grained simulations of the reaction kinetics as it is, for example, done using the PLATO (PLAsmas in Thermodynamic nOn-equilibrium) software.75

👁 Refer to caption
Figure 8: NN-STD Model for final state distributions from QCT simulations using the PIP-PES. Reference amplitudes (Grid) obtained from taking moving averages of the raw QCT data in comparison to STD predictions (NN). The final , , and are plotted as a function of the reaction probability. Three cases for the quality of the NN-trained models are distinguished: initial condition for which the prediction is best (A to C, highest ), is closest to the average (D to F, ), and worst (G to I, ). The corresponding initial conditions are [ eV, , ]; [ eV, , ]; [ eV, , ] for the best, average and worst cases, respectively.

4 Discussion and Conclusions

The present work presents, analyzes and uses a new RKHS-based reactive PES for the O(3P) + O collision system suitable for applications to hypersonics. Validation of the PES is done through extensive QCT simulations and comparison with thermal rates for the exchange and atomization reactions at the respective measurement conditions. In both cases the dependence of the rates is correctly captured. Absolute rates are too low by about 2 orders of magnitude for the atomization reaction. Reasons for this include the level of theory at which the electronic structure calculations were carried out and the degeneracy factor which can be expected to be larger than 1/27 as used here and in previous work.31 The chosen level of theory (MRCI+Q/AVTZ) was primarily motivated by the fact that with the present PES a full set of reactive and validated PESs for the [NNO],48 [NOO],47 and [NNN],76 reactive systems is now available which can be used in more coarse-grained studies of combustion processes and reaction networks. Further improvements, such as the use of larger basis sets, are possible but are unlikely to yield qualitatively different results.

QCT simulations using an earlier PIP-PES are consistent with the findings for when running simulations with the RKHS PES. Although the shapes of the PESs are related, the atom exchange rates differ slightly which underlines the sensitivity of such simulations to the topography and features of the PESs. Given that two related PESs (RKHS and PIP) yield somewhat different thermal rates for the exchange reaction motivates the question how to further improve the PESs for yet better agreement with experiment. This can be achieved for example by using morphing approaches guided by experimental information.77, 78 This has been successfully done for the He–H collision system by using measurements of the H translational spectrum.

👁 Refer to caption
Figure 9: The cumulative probability distributions (black isocontours) from 1000 trajectories of each type (elastic, inelastic, atom exchange, atomization - panels A to D) projected onto the 2d-PES for a0. The distributions are shown as normalized heatmaps with minimum and maximum amplitudes between [0.0, 1.0]. White to dark green regions refer to amplitudes of 0.0 and 1.0, respectively.

To better understand how the process in question is related to sampling the underlying PES, Figure 9 reports the PES for a diatomic bond length a0 together with probability distributions from 1000 QCT-trajectories for each of the four relevant processes: elastic (A), inelastic (B), reactive (C), and dissociative (D) atom-plus-diatom collisions. Elastic scattering (Figure 9A) probes primarily the long-range part of the PES and does not sample the region around the global minimum. Inelastic scattering (panel B), on the other hand, samples the region around the global minimum and migration along the two directions and can be clearly identified. In addition, a high-probability zone is also around . Contrary to that, atom exchange reactions penetrate deeper into the region of the two symmetry-related minima but do not sample the T-shaped structure at all, see Figure 9C. Finally, atomization trajectories approach the diatomic mainly along and primarily probe the repulsive wall (Figure 9D), as expected.

Notably, both PESs used in the present work - RKHS and PIP - feature “reefs” in the entrance channel, see Figure 3. Nevertheless, the QCT-simulations using both PESs reproduce the experimentally observed negative dependence of the thermal rates for the exchange reaction. Hence, the notion that the “reef” in earlier PESs is responsible for the positive dependence of , which is inconsistent with experimental observations, is not supported by the present work. For further validation of the RKHS-PES quantum bound state calculations are envisaged for which experimental data is also available for comparison.79

With the present RKHS-PES all relevant atom + diatom reactions for burning air are now described at a uniformly high level of quantum chemical theory (MRCI+Q/AVTZ) with all PESs represented using RKHS. The PESs for all systems [NOO], [NNO], [NNN], and [OOO] were validated vis-a-vis thermal rates for the exchange and atomization reactions and , respectively. This provides a unified framework to investigate the reaction network using more coarse-grained approaches, for example.

As evidenced here, further improvements can be envisaged, for example through transfer learning80 of the PESs to higher levels of theory, such as MRCI+Q/AVQZ. On the other hand, the present work and earlier investigations of the [NNN] reaction system provide evidence that different levels of quantum chemical theory, e.g. CASPT2 instead of MRCI+Q, or larger basis sets do not lead to fundamental changes in the rates, specifically for high-energy collisions. This may be different, though, for reactions at low temperatures for which additional care needs to be exercised.81

In summary, the thermal rates for the atom exchange and atomization reactions using QCT simulations validate the RKHS representation of the PES. Although both surfaces considered in the present work exhibit more or less pronounced “reef” structure in the long-range part of the PES, QCT simulations using them correctly capture the experimentally observed negative dependence for the exchange reaction. This is observed for two different representations of the PES based on two different quantum chemical levels of theory, which further corroborates the validity of the RKHS- and PIP-PESs.

Data Availability

The codes and data for the present study are available from https://github.com/MMunibas/NN-STD-O3 upon publication.

Acknowledgment

The authors gratefully acknowledge financial support from the AFOSR under award number FA8655-21-1-7048, the Swiss National Science Foundation through grants (MM), (MM), the NCCR-MUST (MM), and the University of Basel (MM). The authors thank Prof. R. Dawes for providing source code and help with the DLLJG-PES and Prof. G. Schatz for scientific correspondence.

References

  • Park 1990 Park, C. Nonequilibrium Hypersonic Aerothermodynamics; Wiley-Interscience: New York, 1990.
  • Boyd 2015 Boyd, I. D. Computation of hypersonic flows using the direct simulation Monte Carlo method. J. Space. Rockets 2015, 52, 38–53.
  • Boyd and Schwartzentruber 2017 Boyd, I. D.; Schwartzentruber, T. E. Nonequilibrium Gas Dynamics and Molecular Simulation; Cambridge University Press, 2017; Vol. 42.
  • Koner et al. 2020 Koner, D.; Bemish, R. J.; Meuwly, M. Dynamics on multiple potential energy surfaces: Quantitative studies of elementary processes relevant to hypersonics. J. Phys. Chem. A 2020, 124, 6255–6269.
  • Schouler et al. 2020 Schouler, M.; Prévereaud, Y.; Mieussens, L. Survey of flight and numerical data of hypersonic rarefied flows encountered in earth orbit and atmospheric reentry. Prog. Aerosp. Sci. 2020, 118, 100638.
  • Koner et al. 2019 Koner, D.; Unke, O. T.; Boe, K.; Bemish, R. J.; Meuwly, M. Exhaustive state-to-state cross sections for reactive molecular collisions from importance sampling simulation and a neural network representation. J. Chem. Phys. 2019, 150, 211101.
  • Houston et al. 2019 Houston, P. L.; Nandi, A.; Bowman, J. M. A machine learning approach for prediction of rate constants. J. Phys. Chem. Lett. 2019, 10, 5250–5258.
  • Arnold et al. 2020 Arnold, J.; Koner, D.; Käser, S.; Singh, N.; Bemish, R. J.; Meuwly, M. Machine Learning for Observables: Reactant to Product State Distributions for Atom–Diatom Collisions. J. Phys. Chem. A 2020, 124, 7177–7190.
  • Arnold et al. 2022 Arnold, J.; San Vicente Veliz, J. C.; Koner, D.; Singh, N.; Bemish, R. J.; Meuwly, M. Machine Learning Product State Distributions from Initial Reactant States for a Reactive Atom-Diatom Collision System. J. Chem. Phys. 2022, 156, 034301.
  • Komp et al. 2022 Komp, E.; Janulaitis, N.; Valleau, S. Progress towards machine learning reaction rate constants. Phys. Chem. Chem. Phys. 2022, 24, 2692–2705.
  • Hong et al. 2023 Hong, Q.; Storchi, L.; Bartolomei, M.; Pirani, F.; Sun, Q.; Coletti, C. Inelastic N2 + H2 collisions and quantum-classical rate coefficients: large datasets and machine learning predictions. Eur. Phys. J. D 2023, 77, 128.
  • Huang and Cheng 2024 Huang, X.; Cheng, X. State-to-state dynamics and machine learning predictions of inelastic and reactive O(3P) + CO() collisions relevant to hypersonic flows. J. Chem. Phys. 2024, 160.
  • Fleurat-Lessard et al. 2003 Fleurat-Lessard, P.; Grebenshchikov, S. Y.; Siebert, R.; Schinke, R.; Halberstadt, N. Theoretical investigation of the temperature dependence of the O + O2 exchange reaction. J. Chem. Phys. 2003, 118, 610–621.
  • Babikov et al. 2003 Babikov, D.; Kendrick, B. K.; Walker, R. B.; T Pack, R.; Fleurat-Lesard, P.; Schinke, R. Metastable states of ozone calculated on an accurate potential energy surface. J. Chem. Phys. 2003, 118, 6298–6308.
  • Ayouz and Babikov 2013 Ayouz, M.; Babikov, D. Global permutationally invariant potential energy surface for ozone forming reaction. J. Chem. Phys. 2013, 138.
  • Dawes et al. 2013 Dawes, R.; Lolur, P.; Li, A.; Jiang, B.; Guo, H. Communication: An accurate global potential energy surface for the ground electronic state of ozone. J. Chem. Phys. 2013, 139, 201103.
  • Li et al. 2014 Li, Y.; Sun, Z.; Jiang, B.; Xie, D.; Dawes, R.; Guo, H. Communication: Rigorous quantum dynamics of O + O2 exchange reactions on an ab initio potential energy surface substantiate the negative temperature dependence of rate coefficients. J. Chem. Phys. 2014, 141, 081102.
  • Tyuterev et al. 2014 Tyuterev, V. G.; Kochanov, R.; Campargue, A.; Kassi, S.; Mondelain, D.; Barbe, A.; Starikova, E.; De Backer, M.; Szalay, P.; Tashkun, S. Does the “reef structure” at the ozone transition state towards the dissociation exist? New insight from calculations and ultrasensitive spectroscopy experiments. Phys. Rev. Lett. 2014, 113, 143002.
  • Hernandez-Lamoneda et al. 2002 Hernandez-Lamoneda, R.; Salazar, M. R.; Pack, R. Does ozone have a barrier to dissociation and recombination? Chem. Phys. Lett. 2002, 355, 478–482.
  • Schinke and Fleurat-Lessard 2004 Schinke, R.; Fleurat-Lessard, P. The transition-state region of the O(3P) + O2() potential energy surface. J. Chem. Phys. 2004, 121, 5789–5793.
  • Holka et al. 2010 Holka, F.; Szalay, P. G.; Müller, T.; Tyuterev, V. G. Toward an improved ground state potential energy surface of ozone. J. Phys. Chem. A 2010, 114, 9927–9935.
  • Siebert et al. 2001 Siebert, R.; Schinke, R.; Bittererová, M. Spectroscopy of ozone at the dissociation threshold: Quantum calculations of bound and resonance states on a new global potential energy surface. Phys. Chem. Chem. Phys. 2001, 3, 1795–1798.
  • Siebert et al. 2002 Siebert, R.; Fleurat-Lessard, P.; Schinke, R.; Bittererová, M.; Farantos, S. The vibrational energies of ozone up to the dissociation threshold: Dynamics calculations on an accurate potential energy surface. J. Chem. Phys. 2002, 116, 9749–9767.
  • Ju-Xiang et al. 2007 Ju-Xiang, S.; Zheng-He, Z.; Duo-Hui, H.; Jun, W.; Xin-Lu, C.; Xiang-Dong, Y. Potential energy surfaces of ozone in the ground state. Chin. Phys. 2007, 16, 2650.
  • Mankodi et al. 2017 Mankodi, T. K.; Bhandarkar, U. V.; Puranik, B. P. Dissociation cross sections for N2 + N 3N and O2 + O 3O using the QCT method. J. Chem. Phys. 2017, 146.
  • Dawes et al. 2011 Dawes, R.; Lolur, P.; Ma, J.; Guo, H. Communication: Highly accurate ozone formation potential and implications for kinetics. J. Chem. Phys. 2011, 135, 081102.
  • Varga et al. 2017 Varga, Z.; Paukku, Y.; Truhlar, D. G. Potential energy surfaces for O + O2 collisions. J. Chem. Phys. 2017, 147, 154312.
  • Shu et al. 2024 Shu, Y.; Akher, F. B.; Guo, H.; Truhlar, D. G. Parametrically Managed Activation Functions for Improved Global Potential Energy Surfaces for Six Coupled 5A′ States and Fourteen Coupled 3A′ States of O + O2. J. Phys. Chem. A 2024, 128, 1207–1217.
  • Xie et al. 2000 Xie, D.; Guo, H.; Peterson, K. A. Accurate ab initio near-equilibrium potential energy and dipole moment functions of the ground electronic state of ozone. J. Chem. Phys. 2000, 112, 8378–8386.
  • Tyuterev et al. 2013 Tyuterev, V. G.; Kochanov, R. V.; Tashkun, S. A.; Holka, F.; Szalay, P. G. New analytical model for the ozone electronic ground state potential surface and accurate ab initio vibrational predictions at high energy range. J. Chem. Phys. 2013, 139.
  • Andrienko and Boyd 2016 Andrienko, D. A.; Boyd, I. D. Rovibrational energy transfer and dissociation in O2–O collisions. J. Chem. Phys. 2016, 144.
  • Cacciatore et al. 1978 Cacciatore, M.; Capitelli, M.; Dilonardo, M. Non equilibrium vibrational population and dissociation rates of oxygen in electrical discharges: the role of atoms and of the recombination process. Beiträge aus der Plasmaphysik 1978, 18, 279–299.
  • Van Wyngarden et al. 2007 Van Wyngarden, A. L.; Mar, K. A.; Boering, K. A.; Lin, J. J.; Lee, Y. T.; Lin, S.-Y.; Guo, H.; Lendvay, G. Nonstatistical behavior of reactive scattering in the 18O + 32O2 isotope exchange reaction. J. Am. Chem. Soc. 2007, 129, 2866–2870.
  • Lahankar et al. 2016 Lahankar, S. A.; Zhang, J.; Minton, T. K.; Guo, H.; Lendvay, G. Dynamics of the O-atom exchange reaction 16O(3P) + 18O18O() 16O18O() + 18O(3P) at hyperthermal energies. J. Phys. Chem. A 2016, 120, 5348–5359.
  • Fleurat-Lessard et al. 2003 Fleurat-Lessard, P.; Grebenshchikov, S. Y.; Schinke, R.; Janssen, C.; Krankowsky, D. Isotope dependence of the O+ O2 exchange reaction: Experiment and theory. J. Chem. Phys. 2003, 119, 4700–4712.
  • Fleurat-Lessard et al. 2004 Fleurat-Lessard, P.; Grebenshchikov, S. Y.; Schinke, R.; Janssen, C.; Krankowsky, D. Erratum:“Isotope dependence of the O + O2 exchange reaction: Experiment and theory”[J. Chem. Phys. 119, 4700 (2003)]. J. Chem. Phys. 2004, 120, 4993–4993.
  • Lendvay 2019 Lendvay, G. Mechanism change in the dynamics of the O′ + O2 O′O + O atom exchange reaction at high collision energies. J. Phys. Chem. A 2019, 123, 10230–10239.
  • Unke and Meuwly 2017 Unke, O. T.; Meuwly, M. Toolkit for the Construction of Reproducing Kernel-based Representations of Data: Application to Multidimensional Potential Energy Surfaces. J. Chem. Inf. Model. 2017, 57, 1923–1931.
  • Langhoff and Davidson 1974 Langhoff, S. R.; Davidson, E. R. Configuration interaction calculations on the nitrogen molecule. Int. J. Quant. Chem. 1974, 8, 61–72.
  • Werner and Knowles 1988 Werner, H.; Knowles, P. J. An efficient internally contracted multiconfiguration–reference configuration interaction method. J. Chem. Phys. 1988, 89, 5803–5814.
  • Dunning 1989 Dunning, J., Thom H. Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen. J. Chem. Phys. 1989, 90, 1007–1023.
  • Lin et al. 2016 Lin, W.; Varga, Z.; Song, G.; Paukku, Y.; Truhlar, D. G. Global triplet potential energy surfaces for the N2(X) + O(3P) → NO(X) + N(4S) reaction. J. Chem. Phys. 2016, 144, 024309.
  • Koner et al. 2018 Koner, D.; Bemish, R. J.; Meuwly, M. The C(3P) + NO(X) O(3P) + CN(X), N(2D)/N(4S) + CO(X) reaction: Rates, branching ratios, and final states from 15 K to 20 000 K. J. Chem. Phys. 2018, 149, 094305.
  • San Vicente Veliz et al. 2020 San Vicente Veliz, J. C.; Koner, D.; Schwilk, M.; Bemish, R. J.; Meuwly, M. The N(4S)+ O2(X) O(3P)+ NO(X) reaction: thermal and vibrational relaxation rates for the 2A′, 4A′ and 2A′′ states. Phys. Chem. Chem. Phys. 2020, 22, 3927–3939.
  • Patra et al. 2022 Patra, S.; San Vicente Veliz, J. C.; Koner, D.; Bieske, E. J.; Meuwly, M. Photodissociation dynamics of N. J. Chem. Phys. 2022, 156.
  • Veliz et al. 2021 Veliz, J. C. S. V.; Koner, D.; Schwilk, M.; Bemish, R. J.; Meuwly, M. The C(3P)+ O2() CO2 CO() + O(1D)/O(3P) Reaction: Thermal and Vibrational Relaxation Rates from 15 K to 20000 K. Phys. Chem. Chem. Phys. 2021, 23, 11251–11263.
  • San Vicente Veliz et al. 2020 San Vicente Veliz, J. C.; Koner, D.; Schwilk, M.; Bemish, R. J.; Meuwly, M. The N(4S) + O2(X) O(3P) + NO(X) Reaction: Thermal and Vibrational Relaxation Rates for the 2A′, 4A′ and 2A′′ States. Phys. Chem. Chem. Phys. 2020, 22, 3927–3939.
  • Koner et al. 2020 Koner, D.; San Vicente Veliz, J. C.; ; Bemish, R. J.; Meuwly, M. Accurate Reproducing Kernel-based Potential Energy Surfaces for the Triplet Ground States of N2O and Dynamics for the N + NO O + N2 Reaction. Phys. Chem. Chem. Phys. 2020, 22, 18488–18498.
  • Werner et al. 2020 Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; Celani, P.; Györffy, W.; Kats, D.; Korona, T.; Lindh, R. et al. MOLPRO, version 2020, a package of ab initio programs. 2020.
  • Rintelman et al. 2005 Rintelman, J. M.; Adamovic, I.; Varganov, S.; Gordon, M. S. Multireference second-order perturbation theory: How size consistent is “almost size consistent”. J. Chem. Phys. 2005, 122, 044105.
  • Shiozaki et al. 2011 Shiozaki, T.; Knizia, G.; Werner, H.-J. Explicitly correlated multireference configuration interaction: MRCI-F12. J. Chem. Phys. 2011, 134, 034113.
  • Gross and Billing 1997 Gross, A.; Billing, G. D. Isotope effects on the rate constants for the processes O2 + O O + O2 and O2 + O + Ar O3 + Ar on a modified ground-state potential energy surface for ozone. Chem. Phys. 1997, 217, 1–18.
  • Andrienko and Boyd 2015 Andrienko, D. A.; Boyd, I. D. High fidelity modeling of thermal relaxation and dissociation of oxygen. Phys. Fluids 2015, 27, 116101.
  • Arnold et al. 2022 Arnold, J.; San Vicente Veliz, J. C.; Koner, D.; Singh, N.; Bemish, R. J.; Meuwly, M. Machine learning product state distributions from initial reactant states for a reactive atom-diatom collision system. J. Chem. Phys. 2022, 156, 034301.
  • LeCun et al. 2012 LeCun, Y. A.; Bottou, L.; Orr, G. B.; Müller, K.-R. Neural networks: Tricks of the trade; Springer, 2012; pp 9–48.
  • Glorot and Bengio 2010 Glorot, X.; Bengio, Y. Understanding the difficulty of training deep feedforward neural networks. Proc. 13th Int. Conf. Artif. Intell. Stat. 2010; pp 249–256.
  • Kingma and Ba 2014 Kingma, D.; Ba, J. Adam: A method for stochastic optimization. arXiv preprint arXiv:1412.6980 2014,
  • Abadi et al. 2016 Abadi, M.; Barham, P.; Chen, J.; Chen, Z.; Davis, A.; Dean, J.; Devin, M.; Ghemawat, S.; Irving, G.; Isard, M. et al. Tensorflow: A system for large-scale machine learning. 12th USENIX symposium on operating systems Design and Implementation (OSDI 16). 2016; pp 265–283.
  • Schinke and McBane 2010 Schinke, R.; McBane, G. C. Photodissociation of ozone in the Hartley band: Potential energy surfaces, nonadiabatic couplings, and singlet/triplet branching ratio. J. Chem. Phys. 2010, 132, 044305.
  • Henkelman et al. 2000 Henkelman, G.; Uberuaga, B.; Jonsson, H. A climbing image nudged elastic band method for finding saddle points and minimum energy paths. J. Chem. Phys. 2000, 113, 9901–9904.
  • Kolsbjerg et al. 2016 Kolsbjerg, E. L.; Groves, M. N.; Hammer, B. An automated nudged elastic band method. J. Chem. Phys. 2016, 145, 094107.
  • Tanaka and Morino 1970 Tanaka, T.; Morino, Y. Coriolis interaction and anharmonic potential function of ozone from the microwave spectra in the excited vibrational states. J. Mol. Struct. 1970, 33, 538–551.
  • Gadzhiev et al. 2013 Gadzhiev, O. B.; Ignatov, S. K.; Kulikov, M. Y.; Feigin, A. M.; Razuvaev, A. G.; Sennikov, P. G.; Schrems, O. Structure, Energy, and Vibrational Frequencies of Oxygen Allotropes On () in the Covalently Bound and van der Waals Forms: Ab Initio Study at the CCSD(T) Level. J. Chem. Theor. Comp. 2013, 9, 247–262.
  • Kalemos and Mavridis 2008 Kalemos, A.; Mavridis, A. Electronic structure and bonding of ozone. J. Chem. Phys. 2008, 129, 054312.
  • Chen and Hu 2011 Chen, J.-L.; Hu, W.-P. Theoretical Prediction on the Thermal Stability of Cyclic Ozone and Strong Oxygen Tunneling. J. Am. Chem. Soc. 2011, 133, 16045–16053.
  • Chase et al. 1985 Chase, M. W.; Davies, C. A.; Downey, J. R.; Frurip, D. J.; McDonald, R. A.; Syverud, A. N. JANAF Thermochemical Tables - 3rd edition .1. Al-Co. J. Phys. Chem.  Ref.  Data 1985, 14, 1–926.
  • Müller et al. 1998 Müller, T.; Xantheas, S. S.; Dachsel, H.; Harrison, R. J.; Nieplocha, J.; Shepard, R.; Kedziora, G. S.; Lischka, H. A systematic ab initio investigation of the open and ring structures of ozone. Chem. Phys. Lett. 1998, 293, 72–80.
  • Wiegell et al. 1997 Wiegell, M. R.; Larsen, N. W.; Pedersen, T.; Egsgaard, H. The temperature dependence of the exchange reaction between oxygen atoms and dioxygen molecules studied by means of isotopes and spectroscopy. Int. J. Chem. Kinet. 1997, 29, 745–753.
  • Byron 1959 Byron, S. R. Measurement of the Rate of Dissociation of Oxygen. J. Chem. Phys. 1959, 30, 1380–1392.
  • Shatalov 1973 Shatalov, O. P. Molecular dissociation of oxygen in the absence of vibrational equilibrium. Combust. Explos. Shock Waves 1973, 9, 610–613.
  • Bytautas and Ruedenberg 2010 Bytautas, L.; Ruedenberg, K. Accurate ab initio potential energy curve of O2. I. Nonrelativistic full configuration interaction valence correlation by the correlation energy extrapolation by intrinsic scaling method. J. Chem. Phys. 2010, 132, 074109.
  • Ruscic et al. 2004 Ruscic, B.; Pinzon, R. E.; Morton, M. L.; von Laszevski, G.; Bittner, S. J.; Nijsure, S. G.; Amin, K. A.; Minkoff, M.; Wagner, A. F. Introduction to active thermochemical tables: Several “key” enthalpies of formation revisited. J. Phys. Chem. A 2004, 108, 9979–9997.
  • Creek and Nicholls 1975 Creek, D.; Nicholls, R. A comprehensive re-analysis of the O2 (B) Schumann-Runge band system. Proc. R. Soc. London, Ser. A 1975, 341, 517–536.
  • Nikitin 1974 Nikitin, E. E. Theory of Elementary Atomic and Molecular Processes in Gases; Oxford University Press: Oxford, 1974.
  • Munafò et al. 2020 Munafò, A.; Alberti, A.; Pantano, C.; Freund, J. B.; Panesi, M. A computational model for nanosecond pulse laser-plasma interactions. J. Comput. Phys. 2020, 406, 109190.
  • Wang et al. 2024 Wang, J.; San Vicente Veliz, J. C.; Meuwly, M. High-Energy Reaction Dynamics of N3. J. Phys. Chem. A 2024, 128, 8322–8332.
  • Meuwly and Hutson 1999 Meuwly, M.; Hutson, J. M. Morphing ab initio potentials: A systematic study of Ne–HF. J. Chem. Phys. 1999, 110, 8338–8347.
  • Horn et al. 2024 Horn, K. P.; Vazquez-Salazar, L. I.; Koch, C. P.; Meuwly, M. Improving potential energy surfaces using measured Feshbach resonance states. Sci. Adv. 2024, 10, eadi6462.
  • Chang et al. 1994 Chang, B.-Y.; Kung, C.-Y.; Kittrell, C.; Hsiao, C.-W.; Johnson, B. R.; Glogover, S. G.; Kinsey, J. L. High-accuracy measurement of vibrational Raman bands of ozone at 266 and 270 nm excitations. J. Chem. Phys. 1994, 101, 1914–1922.
  • Käser and Meuwly 2023 Käser, S.; Meuwly, M. Transfer-learned potential energy surfaces: Toward microsecond-scale molecular dynamics simulations in the gas phase at CCSD (T) quality. J. Chem. Phys. 2023, 158, 214301.
  • Hickson et al. 2023 Hickson, K. M.; Veliz, J. C. S. V.; Koner, D.; Meuwly, M. Low-temperature kinetics for the N + NO reaction: experiment guides the way. Phys. Chem. Chem. Phys. 2023, 25, 13854–13863.

SUPPORTING INFORMATION: High-Energy Reaction Dynamics of O3

👁 Refer to caption
Figure S1: of the 3D RKHS refitted with one extra training data compared to the MRCI training data at and a0.
Figure S2: From left to right: , , .

Contour of between a0 of all three PESs at their corresponding angles of tight minimum in valence coordinates. The subfigures from left to right: (a) RKHS-PES (b) PIP-PES;27 (c) PES3.26, 16

(K) Ratio
1000 0 5002227 0
2000 0 5001402 0
3000 0 4999216 0
4000 6 5000171
5000 134 5000606
6000 816 4999019 0.00016
8000 6627 5001339 0.0013
10000 21310 5000949 0.0043
15000 91174 5000299 0.018
20000 177426 4998451 0.035
Table S1: Statistics of the QCT results for dissociation reaction on RKHS PES
(K) Ratio
1000 0 5000800 0
2000 0 5001996 0
3000 0 5000755 0
4000 2 5000152
5000 85 5000915
6000 514 5000549 0.00010
8000 4529 5002607 0.00091
10000 15821 4998037 0.0032
15000 76869 4997764 0.015
20000 156768 5000173 0.031
Table S2: Statistics of the QCT results for dissociation reaction on PIP PES