Open Access
Issue
A&A
Volume 688, August 2024
Article Number A208
Number of page(s) 7
Section Atomic, molecular, and nuclear data
DOI https://doi.org/10.1051/0004-6361/202450738
Published online 22 August 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1 Introduction

Inelastic collisions between water molecules, H2O + H2O, play an important role in collisional energy transfer in environments such as atmospheres of icy planets (Wirström et al. 2020; Vorburger et al. 2022), in cometary comae (Cochran et al. 2015; Cordiner et al. 2022), and in the atmospheres of water-rich exoplanets (García Muñoz et al. 2024). In nonequilibrium conditions, which often exist in these environments, the collisional energy transfer competes with radiative energy transfer, alters the populations of molecular states, and affects the spectral properties of the emitted radiation (submillimeter, terahertz, and deep infrared). Therefore, an interpretation of the radiation that is emitted or absorbed by molecules in the environments that deviate from local thermodynamic equilibrium (not in local thermodynamic equilibrium, NLTE) requires detailed modeling using codes such as RADEX (Van der Tak et al. 2007), LIME (Brinch & Hogerheijde 2010), or MOLPOP (Ramos & Elitzur 2018), which in turn require as input accurate rate coefficients for collision-induced state-to-state transitions in the molecules of interest. For this purpose, the databases of collisional rate coefficients such as BASECOL (Dubernet et al. 2024) and LAMDA (van der Tak et al. 2020) were developed, and active computational research is ongoing to populate these databases with more rate coefficients to include new molecular collision partners and to broaden the temperature range. The prediction of the rate coefficients requires a quantum mechanical treatment of the internal molecular states, which is implemented in computational chemistry programs such as MOLSCAT (Hutson & Le Sueur 2020, 2019) and HIBRIDON (Alexander et al. 2023), and more recently, MQCT (Mandal et al. 2024).

The H2O + H2O collision process is very difficult to describe with quantum mechanics (Agg & Clary 1991; Boursier et al. 2020; Buffa et al. 2000) because both collision partners must be treated as general asymmetric top rotors, which creates a very large number of individual quantum states of the H2O + H2O system as a whole. Recently, noticeable progress in the modeling of H2O + H2O rotationally inelastic collisions was achieved using a mixed quantum/classical theory approach (Mandal et al. 2023) that is implemented in the code MQCT. This approximate method treats the relative translational motion of collision partners classically (using the mean-field trajectories), while the internal motion of collision partners such as rotations and vibrations are treated with quantum mechanics (time-dependent Schrödinger equation). Based on this approach, first, a database of thermally averaged cross sections (TACS) was reported (Mandal & Babikov 2023a), which was then converted into a database of thermal rate coefficients (Mandal & Babikov 2023b) for 231 transitions in para-H2O and 210 transitions in ortho-H2O due to collisions with a thermalized bath of H2O molecules in the temperature range 5 ≤ T ≤ 1000 K.

It should be stressed, however, that all databases that were developed so far for the H2O + H2O system involve averaging over the rotational states of the projectile H2O assuming a Boltzman distribution of these states at given kinetic temperature, that is, the equilibrium conditions. These data can probably be used for NLTE simulations, but only when the deviation from LTE is relatively small. Therefore, it is still desirable to have a database of state-to-state transitions in H2O + H2O that is specifically built for NLTE simulations.

In principle, a database of rate coefficients might be computed at a certain kinetic temperature Tkin for all possible state-to-state transitions in H2O + H2O system, including all individual states of the target and projectile molecules (i.e., without averaging over the distribution of the internal states of the projectile H2O). Then, in the modeling of radiation transfer, these rate coefficients might be averaged over the states of the projectile H2O using a non-Boltzmann distribution of the internal states (with the appropriate distribution for the NLTE conditions of interest). The practical realization of this approach is nearly impossible for two reasons, however. First, the number of individual transitions in H2O + H2O system is unmanageably large, close to 1.3 million. Second, even if such a database were prepared, it could not be used because the standard radiative transfer codes mentioned above are not equipped to take these state-specific rate coefficients as input for the modeling.

An alternative approach, which is undertaken in this paper, is to assume that although there is no equilibrium between all possible degrees of freedom in the system, the translational and rotational degrees of freedom can be expected to achieve their own equilibria independently, and that they can be approximately characterized by Boltzmann distributions at two different temperatures, Tkin and Trot. Then, Tkin can be used to average over the collision energies of two H2O molecules, while Trot can be used to average over the rotational states of the projectile H2O. This approximation permits us to reduce the dataset to a manageable size, but also it causes the rate coefficients for state-to-state transitions in the target H2O to become dependent on two variables: Tkin and Trot (instead of one common temperature in the case of equilibrium). This method was applied to CO + CO in cometary coma and to H3O+, H2O, and HDO collided with H2 in the interstellar medium (Cordiner et al. 2022; Demes et al. 2023; Faure et al. 2024). A similar two-temperature model is widely used for the description of electrons and ions in nonequilibrium plasma of a gas discharge (Raizer & Allen 1997) because the collisional energy exchange between light and heavy particles is inefficient.

In this paper, we outline a method that permits us to implement this idea in the case of H2O + H2O, and we explore the implications of nonequilibrium conditions for an astrophysical modeling of H2O-rich environments, such as the atmospheres of icy planets and cometary comae.

2 Details of the method

The standard procedure for predicting thermal state-to-state transition rate coefficients kn1n1(T)$\[k_{n_1 \rightarrow n_1^{\prime}}(T)\]$ begins with averaging the state-to-state transition cross sections σn1n2n1n2(E)$\[\sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(E)\]$ over the Boltzman-Maxwell distribution of collision energies E at a certain temperature T, kn1n2n1n2(T)=vave(T)(kBT)2E=0E σn1n2n1n2(E) eEkBT dE.$\[k_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(T)=\frac{v_{\mathrm{ave}}(T)}{\left(k_{\mathrm{B}} T\right)^2} \int_{E=0}^{\infty} E ~\sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(E) ~e^{-\frac{E}{k_{\mathrm{B}} T}} \mathrm{~d} E.\]$(1)

In this expression, vave(T)=8kBT/πμ$\[v_{\mathrm{ave}}(T)=\sqrt{8 k_{\mathrm{B}} T / \pi \mu}\]$ is the average collision velocity, kB is the Boltzmann constant, μ is the reduced mass of the two collision partners, and subscripts n1n2 and n1n2$\[n_1^{\prime} n_2^{\prime}\]$ indicate their initial and final states, respectively. Each n is a composite index that represents a full set of quantum numbers for one molecule. Namely, for rotational states of an asymmetric-top-rotor molecule, such as water, n denotes jkAkC$\[j_{k_A k_C}\]$. If we focus on the first water molecule and treat it as a target, while the second water molecule is treated as a projectile, the next step is the summation over the final states and averaging over the initial states of the projectile, kn1n1(T)=n2wn2(T)n2kn1n2n1n2(T),$\[k_{n_1 \rightarrow n_1^{\prime}}(T)=\sum_{n_2} w_{n_2}(T) \sum_{n_2^{\prime}} k_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(T),\]$(2)

where the thermal populations (weights) of the initial states of the projectile are defined as wn2(T)=(2j2+1)eE2kBTQ2(T).$\[w_{n_2}(T)=\frac{\left(2 j_2+1\right) e^{-\frac{E_2}{k_{\mathrm{B}} T}}}{Q_2(T)}.\]$(3)

Here E2 denotes the energies of the rotational states n2 of the projectile, and Q2(T) is the standard partition function of the projectile, Q2(T)=n2(2j2+1)eE2kBT.$\[Q_2(T)=\sum_{n_2}\left(2 j_2+1\right) e^{-\frac{E_2}{k_{\mathrm{B}} T}}.\]$(4)

Since there in no ortho-to-para conversion (of the projectile states n2) during nonreactive inelastic collisions (Lique et al. 2014), Eqs. (2)(4) are applied separately to the para- and ortho-states of the projectile as if they would be different species, which gives two values of kn1n1(T)$\[k_{n_1 \rightarrow n_1^{\prime}}(T)\]$ for each n1n1$\[n_1 \rightarrow n_1^{\prime}\]$ transition in the target molecule. These need to be added together considering the nuclear spin weights of the ortho- and para-states of the projectile molecule. For example, for H2O water molecules, the ortho-to-para ratio is 3 to 1, and the weights are 1/4 and 3/4 for contributions to kn1n1(T)$\[k_{n_1 \rightarrow n_1^{\prime}}(T)\]$ of the para- and ortho-states of the projectile.

In nonequilibrium conditions, this procedure is modified as follows. It is assumed that the translational (relative) motion of the molecules is thermalized at certain kinetic temperature Tkin, and this temperature enters the analog of Eq. (1), where the averaging over molecule-molecule collision energies is computed, kn1n2n1n2(Tkin )=vave (Tkin )(kBTkin )2E=0E σn1n2n1n2(E) eEkBTkin dE.$\[k_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}\left(T_{\text {kin }}\right)=\frac{v_{\text {ave }}\left(T_{\text {kin }}\right)}{\left(k_{\mathrm{B}} T_{\text {kin }}\right)^2} \int_{E=0}^{\infty} E ~\sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(E) ~e^{-\frac{E}{k_{\mathrm{B}} T_{\text {kin }}}} \mathrm{d} E.\]$(5)

Then, it is assumed that although the rotational motion is also thermalized, it is characterized by a different temperature Trot, and TrotTkin. This rotational temperature describes the distribution of the internal states of the projectile, Q2(Trot )=n2(2j2+1)eE2kBTrot,$\[Q_2\left(T_{\text {rot }}\right)=\sum_{n_2}\left(2 j_2+1\right) e^{-\frac{E_2}{k_B T_{\mathrm{rot}}}},\]$(6) wn2(Trot)=(2j2+1)eE2kBTrotQ2(Trot),$\[w_{n_2}\left(T_{\mathrm{rot}}\right)=\frac{\left(2 j_2+1\right) e^{-\frac{E_2}{\mathrm{k}_{\mathrm{B}} T_{\mathrm{rot}}}}}{Q_2\left(T_{\mathrm{rot}}\right)},\]$(7)

and is used for averaging over this distribution, kn1n1(Tkin ,Trot )=n2wn2(Trot )n2kn1n2n1n2(Tkin ).$\[k_{n_1 \rightarrow n_1^{\prime}}\left(T_{\text {kin }}, T_{\text {rot }}\right)=\sum_{n_2} w_{n_2}\left(T_{\text {rot }}\right) \sum_{n_2^{\prime}} k_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}\left(T_{\text {kin }}\right).\]$(8)

In this way, the value of state-to-state transition rate coefficient depends on two temperatures, Tkin and Trot. At equilibrium, when Trot = Tkin, the result of Eqs. (5)(8) is equivalent to that of Eqs. (1)(4), which can be used to confirm the correctness of the code and the resultant rate coefficients.

The practical implementation of this approach encounters methodological difficulties in the case of a complex molecular system, such as H2O + H2O. The reason is the huge number of individual state-to-state transitions n1n2n1n2$\[n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}\]$ in a model that considers the rotational states of both target and projectile molecules. For example, as previously reported (Mandal & Babikov 2023a), we focused on 231 para-para and 210 ortho-ortho lower-energy transitions in the target H2O, and we used a rotational basis set of 38 para- and 38 ortho-states for the projectile H2O molecule (2888 transitions in the projectile, excluding ortho-to-para conversion), which led to close to 1.3 million transitions overall. For each of these, the energy dependence of the cross section σn1n2n1n2(E)$\[\sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(E)\]$, computed on a grid of points, needs to be smoothly interpolated between these points and extrapolated toward the process threshold at low collision energies and toward the high-energy limit using a smooth function. The result of this fitting must be confirmed in some way (e.g., visually) before this dependence is numerically integrated over the collision energy E. While this is possible for a small number of transitions in simple systems, this procedure becomes unfeasible when the number of transitions is about one million or larger, as is the case here.

Therefore, we developed a different approach. We substituted Eq. (5) into Eq. (8) and exchanged the order of integration over collision energies (taking it outside) with a summation over the states of the projectile (bringing it inside), which gives kn1n1(Tkin ,Trot )=vave (Tkin )(kBTkin )2E=0σn1n1(E,Trot )eEkBTkin EdE.$\[k_{n_1 \rightarrow n_1^{\prime}}\left(T_{\text {kin }}, T_{\text {rot }}\right)=\frac{v_{\text {ave }}\left(T_{\text {kin }}\right)}{\left(k_{\mathrm{B}} T_{\text {kin }}\right)^2} \int_{E=0}^{\infty} \sigma_{n_1 \rightarrow n_1^{\prime}}\left(E, T_{\text {rot }}\right) e^{-\frac{E}{k_{\mathrm{B}} T_{\text {kin }}}} E \mathrm{d} E.\]$(9)

We introduced a thermally averaged cross section (TACS) for the transition n1n1$\[n_1 \rightarrow n_1^{\prime}\]$, defined as σn1n1(E,Trot )=n2wn2(Trot )n2σn1n2n1n2(E).$\[\sigma_{n_1 \rightarrow n_1^{\prime}}\left(E, T_{\text {rot }}\right)=\sum_{n_2} w_{n_2}\left(T_{\text {rot }}\right) \sum_{n_2^{\prime}} \sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(E).\]$(10)

The TACS are computed from the original state-to-state transition cross sections as the sum over the final and average over the initial states of the projectile at chosen values of rotational temperature Trot and for all values of the considered collision energy E, making each TACS a function of these two variables. The number of transitions between the states of a target molecule is relatively small. For example, in our case, it was 231 transitions in para-water and 210 transitions in ortho-water (for the quenching and the corresponding excitation processes), which were computed at six values of the collision energy. Therefore, it is quite manageable to check the behavior of all TACS before they are integrated in Eq. (9).

Equation (9) can be used as is to process the results of exact methods, such as full quantum scattering theory (Demes et al. 2023). We built upon the method developed by Mandal & Babikov (2023b) specifically for the mixed quantum/classical calculations, which permitted us to enforce microscopic reversibility using the cross-section calculations in the excitation and quenching directions. In this case, Eq. (9) is replaced by a slightly more complicated formula (Mandal & Babikov 2023b), kn1n1(Tkin ,Trot )=vave (Tkin )(kBTkin )2(eΔE2kBTkin 2j+1)×Umin σ~n1n1(U,Trot )eUkBTkin [1+(ΔE4U)2][1(ΔE4U)2]UdU.$\[\begin{aligned}& k_{n_1 \rightarrow n_1^{\prime}}\left(T_{\text {kin }}, T_{\text {rot }}\right)=\frac{v_{\text {ave }}\left(T_{\text {kin }}\right)}{\left(k_{\mathrm{B}} T_{\text {kin }}\right)^2}\left(\frac{e^{-\frac{\Delta E}{2k_{\mathrm{B}} T_{\text {kin }}}}}{2 j+1}\right) \\& \times \int_{U_{\text {min }}}^{\infty} \tilde{\sigma}_{n_1 n_1^{\prime}}\left(U, T_{\text {rot }}\right) e^{-\frac{U}{k_{\mathrm{B}} T_{\text {kin }}}\left[1+\left(\frac{\Delta E}{4 U}\right)^2\right]}\left[1-\left(\frac{\Delta E}{4 U}\right)^2\right] U \mathrm{d} U.\end{aligned}\]$(11)

Here, U is the kinetic energy of the mixed quantum/classical trajectories, and Umin = ΔE/4 is its threshold value for a given n1n1$\[n_1 \rightarrow n_1^{\prime}\]$ transition. ΔE=En1En1$\[\Delta E=E_{n_1^{\prime}}-E_{n_1}\]$ is the energy change in the transition, which is negative for quenching and positive for excitation processes. Finally, the TACS in Eq. (11) represents a weighted average of the TACS for the quenching and excitation directions of the same transition (Mandal & Babikov 2023b), σ~n1n1(U,Trot)=12[(2j+1)σn1n1(U,Trot)+(2j+1)σn1n1(U,Trot)].$\[\begin{aligned}\tilde{\sigma}_{n_1 n_1^{\prime}}\left(U, T_{\mathrm{rot}}\right)= & \frac{1}{2}\left[(2 j+1) \sigma_{n_1 \rightarrow n_1^{\prime}}\left(U, T_{\mathrm{rot}}\right)\right. \\& \left.+\left(2 j^{\prime}+1\right) \sigma_{n_1^{\prime} \rightarrow n_1}\left(U, T_{\mathrm{rot}}\right)\right].\end{aligned}\]$(12)

For a detailed derivation of Eqs. (11)(12), we refer to Mandal & Babikov (2023b). Each TACS in Eq. (12) was computed as the sum over the final and average over the initial states of the projectile, as before, σn1n1(U,Trot )=n2wn2(Trot )n2σn1n2n1n2(U).$\[\sigma_{n_1 \rightarrow n_1^{\prime}}\left(U, T_{\text {rot }}\right)=\sum_{n_2} w_{n_2}\left(T_{\text {rot }}\right) \sum_{n_2^{\prime}} \sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}(U).\]$(13)

thumbnail Fig. 1

Function F, which includes all factors under the integral of Eq. (11) for several transitions between the ground state 000 of para-H2O and its several low-lying rotational states (as indicated by different colors in the figure) at Tkin = 100 K. The solid, dashed, and dash-dotted lines correspond to a rotational temperature of Trot =10, 50, and 100 K, respectively.

3 Results and discussion

From previous work (Mandal & Babikov 2023a), the cross-sections σn1n2n1n2$\[\sigma_{n_1 n_2 \rightarrow n_1^{\prime} n_2^{\prime}}\]$ are available at six collision energies: U = 133, 200, 267, 400, 533, and 708 cm−1. These were used to compute six values of σn1n1(U,Trot )$\[\sigma_{n_1 \rightarrow n_1^{\prime}}\left(U, T_{\text {rot }}\right)\]$ at any given input temperature Trot using Eq. (10). The excitation and quenching cross sections were combined in Eq. (12) to give σ~n1n1(U,Trot )$\[\tilde{\sigma}_{n_1 n_1^{\prime}}\left(U, T_{\text {rot }}\right)\]$. These were interpolated and integrated over U for any given input temperature Tkin according to Eq. (11) using a numerical procedure (Mandal & Babikov 2023b) with minor modifications. Figure 1 gives examples of the function under the integral of Eq. (11) for several transitions (which, when divided by kBTkin, represents a scaled cross section and has units of Å2). This figure shows the data points computed using the MQCT method in Mandal & Babikov (2023a) and a smooth fitting function from Mandal & Babikov (2023b) that interpolates these data, expands toward the process threshold at low energies, and extrapolates toward the high-energy limit. The resultant values of kn1n1(Tkin ,Trot )$\[k_{n_1 \rightarrow n_1^{\prime}}\left(T_{\text {kin }}, T_{\text {rot }}\right)\]$ are analyzed below and were used to examine the implications for an astrophysical modeling. Tkin and Trot range from 5 to 1000 K.

First, in order to confirm that the new theory and code are correct, we computed the thermal rate coefficients kn1n1(T)$\[k_{n_1 \rightarrow n_1^{\prime}}(T)\]$ at equilibrium, when Tkin = Trot, and compared them with the data published previously. In Fig. 2 this comparison is presented for para-H2O at three temperatures. In Fig. 3, this information is presented for ortho-H2O. The agreement at lower temperatures is reasonably good, but worsens slightly at higher temperatures. The differences for some rate coefficients reach a factor of two to three. This is expected, because in the previously published data, a more approximate method was used to compute the TACS. Namely, in the Eq. (13) above, the TACS is a function of two variables, σn1n1(U,Trot)$\[\sigma_{n_1 \rightarrow n_1^{\prime}}\left(U, T_{\mathrm{rot}}\right)\]$, and the temperature Trot is an input parameter that is independent of U. Importantly, this temperature defines the values of the state populations wn2(Trot)$\[w_{n_2}\left(T_{\text {rot}}\right)\]$ in Eq. (13), according to Eq. (7). In contrast, in the previous work, we assumed that the average collision energy U(T) can approximately represent all collisions at a given temperature. This in turn defined the values of the state populations wn2(U)$\[w_{n_2}(U)\]$ and effectively meant that the TACS were a function of one variable only, σn1n1(U)$\[\sigma_{n_1 \rightarrow n_1^{\prime}}(U)\]$, without an explicit dependence on temperature. As a result of these assumptions, in the previous work, the weights wn2(U)$\[w_{n_2}(U)\]$ were different for different energy points U, while here, the weights wn2(Trot)$\[w_{n_2}\left(T_{\text {rot}}\right)\]$ are the same for different energy points because they are detrmined by the chosen rotational temperature, not by the collision energy. This explains some of the differences between the results of this work and previous work shown in Figs. 2 and 3, and it suggests that the approximation employed in previous works is fine at lower temperatures when the distribution of the collision energies is narrow (e.g., at 200 K), but becomes less accurate at higher temperatures when the distribution broadens (e.g., at 800 K). The method implemented here is more accurate, and we recommend following it.

We would also like to note that in our previous work (Mandal & Babikov 2023a,b), we assumed equilibrium between all rotational states of the H2O projectile, including the para- and ortho-states. This was done largely for historical reasons, following the earlier work on H2O + H2O by other groups (Buffa et al. 2000; Boursier et al. 2020), which did not distinguish between para-and ortho-states of the H2O projectile. Here, we employed a rigorous approach by considering para- and ortho-H2O projectiles separately, averaging over their states independently, and then computing the weighted average of them using a standard ortho-to-para ratio of 3 (as described in Sec. 2 above). Moreover, we created a user-ready code that allows computing a weighted average of ortho and para projectiles using any other value of the ortho-to-para ratio that is specified as an input parameter. For example, if in a certain model for a certain astrochemical application, we wish to eliminate ortho-H2O completely, this can easily be done using the input value of zero for the ortho-to-para ratio. This code is available for download from the GitHub1.

In Fig. 4, we present for several of the most intense transitions in H2O + H2O the dependence of the rate coefficient on both Tkin and Trot as a 2D function. In each frame of this figure, the points along the diagonal line correspond to the equilibrium condition, Tkin = Trot, while the off-diagonal points describe nonequilibrium situations: Tkin < Trot is the area below the diagonal line, and Tkin > Trot is the area above the diagonal line. The farther away from the diagonal line, the larger the deviation from equilibrium. The temperatures range from 5 to 1000 K. In many cases, the highest and lowest rate coefficients are observed along the horizontal and vertical axis of the frames in Fig. 4, where one of the two temperatures is relatively low (on the order of 10 K), while the other temperature is far higher (on the order of 100 K and above). We found that in extreme nonequilibrium conditions like this, the rate coefficients for many transitions can deviate significantly from their equilibrium values.

thumbnail Fig. 2

Comparison of the rate coefficients computed here for the case of thermodynamic equilibrium Tkin = Trot (vertical) vs. thermal rate coefficients reported previously (Mandal & Babikov 2023b; horizontal) for the quenching of 231 states of para-H2O at three temperatures.

thumbnail Fig. 3

Same as Fig. 2, but for the quenching of 210 states of ortho-H2O.

thumbnail Fig. 4

2D dependence of the nonequilibrium rate coefficients for 16 strong transitions in the H2O + H2O system, as indicated in the figure. The rotational and kinetic temperatures are plotted along the horizontal and vertical axis, respectively. Color is used to show the magnitude of rate coefficients in units of cm3 s−1 in log scale as indicated in the color bar.

4 Nonequilibrium radiative transfer modeling

In order to determine the impact of the rotational temperature of H2O projectiles on the excitation conditions of H2O targets in nonequilibrium conditions, we performed radiative transfer modeling using the new collisional data.

The coupled radiative transfer and statistical equilibrium equations were solved using the escape probability approximation at steady state, implemented in the public version of the RADEX code (Van der Tak et al. 2007). We selected physical conditions corresponding to a typical comet located at Rh = 1 AU from the Sun with a total production rate from the nucleus QH2O=1029 s1$\[Q_{\mathrm{H}_2 \mathrm{O}}=10^{29} \mathrm{~s}^{-1}\]$, an expansion velocity vexp = 0.8 km s−1, and at two different gas kinetic temperatures, Tkin = 50 and 100 K. Water was the only projectile considered, and the H2O density was described by a standard isotropic Haser model (Haser 1957). The H2O column density was set to 1017 cm−2, which corresponds to the typical column density in the part of the coma where NLTE effects are expected to be important. At this column density, the opacities of the H2O lines will be high. We then anticipate that the transition from LTE to NLTE will occur at a lower H2O density than in the optically thin regime. For the background radiation, we incorporated in the calculations the cosmic microwave radiation (2.73 K) along with solar radiation. Regarding the latter, we adopted the same approach as Faure et al. (2020) and used a diluted blackbody radiation of 5770 K with a dilution factor W = Ωs/4π, where Ωs equals 6.79 × 10−5 sr, representing the solid angle of the Sun as seen from 1 AU.

The results of our NLTE models are presented in Figs. 5 and 6, where the populations of several lower-energy rotational levels of H2O (both para- and ortho-H2O) are plotted as a function of the H2O density for Tkin = 50 and 100 K. We selected these levels since they are the most populated levels in typical physical coma conditions, but also because they are all radiatively connected to one or more rotational states and then would be involved in radiative transitions detected by telescopes. The distribution of the H2O rotational levels of the projectiles was performed assuming a rotational temperature of Trot = 10, 30, and 50 K at Tkin = 50 K and of Trot = 10, 50, and 100 K at Tkin = 100 K. When Trot and Tkin are equal, this is equivalent to the thermal equilibrium distribution of the rotational states of the projectiles.

First, three different regimes exist: the fluorescence equilibrium at low densities (nH2O<10 cm3)$\[\left(n_{\mathrm{H}_2 \mathrm{O}}<10 \mathrm{~cm}^{-3}\right)\]$, the NLTE regime in the range 10<nH2O<106 cm3$\[10<n_{\mathrm{H}_2 \mathrm{O}}<10^6 \mathrm{~cm}^{-3}\]$, and the LTE regime at higher densities. Therefore, the NLTE regime is expected for the intermediate part of the coma, and this conclusion is relatively independent of the rotational distribution of the projectiles. We also note that the energy level populations at fluorescence equilibrium do not significantly vary with temperature. In contrast, the LTE population distribution is sensitive to the temperature.

Then, in the intermediate (NLTE) regime, the population of the energy levels is dependent on the rotational temperature of the projectile. At a lower kinetic temperature (Tkin = 50 K, Fig. 5), we found that for para-H2O, large deviations exist for the energy level populations in this regime (rotational temperature of 10 K vs. thermalized projectiles). Namely, near nH2O$\[n_{\mathrm{H}_2 \mathrm{O}}\]$ = 103 cm−3, the deviation of the state populations (solid vs. dash-dotted line) reaches 21% for the ground state 000 and 54% for the first excited state 111. In contrast, for ortho-H2O, these deviations are smaller, about 3 and 9% for states 101 and 110, respectively. However, as the temperature increases, their behavior reverses. At Tkin = 100 K (in Fig. 6), we found that when the rotational temperature is 10 K, the state populations for para-H2O are within 10% of the thermalized case (again, solid vs. dash-dotted line), but for ortho-H2O, some of these deviations are larger: 9% for the ground state 101 and 22% for the first excited state 110. These differences can be explained by the strong dependence of the excitation rate coefficients on the rotational state of H2O.

This comparison highlights the need for accurately determining and constraining the rotational population of the projectiles when modeling a cometary coma. When a thermal distribution of the projectile in the NLTE regime is assumed, this may lead to a significant under- or overestimation of the population of the targets and accordingly, of the abundance of the studied molecules. The case of Trot = 10 K is on the extreme side and could only occur in more distant comets, where Tkin is already low, for instance, about 20 K. The other cases considered here (e.g., Trot = 30 K when Tkin = 50 K) are more typical, however, and they also demonstrate a non-negligible effect of the NLTE conditions on the population of H2O states.

thumbnail Fig. 5

Relative populations of several low-lying rotational states of para- and ortho-H2O (targets) as a function of the H2O (projectile) density at Tkin = 50 K. The solid, dashed, and dash-dotted lines correspond to a rotational temperature of Trot =10, 30, and 50 K, respectively.

thumbnail Fig. 6

Same as in Fig. 5, but at Tkin = 100 K and Trot = 10, 50, and 100 K represented by solid, dashed, and dash-dotted lines, respectively.

5 Conclusions

We implemented a two-temperature model to generate a database of rate coefficients for rotational state-to-state transitions in H2O + H2O collisions that is suitable for the modeling of energy transfer under conditions when the distribution of rotational states of H2O deviates from thermodynamic equilibrium. Individual state-to-state transition cross sections from our previous work, predicted using the mixed quantum/classical theory, were re-employed here, but the procedure for averaging over the distribution of rotational states of the projectile H2O molecule was modified to permit the deviation of the rotational temperature from the kinetic temperature. Thus, the nascent rate coefficients for rotational transitions in the target H2O molecule were obtained as a 2D function of both rotational and translational temperatures.

The results of radiative transfer modeling indicate that for certain conditions that are determined by the column density, significant deviations of the rotational temperature from the kinetic temperature are possible and resultant deviations from the previously used models (which assumed Trot = Tkin) of about 10–20% are typical, but may be even larger for some emission lines of water in the part of the cometary coma where NLTE regime occurs.

These data can be used to model cometary comae and the atmospheres of icy planets where H2O + H2O collisions are important, but the conditions deviate from thermodynamic equilibrium. For the potential users of this database, we created a computer code that generates rate coefficients for 231 transitions between para-states of water and 210 transitions between ortho-states of water (in both excitation and quenching directions) for any given value of kinetic and rotational temperatures in the range from 5 to 1000 K, and for any given value of ortho-to-para ratio of water molecules. This code is available at the GitHub2. The data will be submitted to the BASECOL database.

Acknowledgements

This research was supported by NASA Astrophysics Program, grant number 80NSSC24K0208, and by French INSU-CNRS “programme national de planétologie, PNP”. D.B. acknowledges the support of Way Klingler Fellowship and Haberman research fund. F.L. acknowledges Rennes Metropole and the Institut Universitaire de France for financial support. M.A.C. was supported by the National Science Foundation under Grant No. AST-2009253. We used resources of the National Energy Research Scientific Computing Center, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-5CH11231. This research also used HPC resources at Marquette funded in part by the National Science foundation award CNS-1828649.

References

  1. Agg, P., & Clary, D. 1991, J. Chem. Phys., 95, 1037 [NASA ADS] [CrossRef] [Google Scholar]
  2. Alexander, M., Dagdigian, P. J., Werner, H.-J., et al. 2023, Comp. Phys. Commun., 289, 108761 [NASA ADS] [CrossRef] [Google Scholar]
  3. Boursier, C., Mandal, B., Babikov, D., & Dubernet, M. 2020, MNRAS, 498, 5489 [NASA ADS] [CrossRef] [Google Scholar]
  4. Brinch, C., & Hogerheijde, M. 2010, A&A, 523, A25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Buffa, G., Tarrini, O., Scappini, F., & Cecchi-Pestellini, C. 2000, ApJS, 128, 597 [NASA ADS] [CrossRef] [Google Scholar]
  6. Cochran, A. L., Levasseur-Regourd, A.-C., Cordiner, M., et al. 2015, Space Sci. Rev., 197, 9 [CrossRef] [Google Scholar]
  7. Cordiner, M., Coulson, I., Garcia-Berrios, E., et al. 2022, ApJ, 929, 38 [CrossRef] [Google Scholar]
  8. Demes, S., Lique, F., Faure, A., & van der Tak, F. F. 2023, MNRAS, 518, 3593 [Google Scholar]
  9. Dubernet, M., Boursier, C., Denis-Alpizar, O., et al. 2024, A&A, 683, A40 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Faure, A., Lique, F., & Loreau, J. 2020, MNRAS, 493, 776 [NASA ADS] [CrossRef] [Google Scholar]
  11. Faure, A., Żółtowski, M., Wiesenfeld, L., Lique, F., & Bergeat, A. 2024, MNRAS, 527, 3087 [Google Scholar]
  12. García Muñoz, A., Ramos, A. A., & Faure, A. 2024, Icarus, 415, 116080 [CrossRef] [Google Scholar]
  13. Haser, L. 1957, Bulletins de l'Académie Royale de Belgique, 43, 740 [Google Scholar]
  14. Hutson, J. M., & Le Sueur, C. R. 2019, Comput. Phys. Commun., 241, 9 [NASA ADS] [CrossRef] [Google Scholar]
  15. Hutson, J. M., & Le Sueur, C. R. 2020, MOLSCAT, Bound and Field, version 2020.0, https://github.com/molscat/molscat [Google Scholar]
  16. Lique, F., Honvault, P., & Faure, A. 2014, Int. Rev. Phys. Chem., 33, 125 [CrossRef] [Google Scholar]
  17. Mandal, B., & Babikov, D. 2023a, A&A, 671, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Mandal, B., & Babikov, D. 2023b, A&A, 678, A51 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Mandal, B., Joy, C., Bostan, D., Eng, A., & Babikov, D. 2023, J. Phys. Chem. Lett., 14, 817 [CrossRef] [Google Scholar]
  20. Mandal, B., Bostan, D., Joy, C., & Babikov, D. 2024, Comput. Phys. Commun., 294, 108938 [NASA ADS] [CrossRef] [Google Scholar]
  21. Raizer, Y. P., & Allen, J. E. 1997, Gas Discharge Physics (Berlin: Springer), 2 [Google Scholar]
  22. Ramos, A. A., & Elitzur, M. 2018, A&A, 616, A131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Van der Tak, F., Black, J. H., Schöier, F., Jansen, D., & van Dishoeck, E. F. 2007, A&A, 468, 627 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. van der Tak, F. F., Lique, F., Faure, A., Black, J. H., & van Dishoeck, E. F. 2020, Atoms, 8, 15 [NASA ADS] [CrossRef] [Google Scholar]
  25. Vorburger, A., Fatemi, S., Galli, A., et al. 2022, Icarus, 375, 114810 [NASA ADS] [CrossRef] [Google Scholar]
  26. Wirström, E., Bjerkeli, P., Rezac, L., Brinch, C., & Hartogh, P. 2020, A&A, 637, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

All Figures

thumbnail Fig. 1

Function F, which includes all factors under the integral of Eq. (11) for several transitions between the ground state 000 of para-H2O and its several low-lying rotational states (as indicated by different colors in the figure) at Tkin = 100 K. The solid, dashed, and dash-dotted lines correspond to a rotational temperature of Trot =10, 50, and 100 K, respectively.

In the text
thumbnail Fig. 2

Comparison of the rate coefficients computed here for the case of thermodynamic equilibrium Tkin = Trot (vertical) vs. thermal rate coefficients reported previously (Mandal & Babikov 2023b; horizontal) for the quenching of 231 states of para-H2O at three temperatures.

In the text
thumbnail Fig. 3

Same as Fig. 2, but for the quenching of 210 states of ortho-H2O.

In the text
thumbnail Fig. 4

2D dependence of the nonequilibrium rate coefficients for 16 strong transitions in the H2O + H2O system, as indicated in the figure. The rotational and kinetic temperatures are plotted along the horizontal and vertical axis, respectively. Color is used to show the magnitude of rate coefficients in units of cm3 s−1 in log scale as indicated in the color bar.

In the text
thumbnail Fig. 5

Relative populations of several low-lying rotational states of para- and ortho-H2O (targets) as a function of the H2O (projectile) density at Tkin = 50 K. The solid, dashed, and dash-dotted lines correspond to a rotational temperature of Trot =10, 30, and 50 K, respectively.

In the text
thumbnail Fig. 6

Same as in Fig. 5, but at Tkin = 100 K and Trot = 10, 50, and 100 K represented by solid, dashed, and dash-dotted lines, respectively.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.