# Can astronomical observations be used to constrain crucial chemical reactions? The methoxy case. SOLIS XVIII

Nadia Balucani, Cecilia Ceccarelli, Fanny Vazart, Francois Dulieu, Dimitrios Skouteris, Marzio Rosi, Fernando Pirani, Eleonora Bianchi, Paola Caselli, Claudio Codella

► To cite this version:

Nadia Balucani, Cecilia Ceccarelli, Fanny Vazart, Francois Dulieu, Dimitrios Skouteris, et al. Can astronomical observations be used to constrain crucial chemical reactions? The methoxy case. SOLIS XVIII. Monthly Notices of the Royal Astronomical Society, 2024, 528, pp.6706-6719. 10.1093/mn-ras/stae232. insu-04822477

HAL Id: insu-04822477

https://insu.hal.science/insu-04822477v1

Submitted on 6 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Can astronomical observations be used to constrain crucial chemical reactions? The methoxy case. SOLIS XVIII

Nadia Balucani<sup>1,2,3*</sup> Cecilia Ceccarelli<sup>1,2*</sup> Fanny Vazart,<sup>2</sup> Francois Dulieu,<sup>4</sup> Dimitrios Skouteris,<sup>5</sup> Marzio Rosi,<sup>6</sup> Fernando Pirani,<sup>1,6</sup> Eleonora Bianchi,<sup>7</sup> Paola Caselli<sup>8</sup> and Claudio Codella<sup>1,3</sup>

Accepted 2024 January 18. Received 2024 January 13; in original form 2023 November 10

## ABSTRACT

To understand the origin of interstellar molecules we rely on astrochemical models, the gas-phase networks of which contain ≥7000 reactions. However, just a tiny fraction of them have parameters derived in laboratory experiments. Theoretical quantum mechanical (QM) calculations can also provide this information. Unfortunately, sometimes theoretical predictions and experimental values disagree, as is the case for the paradigmatic reaction CH<sub>3</sub>OH + OH → CH<sub>3</sub>O + H<sub>2</sub>O. Both laboratory experiments and QM calculations found an unexpected increase in the rate coefficients with decreasing temperature. However, experimental and theoretical estimates of the rate coefficients diverge by up to two orders of magnitude at the low temperatures of interest in interstellar chemistry. This work aims to test whether astronomical observations can help untangle this confusing situation. To this end, we first carried out new QM calculations to derive the rate coefficients of the major destruction reaction of the methoxy radical, CH<sub>3</sub>O + H, and then we compared astronomical observations from the IRAM/NOEMA Large Programme SOLIS with astrochemical model predictions. Our new rate coefficient for the CH<sub>3</sub>O + H reaction is 5–10 times larger than that in the astrochemical data base KIDA in the 10–100 K range. When including the new methoxy destruction rate coefficients, the comparison between observations and model predictions favours the rate coefficients of the CH<sub>3</sub>OH + OH reaction from QM calculations. We conclude that QM calculations are an important alternative to laboratory experiments when it comes to the harsh conditions of interstellar objects and that astronomical observations can be used to constrain the rate coefficients of relevant reactions.

**Key words:** astrochemistry – ISM: abundances – ISM: molecules.

## 1 INTRODUCTION

Since the discovery of the first interstellar molecule in 1937 (Swings & Rosenfeld 1937), the number of detected species as well as their degree of complexity have steadily increased with time. So far, more than 300 interstellar molecules have been detected with up to 19 atoms (https://cdms.astro.uni-koeln.de/classic/molecules). At the same time, astrochemists have increased the number of species and reactions in their models to predict molecular abundances and compare them with observations. Nowadays, more than 7000 reactions are contained in the reaction networks of gas-phase models and almost twice that in grain-surface models (e.g. Agúndez & Wakelam 2013; Garrod et al. 2022; Tinacci et al. 2023).

Unfortunately, only a small fraction, less than 15 per cent, of those reactions have been included with rate coefficients and product

branching ratios (BRs) determined in laboratory experiments. This fraction is even lower if we consider reactions that have been characterized in experiments reproducing the low-temperature and/or low-pressure conditions of the interstellar medium (ISM). This is due not only to the large number of reactions involved and the relatively small community able to investigate them: many reactions are almost impossible to reproduce in terrestrial laboratories because, for instance, they involve two transient species (radicals or ions) which are difficult to generate in a controlled manner and in a sufficient amount to run reactive experiments.

A complementary approach is provided by theoretical quantum mechanical (QM) calculations, which can indeed give a deep insight into the reaction mechanisms at work by allowing one to determine the potential energy surface (PES) (see, for instance, Barone et al. (2015); Skouteris et al. (2015, 2017, 2018, 2019); Siebrand et al. (2016); Balucani et al. (2018); Gao et al. (2018); Roncero, Zanchet & Aguado (2018); Nguyen, Ruscic & Stanton (2019); Giani et al. (2023)). Once the PES has been derived, reaction rate coefficients can be calculated under the relevant pressure and temperature conditions.

© 2024 The Author(s)

Published by Oxford University Press on behalf of Royal Astronomical Society. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited.

<footnote n="1"><sup>1</sup>Dipartimento di Chimica, Biologia e Biotecnologie, Università degli Studi di Perugia, Perugia I-06123, Italy</footnote>

<footnote n="2"><sup>2</sup>Univ. Grenoble Alpes, CNRS, IPAG, F-38000 Grenoble, France</footnote>

<footnote n="3"><sup>3</sup>INAF – Osservatorio Astrofisico di Arcetri, largo E. Fermi 5, I-50125 Firenze, Italy</footnote>

<footnote n="4"><sup>4</sup>LERMA, Université de Cergy-Pontoise, F-95000 Cergy Pontoise Cedex, France</footnote>

<footnote n="5"><sup>5</sup>Master-Tec, Via Sicilia 41, I-06128 Perugia, Italy</footnote>

<footnote n="6"><sup>6</sup>Dipartimento di Ingegneria Civile e Ambientale, Università degli Studi di Perugia, Perugia I-06125, Italy</footnote>

<footnote n="7"><sup>7</sup>Excellence Cluster ORIGINS, Boltzmannstraße 2, D-85748 Garching bei München, Germany</footnote>

<footnote n="8"><sup>8</sup>Max-Planck-Institut für extraterrestrische Physik, Giessenbachstrasse 1, D-085748 Garching, Germany</footnote>

<footnote>*E-mail: nadia.balucani@unipg.it (NB); cecilia.ceccarelli@univ-grenoble-alpes.fr (CC)</footnote>
In this way, QM calculations can assist the interpretation of the experimental results (see, for instance, Siebrand et al. (2016); Gao et al. (2018); Roncero et al. (2018); Sleiman et al. (2018); Nguyen et al. (2019); Marchione et al. (2022)) or provide a reliable estimate of the rate coefficients and product BR for those cases in which no experimental data are available (see, for instance, Barone et al. (2015); Skouteris et al. (2015, 2017, 2018); Balucani et al. (2018); Giani et al. (2023))

However, it might happen that experimental and theoretical results are at odds. A notable example of this last case concerns the reaction between methanol, a widely spread interstellar molecule, and the ubiquitous OH radical. Recent kinetics experiments performed at very low temperatures in different CRESU (Cinétique de Réaction en Ecoulement Supersonique Uniforme) apparatuses have revealed an unexpected behaviour. Namely, despite the presence of a significant energy barrier, the measured rate coefficients for the reaction OH + CH₃OH increase from a value of ~10⁻¹³ cm³ mol⁻¹ s⁻¹ at 300 K to ~10⁻¹⁰ cm³ mol⁻¹ s⁻¹ at the temperature of 12 K passing through a minimum value at 220 K (Shannon et al. 2013; Gómez Martín et al. 2014; Antiñolo et al. 2016; Ocaña et al. 2019). Moreover, experimental evidence supported by a first series of theoretical calculations (Shannon et al. 2013) indicates that the main product formed in conjunction with a molecule of water is the methoxy radical (CH₃O, reaction channel 1), instead of the expected hydroxymethyl radical (CH₂OH, reaction channel 2)

$$OH + CH_3OH \rightarrow CH_3O + H_2O \quad (1)$$

$$OH + CH_3OH \rightarrow CH_2OH + H_2O. \quad (2)$$

This is counter-intuitive since channel (2) is more exothermic than channel (1), is characterized by a significantly lower energy barrier, and is dominant under the high-temperature conditions of combustion systems. The explanation for this unusual behaviour has been associated with the role played by a pre-reactive complex which is formed along the minimum energy path (at large distances where long-range forces dominate) and by the much more efficient tunnelling through the energy barrier leading to CH₃O + H₂O rather than that leading to CH₂OH + H₂O (Shannon et al. 2013).

Due to the clamor aroused by this unexpected result, the first experimental and theoretical study by Shannon et al. (2013) has been followed by numerous other experimental and theoretical studies that will be briefly analysed in Section 2. Here, we anticipate that most theoretical studies do not confirm the large effect seen in the CRESU experiments and point to an overestimate of the experimental rate coefficient at low temperatures, possibly caused by the a partial collisional stabilization of the pre-reactive complex in the conditions of CRESU experiments.

Given the abundances of OH and CH₃OH in interstellar objects, reaction (1) is expected to be the dominant formation route of the interstellar methoxy radical. A significant source of methoxy from the icy mantle of interstellar grains can be excluded because its presence in the ice is still debated being its isomer CH₂OH strongly favoured and, as in the case of most radicals, it will rather react with other species before undergoing thermal desorption. A detailed discussion on these two aspects is provided in Appendix A and the conclusion is that ice chemistry is very likely a minor player. Because of the relevance of methoxy as a possible precursor of more complex molecules (e.g. Garrod, Widicus Weaver & Herbst 2008; Vasyunin & Herbst 2013; Balucani, Ceccarelli & Taquet 2015), after its first detection towards a cold object (B1-b, ~10 K) by Cernicharo et al. (2012) it has been searched for in numerous sources. However, successful detections remain sparse and always involve cold objects (e.g. Vastel

et al. 2014; Bacmann & Faure 2016; Jiménez-Serra et al. 2016; Zhou et al. 2022; Megías et al. 2023). Particularly surprising has been the failure to observe it in environments where its CH₃OH precursor is abundant and a large presence of the OH radical is expected, as is the case of the hot corinos observed by the Large Programme NOEMA Seeds Of Life In Space (SOLIS; Ceccarelli et al. 2017).

The main aim of this paper is to use the abundance of methoxy and methanol measured in astronomical objects to help solve the earlier-described discrepancy between laboratory and theoretical estimates of the crucial reaction OH + CH₃OH → CH₃O + H₂O. To this end, we present model predictions of the gas-phase abundance of methoxy and methanol by using the recent theoretical and experimental values of rate coefficients and product BR for reactions in equations (1) and (2).

As we are going to see, to address the missing detection of the methoxy radical in several environments, it has also been crucial to reconsider its main destruction route. As discussed in Vasyunin & Herbst (2013), Balucani et al. (2015), and Antiñolo et al. (2016, hereinafter AAJ2016), the major destruction route of methoxy in interstellar conditions is its reaction with abundant hydrogen atoms

$$CH_3O + H \rightarrow \text{products,} \quad (3)$$

and to a much lesser extent, with atomic oxygen. Unfortunately, the rate coefficients for both reactions are poorly constrained under ISM conditions. In particular, in the case of reaction in equation (3), old kinetics experiments at room temperature and above derived a rate coefficient of ~3 × 10⁻¹¹ cm³ s⁻¹ (Hoylemann, Sievert & Wagner 1981; Dobe, Berces & Szilagyi 1991) with the dominant channel being that leading to H₂CO + H₂. More recent QM calculations by Xu, Xu & Lin (2007) show the presence of two major channels for the CH₃O + H reaction, leading to H₂CO + H₂ and CH₃ + OH, respectively. The first channel presents an activation barrier of 0.11 kcal mol⁻¹ (0.46 kJ mol⁻¹, ~ 55 K), which is easily overcome at room temperature. However, at the cold temperatures of the ISM, this channel is likely inhibited, leading to a rate coefficient substantially different from that measured by experiments at room temperature. Since Xu et al. (2007) only derived rate coefficients at temperatures larger than 200 K, we carried out new QM calculations of the PES and evaluated the rate coefficient of the reaction in equation (3) at the low temperatures of the ISM.

The manuscript is organized as follows. In Section 2, we critically review the available experimental and theoretical studies on the reaction OH + CH₃OH. In Section 3, we review the astronomical observations available in the literature. In Section 4, we describe the model predictions and their comparison with the astronomical observations. Finally, in Section 5 we discuss the implications of what we find, and in Section 6 we summarize our conclusions. The grain-surface role in the methoxy chemistry is discussed in the Appendix Section A, and the details of our theoretical calculations on the methoxy destruction route, CH₃O + H, are reported in the Appendix Section B.

## 2 PREVIOUS EXPERIMENTAL AND THEORETICAL CALCULATIONS ON THE REACTION OH + CH₃OH

As already mentioned, the OH + CH₃OH reaction was never believed to play any role in interstellar chemistry because it is characterized by an energy barrier of 15.2 kJ mol⁻¹ when considering channel (1) and 4.2 kJ mol⁻¹ when considering channel (2) (Xu & Lin 2007). For the same reason, this reaction had not been investigated in the past with the CRESU technique (the only technique able to provide reaction
rate coefficients at the temperature of interest in the ISM) because of the rule of thumb that a reaction with a rate coefficient below $10^{-12} \text{ cm}^3 \text{ s}^{-1}$ at room temperature has a negligible rate coefficient at very low temperatures. Nevertheless, non-Arrhenius reactions can have surprising behaviour (Sims et al. 1993; Sabbah et al. 2007; Berteloite et al. 2010). In particular, the presence of a pre-reactive complex, which does not hold a new covalent bond but is characterized by a relatively deep potential well, can have a profound influence on the outcome of a bimolecular reaction, especially at low temperatures (Georgievskii & Klippenstein 2007; Potapov et al. 2017; Heard 2018; Marchione et al. 2022; Recio et al. 2022). According to the CRESU experiments by Shannon et al. (2013) the rate coefficients for the $\text{OH} + \text{CH}_3\text{OH}$ reaction at temperatures between 22 and 88 K are almost two orders of magnitude larger than that at 200 K, with a sharp negative temperature dependence between 120 and 200 K. Shannon et al. (2013) were also able to detect $\text{CH}_3\text{O}$ as a primary reaction product at 80 K.

To rationalize their experimental finding, Shannon et al. (2013) performed master equation calculations of the total rate coefficient for the $\text{OH} + \text{CH}_3\text{OH}$ reaction starting from the electronic structure calculations of the PES carried out by Xu & Lin (2007). According to their results, the reaction features the formation of an H-bonded pre-reactive complex, characterized by a potential well as deep as $20.5 \text{ kJ mol}^{-1}$. At very low temperatures, the pre-reactive complex is long-lived and the system can reach the product asymptote by tunnelling through the barrier, while at high temperatures the lifetime of the pre-reactive complex is small and its effect on the rate coefficient becomes negligible. By considering tunnelling through the two transition states, Shannon et al. (2013) were able to qualitatively reproduce the temperature dependence of the rate coefficient in the 60–250 K range. Shannon et al. (2013) have also repeated their experiments in different pressure regimes to exclude collisional stabilization and found essentially no dependence of the measured rate coefficients on the pressure for the investigated range. Subsequently, other CRESU experiments were performed in the group of E. Jimenez in Ciudad Real (Ocaña et al. 2018). The new set of measurements extended also to lower temperatures (as low as 11.7 K, see Ocaña et al. 2019). The global picture derived by Shannon et al. (2013) was confirmed, even though the values determined in the Ciudad Real experiments are slightly lower (see Fig. 1).

The peculiar behaviour of the $\text{OH} + \text{CH}_3\text{OH}$ rate coefficients stimulated dedicated theoretical calculations. Siebrand et al. (2016) derived a new PES and performed quasi-classical trajectory (QCT) and Rice–Ramsperger–Kassel–Marcus (RRKM) calculations to clarify the reaction mechanism. Remarkably, also this work confirmed the negative temperature dependence of the rate coefficient below 200 K, but predicted much smaller rate coefficients with respect to the experimental values at low temperatures. Siebrand et al. (2016) suggested that a way to reconcile experimental data and theoretical predictions is to consider the formation of methanol dimers during the supersonic expansion through the Laval nozzle, but later Shannon et al. (2018) demonstrated that the percentage of methanol dimers in the CRESU experiments is too small to account for the observed discrepancies.

Gao et al. (2018) performed benchmark calculations of the stationary points on the $\text{OH} + \text{CH}_3\text{OH}$ PES with coupled cluster CCSD(T)-F12a and multireference CASPT2 theory and used the competitive canonical unified statistical (CCUS) model in both the low-pressure and high-pressure limits to derive rate coefficients and product BRs in the 30–2000 K range for channels (1) and (2). They also confirmed the trend with an inversion of the slope around 150 K, but the derived rate coefficient in the low-pressure limit was about

2 orders of magnitude smaller than the experimental determinations in the range 40–70 K (see Fig. 1). According to their work, the experimental values obtained with CRESU can be reproduced by assuming that the pre-reactive complex is partially stabilized by secondary collisions which can occur to some extent in the CRESU experiments but not in the conditions of interstellar clouds. Finally, they derived a temperature-dependent product BR and channel (2) is never negligible, even at the lowest temperatures investigated where it accounts for about 35 per cent of the global reaction.

Almost at the same time, Roncero et al. (2018) derived a full dimensional PES at the CCSD(T)-F12a level of calculations and performed QCT calculations to elucidate the dynamics of the $\text{OH} + \text{CH}_3\text{OH}$ reaction. In their calculations, Roncero et al. (2018) attempted to accurately describe the long-range interactions that are so important for a system like this. However, quantum effects including tunnelling through the barrier are not considered in their classical treatment. Remarkably, the analysis of the long-range interactions pointed out that the most stable configuration of the pre-reactive complex favours the H-abstraction from the $-\text{OH}$ group of methanol, thus reinforcing the suggestion by Shannon et al. (2013). However, the resulting rate coefficients are significantly lower than the experimental values between 40 and 100 K and almost temperature-independent in the range 20–300 K.

Later on, the ring polymer molecular dynamics (RPMD) method was used by del Mazo-Sevillano et al. (2019) on the same full-dimensional PES to derive the rate coefficient from 20 to 1200 K. In agreement with the experimental data, del Mazo-Sevillano et al. (2019) were able to reproduce the V-shaped dependence on temperature with a minimum at 200–300 K by combining their trapping rate coefficients obtained with RPMD and Transition State Theory (TST) calculations by Ocaña et al. (2019). Further theoretical studies by Ocaña et al. (2019) and Nguyen et al. (2019) provided values of rate coefficients, again smaller than those derived from the experiments when considering the low-pressure limit (see Fig. 1). Naumkin et al. (2019) tried to characterize the effects of secondary collisions and concluded that the experimental results below 100 K can also be explained by assuming the zero-pressure limit. New electronic structure calculations were also performed by Györi & Czakó (2023), while a new PES was calculated by Song & Li (2023) at UCCSD(T)-F12a/AVTZ level and used to derive rate coefficients via QCT calculations, but only for T>200 K.

For completeness, we mention here also the study by Hernandez et al. (2015), where the structure of the hydrogen-bonded intermediate along the H-abstraction path was investigated experimentally, and the recent study by Benitez et al. (2021) on the characteristics of the exit channel.

To summarize, there are several values for the rate coefficients of reaction (1) as a function of the temperature which differ by up to two orders of magnitude in the 10–200 K temperature range (see Fig. 1). Essentially, these values are confined within the experimental determinations by Shannon et al. (2013) plus AAJ2016 (upper limit) and those calculated by Gao et al. (2018) (lower limit). Therefore, for our modelling we will consider the two extreme values of the rate coefficients to verify if astronomical observations can discriminate among them (see Section 4).

## 3 ASTRONOMICAL OBSERVATIONS

Methoxy was first detected towards a cold (~10 K) object, B1-b, by Cernicharo et al. (2012). This detection was followed by a few others, always towards cold (~10 K) objects (Vastel et al. 2014; Bacmann & Faure 2016; Jiménez-Serra et al. 2016; Zhou et al. 2022; Megías
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Rate coefficients for the reaction (1) OH + CH₃OH → CH₃O + H₂O, as a function of temperature (within the 1–300 K range). Symbols: circles and squares are experimental points by Shannon et al. (2013), AAJ2016, Gómez Martín et al. (2014), and Ocaña et al. (2019). Continuum lines: theoretical predictions by Siebrand et al. (2016), Gao et al. (2018), Ocaña et al. (2019), and Nguyen et al. (2019). The green continuum line represents the rate coefficients calculated using the α, β, and γ parameters derived by AAJ2016. Dashed and dotted blue lines represent the rate coefficients for channels (1) and (2), respectively, in the theoretical work by Gao et al. (2018).</figcaption></figure>

et al. 2023). Methoxy was also searched for towards a handful of warm (50–100 K) objects targeted by the sensitive observations of the IRAM-NOEMA Large Programme SOLIS (Ceccarelli et al. 2017). However, the search ended up only in upper limits on the abundance of methoxy in those objects (Ceccarelli et al. 2017; Bianchi et al. 2019). Please note that the upper limits towards the warm objects, listed in Table 1, are up to a factor 10 larger with respect to the column densities estimated for the cold objects. This is due to two reasons (1) the larger temperature and (2) the smaller sizes of IRAS4A and SVS13A with respect to the former (which introduces a beam dilution). We notice that the absence of methoxy detection in the SOLIS observations of warm objects is confirmed by the search of methoxy over the entire 1, 2, and 3 mm bands observed in the warm sources targeted by the IRAM-30m Large Programme ASAI (Astrochemical Surveys At Iram: Lefloch et al. 2018). The detections and upper limits of cold and warm objects are summarized in Table 1.

## 4 MODEL SIMULATIONS AND COMPARISON WITH OBSERVATIONS

The detection of methoxy in cold environments and non-detection in warm ones provide important information that we try to interpret by using model simulations and comparing their predictions with the

observations. Specifically, we run models that simulate the conditions in two objects of Table 1

(i) Warm object, L1157-B1 (Section 4.2): L1157-B1 is a well-studied site of shocked gas (e.g. Bachiller et al. 2001; Codella et al. 2017, 2020; Lefloch et al. 2017). We chose to compare the model predictions with the observations towards L1157-B1 because they provide the most stringent constraint among the warm/hot objects. Indeed, in addition to the standard comparison between predicted and observed abundances, the shocked gas site provides also a constraint on the time the comparison should be made, namely in the time lapse corresponding to the approximate age of the shocked gas, which is known from previous studies (Podio et al. 2016).

(ii) Cold object, B1-b (Section 4.3): B1-b is the cold core towards which AAJ2016 tested model predictions using their experimental values for the rate coefficient of channel (1) Cernicharo et al. (2012). Since AAJ2016 used the old value listed in the KIDA data base for the rate coefficients of the reaction CH₃O + H (see Table 2 and Fig. B3), it is important to see also the effect of the new values we have obtained for this process on their model.

We used the code GRAINOBLE (Taquet, Ceccarelli & Kahane 2012; Ceccarelli et al. 2018) to provide the predicted abundances for the warm and cold objects. Despite GRAINOBLE is a gas-grain time-dependent astrochemical model, in this work we used
<table><caption>**Table 1.** Observed column densities of methanol and methoxy (columns 2 and 3), as well as their ratio (column 4) in several objects (source name in column 1). Column 5 reports the reference to the articles quoting the observed values. References: (1) Öberg et al. (2010); (2) Cernicharo et al. (2012); (3) Vastel et al. (2014); (4) Jiménez-Serra et al. (2016); (5) Bacmann & Faure (2016); (6) Zhou et al. (2022); (7) Megías et al. (2023); (8) Codella et al. (2020); (9) Ceccarelli et al. (2017); (10) Taquet et al. (2015); (11) Bianchi et al. (2019).</caption><thead><tr><th>Source</th><th>N(CH<sub>3</sub>OH)<br>(cm<sup>-2</sup>)</th><th>N(CH<sub>3</sub>O)<br>(cm<sup>-2</sup>)</th><th>CH<sub>3</sub>OH/CH<sub>3</sub>O</th><th>References</th></tr></thead><tbody><tr><td colspan="5">Cold objects: T~10 K</td></tr><tr><td>B1-b</td><td>4 × 10<sup>14</sup></td><td>7 × 10<sup>11</sup></td><td>600</td><td>(1), (2)</td></tr><tr><td>L1544-CH<sub>3</sub>OH</td><td>3 × 10<sup>13</sup></td><td>4 × 10<sup>11</sup></td><td>75</td><td>(3), (4)</td></tr><tr><td>L1689B</td><td>1 × 10<sup>14</sup></td><td>1 × 10<sup>12</sup></td><td>100</td><td>(5)</td></tr><tr><td>L1495A-S</td><td>6 × 10<sup>13</sup></td><td>9 × 10<sup>11</sup></td><td>70</td><td>(5)</td></tr><tr><td>L1709A</td><td>7 × 10<sup>13</sup></td><td>5 × 10<sup>11</sup></td><td>140</td><td>(5)</td></tr><tr><td>L1709A</td><td>4 × 10<sup>13</sup></td><td>6 × 10<sup>11</sup></td><td>70</td><td>(5)</td></tr><tr><td>GCC869</td><td>2 × 10<sup>13</sup></td><td>5 × 10<sup>11</sup></td><td>40</td><td>(6)</td></tr><tr><td>GCC1390</td><td>1.5 × 10<sup>13</sup></td><td>3 × 10<sup>11</sup></td><td>50</td><td>(6)</td></tr><tr><td>GCC4149</td><td>1 × 10<sup>13</sup></td><td>3 × 10<sup>11</sup></td><td>30</td><td>(6)</td></tr><tr><td>L1517B</td><td>9 × 10<sup>12</sup></td><td>3 × 10<sup>11</sup></td><td>30</td><td>(7)</td></tr><tr><td colspan="5">Outflow shock site: T~90 K</td></tr><tr><td>L1157-B1</td><td>1.3 × 10<sup>16</sup></td><td>≤4.5 × 10<sup>12</sup></td><td>≥3000</td><td>(8), (9)</td></tr><tr><td colspan="5">Warm objects: T~100 K</td></tr><tr><td>IRAS4A</td><td>5 × 10<sup>18</sup></td><td>≤2 × 10<sup>12</sup></td><td>≥5000</td><td>(9), (10)</td></tr><tr><td>SVS13A</td><td>6 × 10<sup>18</sup></td><td>≤3 × 10<sup>15</sup></td><td>≥2000</td><td>(9), (11)</td></tr></tbody></table>

<table><caption>**Table 2.** Formation (upper table) and destruction (middle table) rate coefficients of the dominant reactions in the formation and destruction of methoxy, used in the astrochemical modelling of Section 4. For the methoxy formation, two cases are reported, from AAJ2016 and Gao et al. (2018), respectively. For the destruction, we list the values from the KIDA data base, also used by AAJ2016, and those calculated in this work. The bottom table reports the CH<sub>2</sub>OH formation rate coefficients from Gao et al. (2018). The rate coefficients *k* are given as a function of the temperature *T* by the usual expression, in different temperature ranges obtained by fitting the theoretical curves (see text): *k* = *α* (T/300K)<sup>β</sup> exp(-γ/T). Columns 4–6 report *α*, *β*, and *γ* (*γ* has been set equal to 0 because these reactions are all barrierless in the considered temperature range (10–150 K)), while column 7 indicates the fit's validity interval. The last two columns report the rate coefficient values at 10 and 90 K, respectively, the two temperatures used for the astrochemical modelling of Section 4. Note that the values at 10 K have been extrapolated using the fit in the interval of column 7.</caption><thead><tr><th>Reaction: author</th><th>Reactants</th><th>Products</th><th>α (cm<sup>3</sup> s<sup>-1</sup>)</th><th>β</th><th>γ</th><th>T range</th><th>k(10 K)</th><th>k(90 K)</th></tr></thead><tbody><tr><td colspan="9">Methoxy formation</td></tr><tr><td>1a: AAJ2016</td><td>OH + CH<sub>3</sub>OH</td><td>CH<sub>3</sub>O + H<sub>2</sub>O</td><td>3.60 × 10<sup>-12</sup></td><td>-1.0</td><td>0</td><td>20–70</td><td>1.08 × 10<sup>-10</sup></td><td>1.20 × 10<sup>-11</sup></td></tr><tr><td rowspan="3">1b: Gao et al. (2018)</td><td rowspan="3">OH + CH<sub>3</sub>OH</td><td rowspan="3">CH<sub>3</sub>O + H<sub>2</sub>O</td><td>1.46 × 10<sup>-13</sup></td><td>-0.11</td><td>0</td><td>100–150</td><td>1.87 × 10<sup>-12</sup></td><td>1.75 × 10<sup>-13</sup></td></tr><tr><td>6.77 × 10<sup>-14</sup></td><td>-0.79</td><td>0</td><td>45–100</td><td></td><td></td></tr><tr><td>3.27 × 10<sup>-14</sup></td><td>-1.19</td><td>0</td><td>30–45</td><td></td><td></td></tr><tr><td colspan="9">Methoxy destruction</td></tr><tr><td>3a: KIDA</td><td>CH<sub>3</sub>O + H</td><td>H<sub>2</sub>CO + H<sub>2</sub></td><td>3.00 × 10<sup>-11</sup></td><td>0</td><td>0</td><td>10–800</td><td>3.00 × 10<sup>-11</sup></td><td>3.00 × 10<sup>-11</sup></td></tr><tr><td rowspan="2">3b: this work</td><td rowspan="2">CH<sub>3</sub>O + H</td><td rowspan="2">CH<sub>3</sub> + OH</td><td>1.70 × 10<sup>-10</sup></td><td>-0.29</td><td>0</td><td>80–300</td><td>1.59 × 10<sup>-10</sup></td><td>2.40 × 10<sup>-10</sup></td></tr><tr><td>3.25 × 10<sup>-10</sup></td><td>0.21</td><td>0</td><td>≤80</td><td></td><td></td></tr><tr><td colspan="9">CH<sub>2</sub>OH formation</td></tr><tr><td rowspan="3">2: Gao et al. (2018)</td><td rowspan="3">OH + CH<sub>3</sub>OH</td><td rowspan="3">CH<sub>2</sub>OH + H<sub>2</sub>O</td><td>2.11 × 10<sup>-13</sup></td><td>0.40</td><td>0</td><td>100–150</td><td>6.89 × 10<sup>-13</sup></td><td>1.37 × 10<sup>-13</sup></td></tr><tr><td>8.75 × 10<sup>-14</sup></td><td>-0.37</td><td>0</td><td>45–100</td><td></td><td></td></tr><tr><td>3.45 × 10<sup>-14</sup></td><td>-0.88</td><td>0</td><td>30–45</td><td></td><td></td></tr></tbody></table>

the version that only considers gas-phase reactions, GRAINOBLE-GAS, which means that the grain-surface chemistry is not taken into account. The reason for this choice and the specific details used for the modelling of the cold and warm objects are described in the relevant Sections 4.3 and 4.2, respectively. In the present simulation, following AAJ2016, we consider that methoxy is mostly formed and destroyed by gas-phase reactions and that grain-surface reactions play a minor role. A full discussion explaining the reasons for this choice is reported in Appendix A.

We adopted the chemical network described in Section 4.1, where we considered the different values of the major reactions of formation and destruction of methoxy.

## 4.1 Chemical network

The chemical network involving the formation and destruction of methoxy is, obviously, a crucial element of the modelling. We adopted the gas-phase network, based on the KIDA 2014 network (Wakelam et al. 2012), that our group has updated and used in several previous works and that it is fully described and revised in Tinacci et al. (2023).

Concerning methoxy, in addition to the standard destruction by molecular ions (He<sup>+</sup>, H<sub>3</sub><sup>+</sup>, and HCO<sup>+</sup>), we employed the rate coefficients of the major reactions of formation (reaction 1) and destruction (reaction 3), summarized in Table 2 and discussed in the following sections. Note that, for completeness, we also included in
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Formation rate constant $k_{\text{formation}}$ as a function of temperature T by AAJ2016 (blue curve) and Gao et al. (2018) (red curve), respectively. The continuum lines show the values quoted by the above authors while the dashed lines are extrapolated values.</figcaption></figure>

the table the rate coefficients by Gao et al. (2018) for the CH$_2$OH formation from reaction (2).

### 4.1.1 Methoxy formation

As discussed earlier, there is a substantial disagreement between experimental values and theoretical predictions for the rate coefficients of reaction (1). If we consider the values in the range 30–100 K by Gao et al. (2018) and those by AAJ2016 in the 20–70 K interval, they differ by a factor of ca. 40 at 30 K and ca. 70 at 70 K. Therefore, we have run our models considering both the rate coefficients derived by AAJ2016 and those derived by Gao et al. (2018) to verify in which case the model predictions reproduce the observations. Table 2 lists the $\alpha$, $\beta$, and $\gamma$ coefficients for these two cases. Please note that we fitted the data in the 30–150 K range by Gao et al. (2018) to derive the relevant $\alpha$ and $\beta$ coefficients (in the fitting we imposed $\gamma = 0$ for barrierless reactions) and extrapolated the curve down to 10 K. Likewise, we extrapolated the value up to 90 K by using the AAJ2016 curve in the 20–70 K interval. Fig. 2 shows the original values and the extrapolated ones. Also, to be noted that, in the considered temperature range (10–150 K) and for both cases, the OH + CH$_3$OH reaction behaves as a barrierless reaction since it exhibits an increasing value of the rate coefficient with decreasing temperature. From this, the choice of using $\gamma = 0$.

### 4.1.2 Methoxy destruction

As also discussed by AAJ2016 and found by our own modelling a posteriori, the major methoxy destruction route is the reaction CH$_3$O + H. AAJ2016 adopted the KIDA value, introduced by Ruaud et al. (2015), which is based on old kinetics experiments obtained at room temperature (Hoyermann et al. 1981; Dobe et al. 1991, see Appendix B). Those experiments identified H$_2$CO and H$_2$ as the dominant reaction products with the formation of CH$_3$ and OH being minor. Collision-free experiments on the decomposition of excited CH$_3$OH formed by the reaction O($^1$D)+CH$_4$ suggested, instead, that the CH$_3$ + OH are the most important products (Yang 2006).

QM calculations by Xu et al. (2007) suggested that the rate coefficient is significantly larger (in the gas kinetics limit, ca. ×

$10^{-10} \text{ cm}^3 \text{ s}^{-1}$) at room temperature and that the major products are CH$_3$ and OH. Our results (reported in Appendix B), fully confirmed the calculations by Xu et al. (2007) and also in the more extended range of T that we have investigated the reaction remains fast with rate coefficients larger than the KIDA ones by a factor of 50 (see Fig. B3). In our network, we have used the new values (see Appendix B) that are reported in Table 2 (also the KIDA values are shown for comparison).

## 4.2 Warm environments: L1157-B1

As described in Section 3, methoxy is not detected in warm environments. The most stringent constraints are provided by the case of the shock site L1157-B1, where, in addition to the upper limit on its abundance, observations also provide an indication of the time at which methoxy should not be detected.

### 4.2.1 Model description

In order to simulate the L1157-B1 case, we ran the gas-phase version of GRAINOBLE-GAS with the parameters constrained by several previous works, and lately summarized in Codella et al. (2017, 2020) and Tinacci et al. (2023). Briefly, the shocked gas is modelled by a sudden increase of the temperature, density, and abundance of the molecules previously frozen onto the grain mantles. This is obtained by running GRAINOBLE-GAS in a first step to compute the gas abundances of the cold molecular cloud gas and then a second step to simulate the shock. At the passage of the shock, we assume that the gas temperature instantaneously reaches 90 K, the number density of total H nuclei increases to $8 \times 10^5 \text{ cm}^{-3}$, and several (previously frozen) species are injected into the gas phase with the following abundances with respect to the H nuclei: $3 \times 10^{-6}$ for CH$_3$OH, $2 \times 10^{-4}$ for H$_2$O, $8 \times 10^{-8}$ for CH$_3$CH$_2$, $3 \times 10^{-5}$ for CO$_2$, $2 \times 10^{-6}$ for OCS, $1 \times 10^{-6}$ for H$_2$CO, and $2 \times 10^{-5}$ for NH$_3$. Finally, the cosmic rays ionization rate is $6 \times 10^{-16} \text{ s}^{-1}$. We emphasize that all these values have been adopted based on observations, which have been interpreted via radiative transfer non-LTE modelling (e.g. gas temperature and density: Lefloch et al. 2012) or which are compared with astrochemical model predictions (e.g. cosmic rays ionization rate: Podio et al. 2014) or directly (e.g. abundances of injected species: Codella et al. 2017, for details).

The chemical composition of the gas is then computed as a function of time over an interval from one hundred to two thousand years, which is the time-scale of the shocks in L1157-B1, as derived by numerous observations (e.g. Podio et al. 2016; Codella et al. 2017).

### 4.2.2 Model predictions and comparison with observations

The results of the simulations are shown in Fig. 3. The predictions obtained with the different rate coefficients of Table 2 are very different and definitively provide strong constraints.

First, when using the AAJ2016 and the KIDA values, the methoxy predicted abundance is much higher than the measured upper limit by more than one order of magnitude. In addition, the predicted methanol abundance becomes incompatible with the measured one at 1500 yr, so that the CH$_3$OH/CH$_3$O abundance ratio is in total disagreement with the model predictions by more than a factor 30 at time ≥1000 yr, which is the approximate age of the shock (Podio et al. 2016). The situation is slightly better when we use our methoxy destruction rate coefficient with the AAJ2016 rate coefficient for reaction (1), because of the increase of the methoxy destruction rate
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Comparison between astrochemical model predictions (Section 4) and observations (Section 3) towards the shocked gas site L1157-B1 (upper panels) and the pre-stellar core B1-b (bottom panels). The lower part of each panel shows the computed abundance of CH<sub>3</sub>O (red curves) and CH<sub>3</sub>OH (blue curves) as a function of time, whereas the upper part the CH<sub>3</sub>OH/CH<sub>3</sub>O abundance ratio (green curves). The three curves show the model predictions obtained considering: solid lines, (Gao et al. 2018) (coefficients 1b of Table 2) and our new computations (coefficients 2b of Table 2); dashed lines, AAJ2016 (coefficients 1a of Table 2) and the KIDA value (coefficients 2a of Table 2); dotted-dashed lines, AAJ2016 (coefficients 1a of Table 2) and our new computations (coefficients 2b of Table 2). The dashed areas show the range of observed values reported in Table 1. The H nuclei column density used to derive the CH<sub>3</sub>OH and CH<sub>3</sub>O abundances is taken from Lefloch et al. (2012) for L1157-B1 and Öberg et al. (2010) for B1-b, respectively.</figcaption></figure>

by a factor ~10. However, considering simultaneously the methoxy and methanol abundances, the model predictions largely fail to reproduce the observations. On the contrary, the predictions obtained using the Gao et al. (2018) values for the rate coefficient of reaction (1) and our values for reaction (3) are in very good agreement with the observations.

In conclusion, the comparison between the astrochemical model predictions and the observations towards L1157-B1 definitively favours the methoxy formation rate coefficients computed by Gao et al. (2018) against the ones measured by AAJ2016.

## 4.3 Cold environments: B1-b

In this section, we consider the source B1-b, where methoxy was detected for the first time Cernicharo et al. (2012) and for which AAJ2016 carried out a detailed modelling using their experimental rate coefficient for reaction (1). Our goal is to verify the impact of the other values of the rate coefficients, rather than modelling again the chemistry of that source.

### 4.3.1 Model description

In order to simulate the B1-b environment, we ran the gas-phase version of GRAINOBLE-GAS, using exactly the parameters of the model by AAJ2016. We first ran the model with the same reaction rate coefficients used by AAJ2016 and then with the new ones of Table 2. Specifically, for the physical parameters, we used a temperature of 10 K, a number density of total H nuclei of $2 \times 10^4 \text{ cm}^{-3}$, a visual extinction of 30 mag, a cosmic ray ionization rate of $H_2$ equal to $1.3 \times 10^{-17} \text{ s}^{-1}$ and the so-called ‘low-metal’ elemental abundances.<fnref n="1" /> We then impose that the CH<sub>3</sub>OH abundance is $6 \times 10^{-5}$ times the CO abundance, as done by AAJ2016 in order to reproduce the relatively large observed methanol abundance. Since the actual methanol abundance has been determined to be ~$2 \times 10^{-9}$ (Öberg et al. 2010) in B1-b, that is about a factor of ten lower than that predicted by the AAJ2016 model with this assumption (see their fig. 4), we also run a case where CH<sub>3</sub>OH/CO is ten times lower. We anticipate that, while the CH<sub>3</sub>OH/CH<sub>3</sub>O abundance ratio does not change in the two cases, the absolute CH<sub>3</sub>O abundance decreases by a factor ~10 as expected.

### 4.3.2 Model predictions and comparison with observations

The results of the simulations are shown in Fig. 3. First, we verified that we reproduce the results by AAJ2016, obtained when using the same rate coefficients for reactions (1) and (3). As in AAJ2016, the CH<sub>3</sub>O and CH<sub>3</sub>OH abundances, as well as their ratio, are perfectly in line with the observations. However, when introducing our new values for reaction (3), the methoxy abundance predicted by the model decreases, because the new value is a factor of 5 larger than that reported in KIDA. Consequently, the CH<sub>3</sub>OH/CH<sub>3</sub>O abundance ratio is barely compatible with the observed one. When also the much smaller rate coefficient derived by Gao et al. (2018) for reaction (1) is considered, the CH<sub>3</sub>OH/CH<sub>3</sub>O abundance ratio is about a factor of 2 larger than that observed. To be noted that, in this case, we used an extrapolated value rather than the real outcome of Gao et al. (2018) calculations because the lowest temperature considered in their work was 30 K. Given the slope of $k(T)$ in that range, the extrapolation might be wrong.

It is also important to emphasize that the model, which assumes a constant ratio between CO and CH<sub>3</sub>OH as AAJ2016, is rather simplified. When considering these uncertainties the disagreement by a factor 2 may be relatively modest (see Section 5).

<footnote n="1">¹See table 3 in Agúndez & Wakelam (2013), where the C, O, and N abundances are those measured by Jenkins (2009) towards ζ Oph.</footnote>
To summarize, the comparison of the observations towards B1-b with the model predictions obtained using the different rate coefficients of Table 2 does not seem to provide constraints strong enough to discriminate among the two sets of values. Please note that the model has been tailored to B1-b, so that the comparison with the other sources of Table 1 is not significant in this context.

## 5 DISCUSSION

Although methoxy is a relatively small organic molecule with five atoms, in the last few years it has conquered the scene of astrochemistry because it is the unexpected product of the reaction $OH + CH_3OH$, for which a surprisingly large rate coefficient has been measured at low temperatures in spite of the presence of an energy barrier (Shannon et al. 2013, AAJ2016). This indicates that under certain conditions, namely when a pre-reactive complex with a certain stability is formed, neutral-neutral reactions with a sizeable activation energy can still occur because of quantum effects like tunnelling. However, as discussed in Section 2, successive theoretical QM calculations derived much smaller rate coefficients at low temperatures, especially below 70 K (Gao et al. 2018; Roncero et al. 2018; Nguyen et al. 2019; Ocaña et al. 2019). Since $CH_3O$ is also a precursor of important iCOMs (Balucani et al. 2015; Herbst 2021), it is important to understand its chemistry in the ISM.

Astronomical observations have detected methoxy in several cold (~10 K) objects (Cernicharo et al. 2012; Bacmann & Faure 2016; Jiménez-Serra et al. 2016; Zhou et al. 2022; Megías et al. 2023) but failed to detect it in warm (~70–100 K) objects (Ceccarelli et al. 2017; Bianchi et al. 2019), against expectations, as methanol is much more abundant in the latter than in the former. A first glance, this is puzzling but, exactly because of that, solving the puzzle can provide precious information on methoxy formation and cold chemistry in general. In the past, astronomical observations have been successful in identifying or even correcting frequencies of molecular rotational transitions (e.g. Cernicharo et al. 2018). Inspired by those works, we might try to use astronomical observations to constrain the methoxy formation rate coefficients.

The method, in this case, is to compare the observed abundance ratio of methoxy and its parent molecule, methanol, with astrochemical model predictions, since absolute abundances are always difficult to obtain and are intrinsically uncertain by about a factor of 2–3. However, to have reliable model predictions one needs not only reliable rate coefficients for the formation routes, but also for the destruction ones. Since the rate coefficients for the major methoxy destruction reaction, $CH_3O + H$, were not available in the literature for temperatures lower than 200 K, we carried out dedicated ab initio calculations of the PES and kinetics calculations to derive them in the temperature range of interest (Appendix B).

Armed with the new rate coefficients for the $CH_3O + H$ reaction, we ran the astrochemical model GRAINOBLE-GAS and compared the model results with the observations towards a warm object, L1157-B1, and a cold object, B1-b. We considered the methoxy formation rate coefficients derived from the AAJ2016 laboratory experiments and from the QM calculations of Gao et al. (2018). The experimental and theoretical rate coefficients differ by a factor ~100 at the temperature of the warm (90 K) and cold (10 K) objects (Table 2).

The comparison between the astrochemical model predictions and the observations towards L1157-B1 strongly favours the values derived by Gao et al. (2018) (Fig. 3). The model using the AAJ2016 rate coefficients fails to reproduce the observations, namely the methanol and methoxy abundances as well as their ratio, by more than a factor of 7, whereas that using the Gao et al. (2018) rate coefficients provides

predictions in very good agreement with the observed values. We emphasize that these results apply even in the case in which some methoxy is synthesized on the grain surfaces (see the discussion in Section A) and injected into the gas. Indeed, the constraint comes from the too-large abundance of methoxy predicted by the $OH + CH_3OH$ reaction compared to the observed one. Since this is a direct consequence of the methanol abundance in the gas phase, the abundance of the methoxy synthesized in the gas phase cannot be lower. The methoxy formation rates have been derived up to 70 K by AAJ2016, but the successive experiments by Ocaña et al. (2019) and Gómez Martín et al. (2014) performed at temperatures between 70 and 90 K (the temperature of relevance in the modelling) are even larger than the extrapolated value and, therefore, model predictions will be even worse in this case. Therefore, even considering the uncertainty associated with both sets of values, the rate coefficients derived by Gao et al. (2018) are much more in line with the observations.

The comparison of the model predictions with the observations of the cold object, B1-b, is much less constraining. Adopting the AAJ2016 rate coefficients instead of those proposed by us provides better results, as the $CH_3OH/CH_3O$ abundance ratio is predicted to be a factor of 2 too large with respect to the observed one in the second case. However, when considering the uncertainty associated with the extrapolation (the calculations by Gao et al. (2018) stop at 30 K) and the simplified model used for that source, we cannot consider this a strong constraint. Finally, it is important to emphasize that, while the temperature and density, as well as the initial species abundances of L1157-B1 are well-established (e.g. Codella et al. 2020; Tinacci et al. 2023, and the herein references), the same cannot be said of the B1-b model, as it is a much less studied source.

In summary, the comparison between astronomical observations and astrochemical model predictions obtained using experimental and theoretical rate coefficients for reaction (1) suggests that the latter are more compatible with astronomical observations. These results clearly indicate that astronomical observations can discriminate between chemical data that are in disagreement.

## 6 CONCLUSIONS

Methoxy has recently become the focus of many astrochemical experimental, theoretical, and observational studies after its first discovery in a cold object by Cernicharo et al. (2012) and the unexpected trend with the temperature of the rate coefficient of the $CH_3OH + OH \rightarrow CH_3O + H_2O$ reaction. In this work, we have reported a comprehensive review of the experimental and theoretical studies. Significant differences emerge among the set of experimental values and most of the theoretical predictions.

To verify whether astronomical observations can contribute to disentangling this confusing situation, we ran astrochemical models using the two extreme values of the rate coefficients for the dominant formation route of methoxy: the experimental values by AAJ2016 and the theoretical values by Gao et al. (2018). Since the $CH_3O$ major destruction route, $CH_3O + H$, was never studied at low (≤200 K) temperatures, we carried out new ab initio calculations to evaluate the rate coefficients in the 5–200 K range.

To compare model predictions with observations, we considered the $CH_3O$ abundances in cold objects (Cernicharo et al. 2012; Bacmann & Faure 2016) and its upper limits in warm regions (Ceccarelli et al. 2017). In the case of warm objects, the model predictions obtained using the theoretical rate coefficients by Gao et al. (2018) are in agreement with the observations, whereas those obtained using the value by AAJ2016 fail to reproduce the observations by about a factor of 10. In the case of cold objects, the situation is less clear, as
the predictions with the values derived by Gao et al. (2018) are only a factor of 2 smaller than the observational values. Considering that the predictions are obtained by extrapolating the rate coefficients by Gao et al. (2018) from 30 to 10 K, in a temperature range where the extrapolation could be very imprecise because of the steep slope, and considering that the model is quite simplistic, a factor of 2 is probably not enough to discriminate between the AAJ2016 and Gao et al. (2018) values. In conclusion, the comparison between the astrochemical model predictions and observations seems to favour the theoretical calculations by Gao et al. (2018).

## ACKNOWLEDGEMENTS

This project has received funding within the European Union's Horizon 2020 research and innovation programme from the European Research Council (ERC) for the project 'The Dawn of Organic Chemistry' (DOC), grant agreement number 741002, and from the Marie Skłodowska-Curie for the project 'Astro-Chemical Origins' (ACO), grant agreement number 811312, and the Agence Nationale de la recherche (ANR) for the SIRC project (grant no. ANR-SPV202448 20202024). The authors thank the Herla Project (http://hscw.herla.unipg.it) – 'Università degli Studi di Perugia' and the 'Dipartimento di Ingegneria Civile e Ambientale' of the University of Perugia within the project 'Dipartimenti di Eccellenza 2018-2022' for allocated computing time. Some computations presented in this paper were performed using the GRICAD infrastructure (https://gricad.univ-grenoble-alpes.fr), which is partly supported by the Equip@Meso project (reference ANR-10-EQPX-29-01) of the programme Investissements d'Avenir supervised by the Agence Nationale pour la Recherche. EB acknowledges the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany's Excellence Strategy – EXC 2094 – 390783311. C. Codella and M. Rosi also acknowledge co-funding by EU - Next Generation EU via Italian MUR PRIN2022 (2022JC2Y93 - Chemical Origins: linking the fossil composition of the Solar System with the chemistry of protoplanetary disks). Finally, we thank the anonymous referee whose comments helped to improve the clarity of this article.

## DATA AVAILABILITY

The data underlying this article are available in the article and in its online supplementary material.

## REFERENCES

* Agúndez M., Wakelam V., 2013, Chem. Rev., 113, 8710

* Antiñolo M., Agúndez M., Jiménez E., Ballesteros B., Canosa A., El Dib G., Albaladejo J., Cernicharo J., 2016, ApJ, 823, L25

* Bachiller R., Pérez Gutiérrez M., Kumar M. S. N., Tafalla M., 2001, A&A, 372, 899

* Bacmann A., Faure A., 2016, A&A, 587, 130

* Balucani N., Ceccarelli C., Taquet V., 2015, MNRAS, 449, 16

* Balucani N., Skouteris D., Ceccarelli C., Codella C., Falcinelli S., Rosi M., 2018, Mol. Astr., 13, 30

* Barone V., Latouche C., Skouteris D., Vazart F., Balucani N., Ceccarelli C., Lefloch B., 2015, MNRAS, 453, 31

* Becke A. D., 1993, The J. Chem. Phys., 98, 1372

* Benitez Y., Nguyen T. L., Parsons A. J., Stanton J. F., Continetti R. E., 2021, The J. Physical Chemistry Lett., 13, 142

* Berteloite C., Le Picard S. D., Balucani N., Canosa A., Sims I. R., 2010, Phys. Chem. Chem. Phys. (Incorporating Faraday Transactions), 12, 3666

* Bertin M. et al., 2016, ApJL, 817, L12

* Bianchi E. et al., 2019, MNRAS, 483, 1850

* Butscher T., Duvernay F., Theule P., Danger G., Carissan Y., Hagebaum-Reignier D., Chiavassa T., 2015, MNRAS, 453, 1587

* Ceccarelli C. et al., 2017, ApJ, 850, L176

* Ceccarelli C., Viti S., Balucani N., Taquet V., 2018, MNRAS, 476, 1371

* Cernicharo J., Marcelino N., Roueff E., Gerin M., Jiménez-Escobar A., Muñoz Caro G. M., 2012, ApJ, 759, L43

* Cernicharo J. et al., 2018, ApJ, 853, L22

* Cernuto A., Tosi P., Martini L. M., Pirani F., Ascenzi D., 2017, Phys. Chem. Chem. Phys., 19, 19554

* Codella C. et al., 2017, A&A, 605, 3

* Codella C. et al., 2020, A&A, 635, 17

* Dartois E., Chabot M., Id Barkach T., Rothard H., Augé B., Agnihotri A. N., Domaracka A., Boduch P., 2019, A&A, 627, 55

* Dulieu F. et al., 2013, Sci. Rep., 3, 1338

* Dobe S., Berces T., Szilagyi I., 1991, J. Ch. Soc. Faraday Trans., 87, 2331

* Fraser H. J., Collings M. P., McCoustra M. R. S., Williams D. A., 2001, MNRAS, 327, 1165

* Frisch M. et al., 2016, Inc., Wallingford Conn., 3

* Gao L., Zheng J., Fernandez-Ramos A., Truhlar D., Xu X., 2018, J. Am. Chem. Soc., 140, 2906

* Garrod R. T., Wakelam V., Herbst E., 2007, A&A 467, 1103

* Garrod R. T., Widicus Weaver S. L., Herbst E., 2008, ApJ, 682, L283

* Garrod R. T., Jin M., Matis K. A., Jones D., Willis E. R., Herbst E., 2022, ApJS, 259, 1

* Georgievskii Y., Klippenstein S. J., 2007, J. Phys. Chem. A, 111, 3802

* Giani L., Ceccarelli C., Mancini L., Bianchi E., Pirani F., Rosi M., Balucani N., 2023, MNRAS, 526, 4535

* Goerigk L., Grimme S., 2011, J. Chem. Theo. Comput., 7, 291

* Gómez Martín J. C., Caravan R. L., Blitz M. A., Heard D. E., Plane J. M. C., 2014, J. Phys. Chem. A, 118, 2693

* Grimme S., 2006, The J. Chem. Phys., 124, 034108

* Grimme S., Ehrlich S., Goerigk L., 2011, J. Comput. Chem., 32, 1456

* Györi T., Czakó G., 2023, The J. Chem. Phys., 158, 034301

* Heard D. E., 2018, Acc. Chem. Reas., 51, 2620

* Herbst E., 2021, Fr. Astr. Sp. Sci., 8, 207

* Hernandez F. J., Brice J. T., Leavitt C. M., Pino G. A., Douberly G. E., 2015, The J. Phys. Chem. A, 119, 8125

* Hoyermann K., Sievert R., Wagner H., 1981, Ber. Bunsen Phys. Chem., 85, 149

* Jenkins E. B., 2009, ApJ, 700, L1299

* Jiménez-Serra I. et al., 2016, ApJ, 830, L6

* Karssemeijer L. J., Cuppen H. M., 2014, A&A, 569, A107

* Kendall R. A., Dunning T. H., Jr, Harrison R. J., 1992, The J. Chem. Phys., 96, 6796

* Lefloch B. et al., 2012, ApJ, 757, L25

* Lefloch B., Ceccarelli C., Codella C., Favre C., Podio L., Vastel C., Viti S., Bachiller R., 2017, MNRAS, 469, 73

* Lefloch B. et al., 2018, MNRAS, 477, 4792

* Marchione D. et al., 2022, J. Phys. Chem. A, 126, 3569

* del Mazo-Sevillano P., Aguado A., Jiménez E., Suleimanov Y. V., Roncero O., 2019, The J. Phys. Chem. Lett., 10, 1900

* Megías A., Jiménez-Serra I., Martín-Pintado J., Vasyunin A. I., Spezzano S., Caselli P., Cosentino G., Viti S., 2023, MNRAS, 519, 1601

* Minissale M., Dulieu F., 2014, The J. Chem. Phys., 141, 014304

* Minissale M., Congiu E., Dulieu F., 2016a, A&A, 585, A146

* Minissale M., Dulieu F., Cazaux S., Hocuk S., 2016b, A&A, 585, 24

* Minissale M., Moudens A., Baouche S., Chaabouni H., Dulieu F., 2016c, MNRAS, 458, 2953

* Naumkin F., del Mazo-Sevillano P., Aguado A., Suleimanov Y. V., Roncero O., 2019, ACS ESC, 3, 1158

* Nguyen T. L., Ruscic B., Stanton J. F., 2019, J. Chem. Phys., 150, 084105

* Öberg K. I., Bottinelli S., Jørgensen J. K., van Dishoeck E. F., 2010, ApJ, 716, L825

* Ocaña A. J., Blázquez S., Ballesteros B., Canosa A., Antiñolo M., Albaladejo J., Jiménez E., 2018, PCCP, 20, 5865

* Ocaña A. J. et al., 2019, PCCP, 21, 6942

* Petersson G., Al-Laham M. A., 1991, The J. Chem. Phys., 94, 6081
Petersson A., Bennett A., Tensfeldt T. G., Al-Laham M. A., Shirley W. A.,
Mantzaris J., 1988, *The J. Chem. Phys.*, **89**, 2193

Pirani F., Maciel G. S., Cappelletti D., Aquilanti V., 2006, *Intern. Rev. Phys. Chem.*, **25**, 165

Pirani F., Brizi S., Roncaratti L. F., Casavecchia P., Cappelletti D., Vecchio-cattivi F., 2008, *Phys. Chem. Chem. Phys.*, **10**, 5489

Podio L., Lefloch B., Ceccarelli C., Codella C., Bachiller R., 2014, *A&A*, **565**, 64

Podio L. et al., 2016, *A&A*, **593**, 4

Pople J. A., Head-Gordon M., Raghavachari K., 1987, *The J. Chem. Phys.*, **87**, 5968

Potapov A., Canosa A., Jiménez E., Rowe B., 2017, *Angewandte Chemie*, **129**, 8742

Purvis III G. D., Bartlett R. J., 1982, *The J. Chem. Phys.*, **76**, 1910

Recio P. et al., 2022, *Nat. Chem.*, **14**, 1405

Roncero O., Zanchet A., Aguado A., 2018, *PCCP*, **20**, 25951

Ruaud M., Loison J. C., Hickson K. M., Gratier P., Hersant F., Wakelam V., 2015, *MNRAS*, **447**, 4004

Sabbah H., Biennier L., Sims I. R., Georgievskii Y., Klippenstein S. J., Smith I. W. M., 2007, *Science*, **317**, 102

Shannon R. J., Blitz M. A., Goddard A., Heard D. E., 2013, *Nature Chem.*, **5**, 745

Shannon R. J. et al., 2018, *Phys. Chem. Chem. Phys. (Incorporating Faraday Transactions)*, **20**, 8349

Siebrand W., Smedarchina Z., Martínez-Núñez E., Fernández-Ramos A., 2016, *Phys. Chem. Chem. Phys.*, **18**, 22712

Sims I. R., Queffelec J.-L., Travers D., Rowe B. R., Herbert L. B., Karthäuser J., Smith I. W. M., 1993, *Chem. Phys. Lett.*, **211**, 461

Skouteris D., Balucani N., Faginas-Lago N., Falcinelli S., Rosi M., 2015, *A&A*, **584**, 76

Skouteris D., Vazart F., Ceccarelli C., Balucani N., Puzzarini C., Barone V., 2017, *MNRAS*, **468**, 1

Skouteris D., Balucani N., Ceccarelli C., Vazart F., Puzzarini C., Barone V., Codella C., Lefloch B., 2018, *ApJ*, **854**, L135

Skouteris D., Balucani N., Ceccarelli C., Faginas Lago N., Codella C., Falcinelli S., Rosi M., 2019, *MNRAS*, **482**, 3567

Sleiman C., El Dib G., Rosi M., Skouteris D., Balucani N., Canosa A., 2018, *PCCP*, **20**, 5478

Song K., Li J., 2023, *Phys. Chem. Chem. Phys.*, **25**, 11192

Swings P., Rosenfeld L., 1937, *ApJ*, **86**, L483

Taquet V., Ceccarelli C., Kahane C., 2012, *A&A*, **538**, 42

Taquet V., López-Sepulcre A., Ceccarelli C., Neri R., Kahane C., Charnley S. B., 2015, *ApJ*, **804**, L81

Tinacci L., Ferrada-Chamorro S., Ceccarelli C., Pantaleone S., Ascenzi D., Maranzana A., Balucani N., Ugliengo P., 2023, *ApJS*, **266**, 38

Vastel C., Ceccarelli C., Lefloch B., Bachiller R., 2014, *ApJ*, **795**, L2

Vasyunin A. I., Herbst E., 2013, *ApJ*, **769**, L34

Vasyunin A. I., Caselli P., Dulieu F., Jiménez-Serra I., 2017, *ApJ*, **842**, 33

Vazart F., Latouche C., Skouteris D., Balucani N., Barone V., 2015, *ApJ*, **810**, L111

Vazart F., Ceccarelli C., Balucani N., Skouteris D., 2022, *ApJ*, **941**, L196

Wakelam V. et al., 2012, *ApJS*, **199**, 21

Wakelam V. et al., 2015, *ApJS*, **217**, 20

Woon D. E., Dunning Jr T. H., 1993, *The J. Chem. Phys.*, **99**, 1914

Xu S., Lin M. C., 2007, *J. Phys. Chem. A*, **111**, 6730

Xu K., Xu Z. F., Lin M. C., 2007, *Molecular Physics*, **105**, 2763

Yang X., 2006, *Phys. Chem. Chem. Phys.*, **8**, 205

Zhou C. et al., 2022, *A&A*, **658**, 131

## APPENDIX A: CONSIDERATIONS ON THE GRAIN-SURFACE ROLE IN THE METHOXY CHEMISTRY

We discuss here the possible role of chemistry on grains. Since methoxy is only observed in the gas phase, we need to consider the combination of both formation and desorption processes. We will therefore first discuss the thermal desorption of methoxy and

then its chemical desorption. Finally, we will focus on sputtering and photodesorption. From the analysis of these processes, it will be clear that the question of whether methoxy forms on the grain surfaces or not is completely irrelevant for the methoxy observed in the gas.

In warm environments, as soon as the grain temperature reaches 100 K, the grain ice mantles (thermally) desorb in a short time-scale: as an order of magnitude, one monolayer of water desorbs in about one year, lifting off in the gas the other species trapped in the ice matrix, including methoxy. However, before desorption, any adsorbate will be able to scan many times the surface of the grain. Indeed the desorption time-scale is $\tau_{des} = 1/v_{des} e^{-E_{binding}/T_s}$, where $E_{binding}$ is the binding energy (expressed in $K/k_B$, $k_B$ being the Boltzmann constant), $v_{des}$ is the frequency of attempts, and $T_s$ is the surface temperature of the grains. The diffusion time-scale, on the other hand, is $\tau_{diff} = 1/v_{diff} e^{-E_{diff}/T_s}$. Therefore, the number of steps $N_{steps}$ on the surface before desorption is $N_{steps} = \tau_{diff}/\tau_{des}$. It is commonly assumed that $v_{des} \sim v_{diff}$ and $E_{diff}/E_{binding}$ is in the range 0.3–0.7 (e.g. Karssemeijer et al. 2014; Minissale et al. 2016a). To have an order of magnitude, we adopt 0.5 (see e.g. Taquet et al. 2012). Adopting a water–ice desorption energy of 5600 K (Fraser et al. 2001), we find that the number of diffusion jumps at 100 K is $N_{steps} = e^{0.5 \times E_{binding}/T_s} \sim 3 \times 10^9$. This number is comparable to that deduced by Amiaud et al. (2006), who discussed the case of $D_2$ on the surface of amorphous ice and should be compared to the average number of surface adsorption sites of a typical grain (radius $\sim 0.1$ µm) of $10^6$. Therefore, in general, any physisorbed molecule sweeps widely over the surface before desorbing. In the case of a radical, such as methoxy, if it encounters another radical, it will react before desorbing. Thus, the desorption of methoxy is very unlikely. The most favourable case is that it is the only radical at the surface, that is, with a concentration of $10^{-6}$. Methanol being typically about a few $10^{-2}$ (compared to water, the major compound of ice mantles), this leaves a radical such as methoxy with a ratio of less than $10^{-4}$ to methanol. If radicals come from the solid phase, it is not via the thermal desorption of physisorbed species.

In cold ($\le 100$ K) environments, thermal desorption is inefficient and processes such as chemical or reactive desorption could be at play (e.g. Garrod et al. 2007; Dulieu et al. 2013). The underlying idea is that excess energy released by a reaction on the surface of the grains may induce the desorption of a formed species. The efficiency of this process primarily depends on the energy released by the reaction and the binding energy, but also on the size of the molecules and the type of surface (e.g. Minissale & Dulieu 2014; Minissale et al. 2016b). Vasyunin et al. (2017) calculated that the probability of chemical desorption of $CH_3O$ formed on the grain surfaces by the hydrogenation of formaldehyde is $\sim 3.4 \times 10^{-5}$. This is low mainly because the $H_2CO + H \rightarrow CH_3O$ reaction is weakly exothermic. Such low probability would inject methoxy into the gas phase in a much smaller fraction than the destruction of methanol in the gas phase, as discussed in Section 4.3. Besides, we note that the hydrogenation of $H_2CO$ leads to $CH_3O$ and $CH_3OH$, so that the methanol to methoxy abundance ratio from chemical desorption would be even lower (see references and discussion in e.g. Minissale et al. 2016c). On top of this, hydrogenation of methanol would rather lead to $CH_3OH (+ H_2)$, because the hydrogen abstraction is easier in the methyl group. In practice, an additional fraction of $CH_3OH$ released from grains upon chemical desorption can also come from the chemical loop of abstraction-addition of methanol, without $CH_3O$ being an intermediary. Therefore, the estimate of the methanol and methoxy abundance ratio from chemical desorption is a lower limit, even if great caution must be taken when calculating the efficiency of chemical desorption.
Other ways to release molecules from ice mantles include photodesorption and cosmic ray bombardment. Bertin et al. (2016) estimated that the photodesorption of methanol (mixed or not with CO) can release a maximum of $3 \times 10^{-6}$ of CH$_3$O or CH$_2$OH per UV photon. Butscher et al. (2015) performed IR spectroscopy of photoinduced radicals as well as quantum calculations and demonstrated that CH$_2$OH is the favoured isomer to form CH$_3$OH. As direct evidence, no IR spectra of CH$_3$O could be taken whereas CH$_2$OH was observed. Therefore, the presence of CH$_3$O in ice is doubtful, and its photodesorption is likely very inefficient. Finally, cosmic rays may contribute to the release of molecules from ice mantles. Contrarily to UV photons, a large part of the molecular mantle desorbs as intact molecules even if over time it will evolve under radiolysis. In particular, Dartois et al. (2019) measured no extra CH$_3$O desorbing from fast heavy ions bombardment of methanol ice. Therefore, even in the case that cosmic rays release stable iCOMs in dark clouds, they do not enrich the gas phase with radicals.

In summary, the presence of CH$_3$O in the solid phase is still debated, and therefore, it is unlikely that the processing (by UV or ions) and sputtering of ices would inject gaseous methoxy in cold and dark regions.

## APPENDIX B: QUANTUM CHEMISTRY CALCULATIONS ON THE REACTION CH$_3$O + H

After old kinetics experiments at room temperature (and above) performed by Hoyermann et al. (1981) and Dobe et al. (1991), Xu et al. (2007) performed electronic structure calculations of the relevant PES for the CH$_3$O + H reaction and Variational Transition State Theory (VTST) coupled with RRKM calculations of its rate coefficient and BR as a function of the temperature and pressure. These calculations provided a deep insight into the CH$_3$O + H reaction and established that there are two possible reaction mechanisms (1) in one case, the CH$_3$O and H radicals recombine to form a methanol molecule in a barrierless addition reaction – the excess energy with which methanol is formed causes its fragmentation mainly into CH$_3$ + OH and, to a much lesser extent, into singlet CH$_2$ + H$_2$O; (2) in the other case, the approaching H atom is able to abstract one of the H atoms of CH$_3$O and to form, in a direct process, H$_2$CO and H$_2$, as experimentally suggested by Hoyermann et al. (1981) and Dobe et al. (1991) – the abstraction process is characterized by a small energy barrier of 0.11 kcal mol$^{-1}$ (0.46 kJ mol$^{-1}$, ~ 55 K) at the level of calculations employed (CCSD(T)/6-311+G(3df,2p)//CCSD/6-311+G(3df,2p)).

In light of the Xu et al. (2007) theoretical calculations, it is plausible that the reaction rate coefficients and product BR derived by Hoyermann et al. (1981) and Dobe et al. (1991) are referring mostly to the H-abstraction channel since the 0.11 kcal mol$^{-1}$ energy barrier can easily be surmounted in room temperature experiments. However, in the cold environments of ISM objects, the presence of an energy barrier could prevent the direct H-abstraction reaction from playing a significant role, and the barrier-less recombination pathway is expected to largely dominate. Two aspects, however, stopped us from using the rate coefficients calculated by Xu et al. (2007) as such, namely: (1) the lowest value of temperature considered is 200 K (that is much higher than the typical one, ~10 K, of the objects where methoxy has been observed) and (2) the rate coefficients were derived for a total pressure between 1 torr and 1000 atm using He as a bath gas. To fill in this gap, we have performed dedicated electronic structure and kinetics calculations to provide rate coefficients and product BRs for the CH$_3$O + H reaction in the 10–200 K range and

at the zero pressure limit, that is, the conditions that best simulate the environments where methoxy has been searched for.

In this section, we present the electronic structure and kinetics calculations for the reaction CH$_3$O + H that we have performed ad hoc to estimate the rate coefficients at the temperatures of interest for astrochemical models. Here, we briefly describe the used methods and the main results.

### B1 Computational methods

#### B1.1 Electronic structure calculations

Electronic structure calculations have been carried out mostly using the Gaussian 16 suite of programmes (Frisch et al. (2016)). Geometry optimizations have been performed followed by harmonic vibrational calculations, to verify that all reactants, intermediates, and products were true minima on PES and that all transition states (TSs) exhibited a single imaginary frequency. This step was done using the B2PLYP double-hybrid functional (Grimme (2006)) in conjunction with the aug-cc-pVTZ triple-ζ basis set (Kendall, Dunning Jr & Harrison (1992); (Woon & Dunning Jr 1993)). Semi-empirical dispersion effects were also included using the D3BJ model of Grimme Grimme, Ehrlich & Goerigk (2011); Goerigk & Grimme (2011)) leading to the so-called B2PLYP-D3 functional. cit The B3LYP hybrid functional (Becke (1993)), in conjunction with the 6-31G(d,p) basis set (Petersson et al. (1988); Petersson & Al-Laham (1991)), was instead used for the variational treatment of the barrierless addition channel. Electronic energies were then re-evaluated by performing single points with the CCSD(T) coupled-cluster method (Purvis III & Bartlett (1982); Pople, Head-Gordon & Raghavachari (1987)), in conjunction with the same aug-cc-pVTZ basis set. In this work, this level of theory will be designated as CCSD(T)/aug-cc-pVTZ//B2PLYP-D3/aug-cc-pVTZ.

With the approach described earlier, however, we failed to characterize the H-abstraction channel TS identified by (Xu et al. 2007). Therefore, we have used a different method to characterize the interaction between the approaching H atom and one of the H atoms of the methoxy radical. Indeed, it was possible to localize the stationary points only by performing a geometry optimization at CCSD(T) level with a basis set of aug-cc-pVTZ quality. This search was helped and corroborated by resorting to a semi-empirical model developed by (Pirani et al. 2008) already used in other cases (see, for instance, Marchione et al. (2022) and (Giani et al. 2023)).

#### B1.2 Kinetics calculations

Kinetic calculations for the addition reaction have been performed by using a code described in previous papers (Skouteris et al. 2015; Vazart et al. 2015; Skouteris et al. 2017, 2018, 2019; Vazart et al. 2022; Giani et al. 2023). The initial bimolecular rate coefficient leading from the reactants to the intermediate, together with the barrierless decomposition process, are evaluated using variational transition state theory (VTST). In the case of the back-dissociation to the reactants, a detailed balance argument has been used, where the unimolecular rate constant is given by the equation

$$k_{\text{back}} = k_{\text{capt}} \frac{\rho(R)}{\rho(I)},$$

where $k_{\text{capt}}$ is the capture rate constant, $\rho(R)$ is the density of states per unit volume for the reactants and $\rho(I)$ is the density of states for the intermediate.
For the steps where a well-defined TS exists, RRKM calculations have been performed. The microcanonical rate constant $k(E)$ is calculated using the formula

$$k(E) = \frac{N(E)}{h\rho(E)},$$

where $N(E)$ denotes the sum of states in the TS at energy $E$, $\rho(E)$ is the reactant density of states at energy $E$ and $h$ is Planck's constant. $N(E)$ is obtained by integrating the relevant density of states up to energy $E$ and the rigid rotor/harmonic oscillator model is assumed. Both densities of states (reactant and TS) are appropriately symmetrized with respect to the number of identical configurations of the reactants and/or TS. Tunnelling and quantum reflection have been taken into account by computing the tunnelling probability for an Eckart barrier having the same negative second derivative at the maximum of the pertinent saddle point. After all calculations have been performed, the BR between products and back-dissociation is determined for each energy and the corresponding capture rate constant is multiplied by this ratio to give the rate constant for the formation of the products. Finally, the rate constants are Boltzmann-averaged in order to provide the rate constants as a function of temperature.

Please note that, given the peculiar characteristics of the abstraction reaction, we have been unable to calculate its rate coefficient. We can only suggest that it is close to the gas kinetics limit (see further).

## B2 Results

### B2.1 Electronic structure: addition mechanism

Fig. B1 shows the stationary points along the PES starting with the addition of atomic hydrogen to methoxy. In line with the calculations by Xu et al. (2007), this is a barrierless process. A new covalent bond between the approaching H atom and the O atom of methoxy is formed and methanol is obtained as an addition intermediate. Starting from CH₃OH, which is ca. 430 kJ mol⁻¹ more stable than the reactants, three different exit channels are possible (i) dissociation into CH₃ and OH without an exit barrier; (ii) H-migration via the TS1 transition state, leading to the LM complex (see Fig. B1) that then dissociates into singlet CH₂ and H₂O; (iii) molecular H₂-elimination accompanied by the formation of either c-HCOH/t-HCOH or H₂CO.

The energy content of all intermediates and TSs are below the initial energy contents of the reactants, which makes all pathways viable in the ISM conditions.

### B2.2 Electronic structure: abstraction mechanism

Alternatively, a direct H-abstraction mechanism (to form H₂CO and H₂) is possible. Xu et al. (2007) located a small entrance barrier for this mechanism. According to our calculations, the interaction between CH₃O and the reacting H-atom pointing to one of the three equivalent H atoms of the methyl group leads first to the formation of a complex that is more stable than the reactants (optimized at the same level of calculation) only by 1.02 kJ mol⁻¹. This small stabilization is in line with the structure of the complex (see Fig. B2). The distance between the two interacting hydrogen atoms is very long (2.403 Å), as we can see from the optimized geometry of the complex. Once formed, the complex evolves towards the products H₂CO + H₂ through a TS whose geometry has been optimized at the same level of calculations (CCSD(T)/aug-cc-pVTZ). The TS is above the adduct by only 2.03 kJ mol⁻¹, that is, it is above the reactants by only 1.01 kJ mol⁻¹ at this level of calculations. From the figure, we can

notice that the distance between the reacting hydrogen and the H atom of the methyl group of the methoxy radical is much shorter in the TS with respect to the complex (1.746 with respect to 2.403 Å). At the same time, the C-H bond distance is increasing (from 1.106 to 1.122 Å).

We have not included the ZPE correction because the energy differences among the reactants, the complex, and the TS are very small (smaller than the predicted error in the calculations). The error in the zero-point energy calculations should be slightly larger for the reactants since they are considered not as a supermolecule but as separate fragments and for the TS with respect to the complex due to the imaginary frequency. For these reasons, we suppose that the ZPE corrections should be comparable for the reactants, complex, and TS. Since the energy differences are very small we check the possible presence of basis set superposition errors (BSSEs). To do that we have performed CCSD(T) calculations at the CCSD(T)/aug-cc-pVTZ optimized geometries, using both the cc-pVnZ and aug-cc-pVnZ with n = Q, 5, 6 basis sets.

The results are reported in Table B1. We can notice that the energy differences do not change significantly, suggesting that the BSSE error is very small also with the aug-cc-pVTZ basis set. At the best level of calculations employed, the complex is 0.9 kJ mol⁻¹ more stable than the reactants, while the TS is 0.8 kJ mol⁻¹ above the complex. Therefore, the TS is submerged by only 0.1 kJ mol⁻¹. Given the uncertainty of these calculations, we performed also an indirect phenomenological validity test of the complex stability and geometry. More specifically, the electronic polarizability $\alpha$ of the methoxy radical, which ca. 3.0 Å³, has been estimated by summing the CH₃ group (2.2 Å³) with the O atom (0.76 Å³) components (Cernuto et al. 2017) and the obtained value is consistent with data reported in the literature (NIST Computational Chemistry Comparison Benchmark Data base).

The $\alpha$ result suggests that, if the long-range interaction in CH₃O-H system is of van der Waals type, that is, determined by the combination of size repulsion with dispersion attraction and then exclusively dependent on $\alpha$ (Pirani et al. 2006), it must be consistent with the one operating in homologous systems involving Kr and Xe. Indeed, the methoxy polarizability value is intermediate between that of Kr ($\alpha$ = 2.49 Å³) and that of Xe ($\alpha$ = 4.04 Å³). The stabilization energies, determined exclusively by van der Waals forces, amount to -0.57 and -0.66 kJ mol⁻¹ for H-Kr and H-Xe, respectively. Related equilibrium distances are 3.57 and 3.82 Å for H-Kr and H-Xe, respectively. Since stability and centre-of-mass distance of the present complex are consistent with the interaction features obtained for H-Kr and H-Xe, this comparison not only confirms the reliability of calculated interaction energies given above but also suggests that the present complex exhibits a pure van der Waals nature.

### B2.3 Rate coefficients for the CH₃O + H → CH₃ + OH

We calculated the rate coefficients of the addition reaction for all possible channels described in Section B2.1 and shown in Fig. B3, from 1 to 300 K. As Xu et al. (2007) in the 200-300 K range, we found that the dominant channel is the one leading to CH₃ + OH, the other ones being negligible, even at the lowest temperatures considered in our study. In Fig. B1 the previous estimates of the rate coefficients are also shown, including those by Xu et al. (2007) (with the extrapolated values at temperatures lower than 200 K) and those adopted in previous astrochemical models (e.g. Ruaud et al. (2015)).

In the range of temperatures between 200-300 K, our rate coefficients are slightly lower than those derived by Xu et al. (2007). This
<figure><img src="image_5.png" /><figcaption>**Figure B1.** Simplified scheme of the PES of the reaction between methoxy and atomic hydrogen. The energies reported in blue (in kJ mol⁻¹) were obtained at the CCSD(T)/aug-cc-pVTZ//B2PLYP-D3/aug-cc-pVTZ level of theory and include the zero-point energy corrections.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure B2.** Optimized geometries of the initial adduct and the TS for the reaction CH₃O + H → CH₂O + H₂ at the CCSD(T)/aug-cc-pVTZ level. Bond distances are in Å and angles in °.</figcaption></figure>

<table><caption>**Table B1.** Relative energies (kJ mol⁻¹) computed at CCSD(T)/aug-cc-pVnZ and CCSD(T)/cc-pVnZ levels with n=T, Q, 5, 6 of the initial adduct and the TS for the reaction CH₃O + H → CH₂O + H₂ with respect to the reactants.</caption><thead><tr><th></th><th>Adduct</th><th>TS</th><th>Barrier height</th></tr></thead><tbody><tr><td>aug-cc-pVTZ</td><td>-1.02</td><td>+1.01</td><td>2.03</td></tr><tr><td>cc-pVQZ</td><td>-0.64</td><td>+0.37</td><td>1.01</td></tr><tr><td>cc-pV5Z</td><td>-0.78</td><td>+0.09</td><td>0.87</td></tr><tr><td>cc-pV6Z</td><td>-0.86</td><td>-0.05</td><td>0.81</td></tr><tr><td>aug-cc-pVQZ</td><td>-0.94</td><td>-0.09</td><td>0.85</td></tr><tr><td>aug-cc-pV5Z</td><td>-0.91</td><td>-0.11</td><td>0.80</td></tr><tr><td>aug-ccpV6Z</td><td>-0.92</td><td>-0.12</td><td>0.80</td></tr></tbody></table>

<table><caption>**Table B2.** The parameters α and β of the function describing the new computed rate coefficients of the reaction CH₃O + H → CH₃ + OH as a function of temperature. The curve has been split into two temperature intervals to obtain a good fit (see text).</caption><thead><tr><th>Temperature (K)</th><th>α (cm³ s⁻¹)</th><th>β</th></tr></thead><tbody><tr><td>1–80</td><td>3.25 × 10⁻¹⁰</td><td>0.21</td></tr><tr><td>80–300</td><td>1.70 × 10⁻¹⁰</td><td>-0.29</td></tr></tbody></table>
<figure><img src="image_7.png" /><figcaption>**Figure B3.** Rate coefficients of the reaction CH₃O + H as a function of temperature as computed by the present work (black curve), Xu et al. (2007) (red curves: the continuum curve refers to the measured values while the dashed curve refers to the extrapolation to temperatures lower than 200 K; see text) and as reported in the Ruaud et al. (2015) model and the KIDA data base (blue curves: the continuum curve refers to the CH₃ + OH channel while the dashed one to the H₂CO + H₂ channel).</figcaption></figure>

might be due to either small differences in the obtained energetics of the stationary points of the PES and/or to the fact that Xu et al. (2007) performed their calculations for a global pressure of 1 torr while we have considered the zero pressure limit, which better resembles the very low ISM number density. Interestingly, the extrapolation to low (≤20 K) temperatures of the trend with the temperature derived by Xu et al. (2007) provides rate coefficients that are almost one order of magnitude different than the values calculated here: for example, at 10 K, our rate coefficient is 1.7 × 10⁻¹⁰ cm³ s⁻¹ while the extrapolated value is 3.7 × 10⁻¹¹ cm³ s⁻¹.

In Fig. B3, we also show the values reported in the KIDA data base<fnref n="2" />, and currently used in many astrochemical models (e.g. Ruaud et al. (2015)). In KIDA, H₂CO + H₂ is reported to be the major channel of the CH₃O + H reaction, with a constant rate coefficient equal to 3 × 10⁻¹¹: at 10 K this is about five times lower than our value and more than ten times lower at ≥100 K. The channel producing CH₃ + OH is also assumed to be constant and its rate coefficient is taken to be a factor ten lower than that of the H₂CO + H₂ channel.

Finally, our computed rate coefficients, shown in Fig. B3, have been fitted with the function used in astrochemical networks, namely $k = \alpha(\frac{T}{300K})^\beta \exp[-\gamma/T]$, where we used $\gamma = 0$. To have a good fit, we fitted the curve in two separate temperature intervals, 1–80

and 80–300 K, respectively. The derived values of $\alpha$, $\beta$, and $\gamma$ are reported in Table B2.

B2.4 The rate for the CH₃O + H → H₂CO + H₂

Given the peculiarity of this reaction pathway in the entrance channel, a theoretical calculation of the rate coefficient challenges the methods available in the literature. However, we can speculate that, even though it is an open channel under the ISM conditions, direct H-abstraction is not expected to be competitive with the barrierless addition mechanism at low temperatures because of the attractive nature of the deep potential well associated with methanol formation. Indeed, the presence of a van der Waals complex and a submerged TS for the H-abstraction mechanism will reduce the reactive flux along this pathway because of the recrossing back towards the reactants. Therefore, it is reasonable to expect this channel to be minor under the low T conditions of relevance to the cases considered here. In all cases, we have tested in the astrochemical models different values for the BR of this reaction, but no significant changes have been observed. This is because there are many other more important formation routes for all the species involved, that is H₂CO, CH₃, and OH.

This paper has been typeset from a TEX/LATEX file prepared by the author.

<footnote n="2">²http://kida.astrophy.u-bordeaux.fr/; Wakelam et al. (2015).</footnote>