Vibrational energy relaxation (VER) of diatomics following collisions with the surrounding medium is an important elementary process for modeling high-temperature gas flow. VER is characterized by two parameters: the vibrational relaxation time $\tau vib$ and the state relaxation rates. Here the vibrational relaxation of CO$(\nu =0\u2190\nu =1)$ in Ar is considered for validating a computational approach to determine the vibrational relaxation time parameter ($p\tau vib$) using an accurate, fully dimensional potential energy surface. For lower temperatures, comparison with experimental data shows very good agreement whereas at higher temperatures (up to 25 000 K), comparisons with an empirically modified model due to Park confirm its validity for CO in Ar. Additionally, the calculations provide insight into the importance of $\Delta \nu >1$ transitions that are ignored in typical applications of the Landau–Teller framework.

The relaxation of vibrational energy from an energized molecule to its environment is a key process for energy redistribution in gas and condensed phase systems.^{1–9} Specifically, the vibrational relaxation of diatomic molecules during collisions with the surrounding gas is an important elementary process in hypersonic flow typical of spacecraft entering an atmosphere.^{10–16} The key parameters for modeling these systems are the rate coefficients describing the collision of diatomic molecules with other atoms and molecules in the atmosphere. The vibrational relaxation process is essential for characterizing the energy flow from vibrational to rotational and translational degrees of freedom and chemical kinetics can depend on the distribution of energy between these modes.

CO in Ar is a suitable proxy for a generic investigation of vibrational energy relaxation (VER) at high temperatures because experimental data up to 3000 K are available for comparison^{17} from which an extrapolation model due to Millikan and White (MW) based on Landau-Teller (LT) theory was developed.^{10} It is also of interest to note that the atmosphere of Mars contains 3% argon and significant amounts of CO which make Ar–CO relevant for space exploration missions.^{12} For Ar–CO, a recent high-level spectroscopically accurate (CCSD(T)-F12b/aug-cc-pV5Z) potential energy surface (PES) is available.^{18} The *ab initio* grid covers CO diatomic distances *r* from 1.0 Å to 1.35 Å. This interval includes the classical turning points up to $\nu =4$ for CO.^{19} In the present work, this PES was employed to generate a comprehensive grid of energies (3861 points) using a reproducing kernel Hilbert space (RKHS)^{20} representation, following the description of Ref. 21. One advantage of an RKHS is that the necessary derivatives can be computed efficiently and analytically.^{22,23} The general approach for quantifying the vibrational energy relaxation *E*_{vib} is to start from the LT formalism,^{24}

where $\u27e8Evib\u27e9$ is the equilibrium thermal energy of the oscillator in the bath and $\tau vib$ is the vibrational relaxation time. This quantity is commonly reported as the vibrational relaxation time parameter $p\tau vib$, where *p* is the pressure and

In Eq. (2), *k*_{10}(*T*) corresponds to the rate coefficient for the vibrational transition $1\u21920$ at temperature *T*, *k*_{B} is the Boltzmann constant, *ℏ* is the reduced Planck constant, and *ω* is the harmonic frequency of the diatomic molecule. The function $p\tau vib(T)$ initially decreases with temperature because the population of Ar atoms with sufficient energy to efficiently couple to the CO molecule increases with increasing *T*. At some temperature, the efficiency to move vibrational energy in a collision is maximal and $p\tau vib$ increases as it is limited by the time between collisions which turns out to be at $T\u224815\u2009000$ K. However, for Ar–CO the experimental data only extend up to 3000 K. Computing $p\tau vib$ requires the calculation of *k*_{10}(*T*), and at high temperatures quasi-classical trajectory (QCT) calculations are a suitable and meaningful way to do so.

QCT simulations have provided valuable insight into processes such as rotational excitation in the charge exchange reaction $N2+$ + N_{2} → N_{2} + $N2+$^{25} or the state-to-state dynamics of the Br+H_{2} reaction.^{26} Here, QCT calculations are used to determine thermal rate coefficients. The dynamics of the system is followed by propagating the equations of motions for given initial conditions using the Velocity-Verlet algorithm^{27} with a time step of 0.02 fs to ensure the conservation of total energy. Initial conditions for CO were generated from a WKB-quantized periodic orbit^{28} of the corresponding rotating Morse oscillator.^{29} The symmetry axis of CO, the axis of its rotation, and the angular momentum were randomly oriented and the importance sampling Monte Carlo scheme was used.^{30} For the impact parameter *b*, stratified sampling was employed.^{31,32} At a given temperature, *b* was sampled on the interval [0;*b*_{max}] for five different values of *b*_{max}. The $2bdbbmax2$ probability function was employed for sampling *b* in each interval. Finally, all trajectories were analysed jointly by grouping them in non-overlapping intervals of *b* with a weight

where *b*^{−} and *b*^{+} are the minimum and maximum values of *b*, respectively. Then the reaction probability is

where *k* labels one particular stratum (interval), and *N*_{k} and *N*_{k,tot} are the number of reactive trajectories (or event of interest, e.g., $\Delta \nu =\u22121$ for vibrational relaxation) and the total number of trajectories in stratum *k*. The complex is considered to be dissociated when the Ar and CO moieties are separated by more than the initial separation and the final state is determined from subtracting the translational and rotational energies from the total energy and using semiclassical quantization to determine an approximate vibrational quantum number *ν*. Starting from $\nu =\u20091$ relaxation to $\nu =\u20090$ was considered to occur if $\nu <0.5$, and similarly when analyzing relaxation from higher excited states. Using such a criterion should be sufficient here as no final state analysis is sought but only the information whether in a particular trajectory vibrational relaxation has occurred is of interest. This conventional binning procedure has recently been compared^{33} for vibrationally induced photodissociation of HSO_{3}Cl with more refined procedures based on Gaussian binning^{34} and found to yield satisfactory agreement. The rate coefficient from state *i* to state *l* is then determined from

For most of the present work, the initial vibrational quantum number was $\nu =1$, and therefore the vibrational temperature is undefined. Several tests were performed by choosing *b*_{max} at a given *T*. The *N*_{k,tot} at each temperature varies from 9000 to 25 000. QCT calculations between 4000 K and 25 000 K were performed for nine temperatures. At $T<4000$ K, almost no trajectory showed vibrational relaxation.

Figure 1 reports the dependence of $\Delta \nu $ on the range of impact parameter sampled. For a given $\Delta \nu $, the interval of *b* increases with temperature. Hence, the effective cross section for changing the vibrational state of CO through collision with Ar increases with *T*. Transitions with $\Delta \nu =\xb11$ can occur over a wide range of *b*, with 889 $(\Delta \nu =\u22121)$ and 897 $(\Delta \nu =+1)$ events, respectively, at 25 000 K (see Figure 1). Furthermore, the probability of a vibrational transition has no angular preference due to the isotropic long range interaction potential.

Figure 2 shows the rate coefficients for vibrational transitions of CO ($\nu =1$) after collision with Ar between 4000 and 25 000 K. Below 10 000 K only processes with $\Delta \nu =\xb11$ are important and even at 25 000 K, $\Delta \nu =\xb11$ is the most probable process. Over the temperature range of interest, the rates increase by three orders of magnitude for $\Delta \nu =\xb11$. Another process that could become important at high *T* is full atomization (Ar+C+O). However, the number of trajectories with dissociated CO ($DeCO=11.2$ eV) even at 20 000 K is minute ($<0.03%$). At this temperature, very high rotational levels are populated which create a centrifugal barrier and largely prevent CO dissociation.

At high temperatures, excited vibrational states of CO are populated (e.g., at 10 000 K, $p(\nu CO=2)=13$%). Running 10 000 trajectories from the initial state $\nu =2$ at 10 000 K showed that the rate coefficients for $\Delta \nu >1$ are one order of magnitude smaller than $\Delta \nu =\xb11$ (e.g., $k20=1.25\xd7\u200910\u221212$ cm^{3} molecule^{−1} s^{−1} and $k23=1.68\xd710\u221211$ cm^{3} molecule^{−1} s^{−1}). The results for $\nu =2$ show that the calculations satisfy detailed balance: from $k21=1.36\xd710\u221211$ cm^{3}molecule^{−1} s^{−1}, $k12=1.01\xd710\u221211$ cm^{3} molecule^{−1} s^{−1} is obtained which is in excellent agreement with the value for $k12=1.09\xd710\u221211$ cm^{3} molecule^{−1} s^{−1} computed directly from the QCT calculations from initial $\nu =1$. For testing the validity of the LT formalism, 10 000 trajectories for $\nu =6$ were run at 10 000 K and 25 000 K. At the lower temperature, the sum of de-excitation coefficients for $\Delta \nu \u2264\u22122$ was less than 20% of *k*_{65} but as high as 72% at the higher temperature. Hence, $\Delta \nu \xb11$ is an appropriate assumption at 10 000 K but not at 25 000 K.

Table I shows the rate coefficients *k*_{10} and the vibrational relaxation parameter $p\tau vib$. Furthermore, at each *T* the *b*_{max} used in the stratified sampling is reported. The computed $p\tau vib$ from 4000 to 25 000 K are compared with the available data and fits to empirical expressions, see Figure 3. The “MW model” is the expression fitted to the experimental data (measured from 1000 to 3000 K) for Ar–CO. For O_{2}–O_{2}, the MW model faithfully describes available experiments up to $\u22488000$ K. With this in mind and accounting for different molecular parameters for the systems, we expect the MW model to be valid for Ar–CO well beyond the experimentally accessible temperature of 3000 K.

T
. | b_{max}
. | N_{tot}
. | k_{10}
. | $\Delta k10$ . | $p\tau vib$ . | $\u2211\nu \u2032\u22601k1\nu \u2032$ . | $p\tau vibrem$ . |
---|---|---|---|---|---|---|---|

4 000 | 4 | 25 000 | $9.36\xd710\u221214$ | $4.01\xd710\u221214$ | $1.09\xd710\u22125$ | $1.36\xd710\u221213$ | $4.07\xd710\u22126$ |

5 000 | 6 | 25 000 | $2.38\xd710\u221213$ | $7.49\xd710\u221214$ | $6.23\xd710\u22126$ | $6.73\xd710\u221213$ | $1.03\xd710\u22126$ |

6 000 | 6 | 10 000 | $8.02\xd710\u221213$ | $1.85\xd710\u221213$ | $2.55\xd710\u22126$ | $1.96\xd710\u221212$ | $4.22\xd710\u22127$ |

7 000 | 8 | 10 000 | $2.46\xd710\u221212$ | $4.28\xd710\u221213$ | $1.09\xd710\u22126$ | $4.91\xd710\u221212$ | $1.97\xd710\u22127$ |

10 000 | 8 | 9 000 | $1.15\xd710\u221211$ | $1.19\xd710\u221212$ | $4.46\xd710\u22127$ | $2.43\xd710\u221211$ | $5.69\xd710\u22128$ |

13 000 | 6 | 10 000 | $2.54\xd710\u221211$ | $1.72\xd710\u221212$ | $3.31\xd710\u22127$ | $5.75\xd710\u221211$ | $3.12\xd710\u22128$ |

15 000 | 8 | 10 000 | $3.48\xd710\u221211$ | $2.10\xd710\u221212$ | $3.17\xd710\u22127$ | $9.42\xd710\u221211$ | $2.20\xd710\u22128$ |

20 000 | 8 | 9 000 | $8.30\xd710\u221211$ | $4.37\xd710\u221212$ | $2.30\xd710\u22127$ | $2.13\xd710\u221210$ | $1.30\xd710\u22128$ |

25 000 | 8 | 10 000 | $1.07\xd710\u221210$ | $4.61\xd710\u221212$ | $2.74\xd710\u22127$ | $3.10\xd710\u221210$ | $1.11\xd710\u22128$ |

T
. | b_{max}
. | N_{tot}
. | k_{10}
. | $\Delta k10$ . | $p\tau vib$ . | $\u2211\nu \u2032\u22601k1\nu \u2032$ . | $p\tau vibrem$ . |
---|---|---|---|---|---|---|---|

4 000 | 4 | 25 000 | $9.36\xd710\u221214$ | $4.01\xd710\u221214$ | $1.09\xd710\u22125$ | $1.36\xd710\u221213$ | $4.07\xd710\u22126$ |

5 000 | 6 | 25 000 | $2.38\xd710\u221213$ | $7.49\xd710\u221214$ | $6.23\xd710\u22126$ | $6.73\xd710\u221213$ | $1.03\xd710\u22126$ |

6 000 | 6 | 10 000 | $8.02\xd710\u221213$ | $1.85\xd710\u221213$ | $2.55\xd710\u22126$ | $1.96\xd710\u221212$ | $4.22\xd710\u22127$ |

7 000 | 8 | 10 000 | $2.46\xd710\u221212$ | $4.28\xd710\u221213$ | $1.09\xd710\u22126$ | $4.91\xd710\u221212$ | $1.97\xd710\u22127$ |

10 000 | 8 | 9 000 | $1.15\xd710\u221211$ | $1.19\xd710\u221212$ | $4.46\xd710\u22127$ | $2.43\xd710\u221211$ | $5.69\xd710\u22128$ |

13 000 | 6 | 10 000 | $2.54\xd710\u221211$ | $1.72\xd710\u221212$ | $3.31\xd710\u22127$ | $5.75\xd710\u221211$ | $3.12\xd710\u22128$ |

15 000 | 8 | 10 000 | $3.48\xd710\u221211$ | $2.10\xd710\u221212$ | $3.17\xd710\u22127$ | $9.42\xd710\u221211$ | $2.20\xd710\u22128$ |

20 000 | 8 | 9 000 | $8.30\xd710\u221211$ | $4.37\xd710\u221212$ | $2.30\xd710\u22127$ | $2.13\xd710\u221210$ | $1.30\xd710\u22128$ |

25 000 | 8 | 10 000 | $1.07\xd710\u221210$ | $4.61\xd710\u221212$ | $2.74\xd710\u22127$ | $3.10\xd710\u221210$ | $1.11\xd710\u22128$ |

Up to 10 000 K, the agreement between the QCT-computed and the experimental^{17} data is very good. For $T>\u200910\u2009000$ K, the MW model (blue trace, “experiment”) and the QCT simulation results (open red symbols) differ. For high temperatures, the computed data decay more slowly while the MW model . The collision time, $p\tau col=(\pi \mu kBT)1/2/22\pi \sigma p,q2$ where $\sigma p,q=0.5(\sigma p+\sigma q)$ and $\sigma p$ and $\sigma q$ are the hard sphere diameters of the CO and Ar, respectively,^{13} was also included (green line Figure 3). For $T>\u200915\u2009000$ K, the MW model predicts $p\tau vib<p\tau col$.

Avoiding such unphysical behaviour (vibrational equilibrium is reached before collisions can occur) was the motivation to include a correction term that maintained $p\tau vib>p\tau col$ up to an arbitrarily high temperature of 50 000 K.^{12} For Ar–CO, $p\tau vib$ assumes the form^{12}

The $p\tau vib$ from the QCT simulations and Park’s parametrized correction (black line) agree very well at high *T*. In fact, Park’s parametrized model is close to a fit up to 15 000 K, see Figure 3. Finally, the QCT-computed $p\tau vib$ data were also fitted (red solid line) to Park’s model^{12} $p\tau vib=Ae\u2212b/T1/3+cT5/2$ and confirm the *T*^{5/2}-dependence at high temperature. The final fitting parameters are $A=5.53\xd710\u221211$ atm s, *b* = −193.53 K^{1/3}, and $c=2.72\xd710\u221218$ atm s K^{−5/2}.

Conventional application of LT theory assumes that collisions behave in a manner analogous to an optical vibrational transition. Within this approximation (truncation of the system-bath coupling at first order), it is considered that including transitions with $\Delta \nu =1$ is sufficient and those with $\Delta \nu >1$ can be ignored. Since simulations track all transitions, limitations of the $\Delta \nu =1$ assumption can be at least qualitatively tested. For CO the vibrational matrix elements ($|R\nu \nu \u2032|2$)—i.e., the integral of the initial and final vibrational wavefunctions over the infrared dipole operator—are available.^{35} For $\nu =1\u2192\nu \u2032=0$, $|R10|2=1.074\xd710\u22122$ D^{2} and for $\nu =3\u2192\nu \u2032=1$, $|R31|2=1.134\xd710\u22124$ D^{2} which is almost 2 orders of magnitude lower.

In the QCT simulation however, the raw probabilities can not be directly compared across $\Delta \nu $. The probabilities include a bias for the availability of energy to drive the transition. In other words, the probability for $\Delta \nu =2$ may be low at temperatures for which the average collision energy is below two vibrational quanta. A simplistic way to remove the bias is to rescale the probabilities against an effective temperature, $T\u2212T\nu $, where $T\nu $ is the characteristic vibrational temperature, i.e., $T\nu =\nu (CO)\u22c5\Delta \nu /0.8\u22c52/3$. Plotting $p\tau 1,\nu $ against $(T\u2212T\nu )$ and fitting it to $p\tau 1\nu =A(\nu )[T\u2212T\nu ]d$, where $d\u2248\u22121$ provides a means to extract $A(\nu )$ which is the likelihood of a collision leading to a change by $\Delta \nu $. Table II summarizes this data and reports the fitted $A(\nu )$ values. For the QCT simulations, the $A(\nu )$ indicate that $\Delta \nu =2$ is actually an order of magnitude more efficient than $\Delta \nu =1$ when the temperature bias is removed. Hence, for accurate $p\tau $ calculations, vibrational matrix elements are not necessarily useful predictors for the collisional vibrational relaxation efficiency.

$\Delta \nu $ . | $|R1\Delta \nu +1|2$ (D^{2})^{35}
. | $A(\nu )$ (atm s K) . | d
. |
---|---|---|---|

1 | $2.13\xd710\u22122$ | $1.76\xd710\u22122$ | −1.17 |

2 | $1.21\xd710\u22124$ | $2.75\xd710\u22121$ | −1.38 |

3 | $6.65\xd710\u22127$ | $5.78\xd710\u22122$ | −1.16 |

4 | $4.27\xd710\u22123$ | −0.88 | |

5 | $9.80\xd710\u22123$ | −0.96 |

$\Delta \nu $ . | $|R1\Delta \nu +1|2$ (D^{2})^{35}
. | $A(\nu )$ (atm s K) . | d
. |
---|---|---|---|

1 | $2.13\xd710\u22122$ | $1.76\xd710\u22122$ | −1.17 |

2 | $1.21\xd710\u22124$ | $2.75\xd710\u22121$ | −1.38 |

3 | $6.65\xd710\u22127$ | $5.78\xd710\u22122$ | −1.16 |

4 | $4.27\xd710\u22123$ | −0.88 | |

5 | $9.80\xd710\u22123$ | −0.96 |

Recently, Andrienko and Boyd^{14} studied the O_{2}–O vibrational relaxation and characterized the vibrational relaxation time by the total removal of population from $\nu =1$ according to

.These values are also reported in Table I. Eqs. (2) and (7) are equivalent at low temperature ($kBT\u226a\u210f\omega $), where the exponential term in Eq. (2) can be neglected and the sum in Eq. (7) only includes the de-excitation rates. However, at high temperatures all terms are required. Therefore, both $p\tau vib$ and $p\tau vibrem$ are useful for modelling conditions at the hypersonic flight regime.

In summary, extrapolation of the QCT calculations using a high-quality PES correctly describe experimental VER data for Ar–CO for $T\u2265\u20093000$ K. Up to $\u224810\u2009000$ K, the empirical MW analysis follows the simulation data which supports the validity of conventional LT theory $(\Delta \nu =1)$. Above 10 000 K, Park’s correction is required and the QCT-computed data confirm the *T*^{5/2}-dependence. The present simulations are consistent with all the available data, and extending the approach to investigate reactive systems such as N + NO $\u2192$ N_{2} + O is the next step. In these systems,^{36} collisions and subsequent reactions can lead to non-Boltzmann distributions in the product state distributions which will profoundly affect the vibrational relaxation rates and other chemical processes relevant in high-temperature gas flows.

The authors acknowledge stimulating discussions with Professor Tom Schwartzentruber. Part of this work was supported by the United State Department of the Air Force which is gratefully acknowledged (to M.M.). Support by the Swiss National Science Foundation through Grant No. 200021-117810, the NCCR MUST (to M.M.), and the University of Basel is also acknowledged.