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

HAL
open science

# Chemical constraints on the dynamical evolution of the cold core L694

A Taillard, V Wakelam, P Gratier, E Dartois, M Chabot, J. A. Noble, L Chu

► To cite this version:

A Taillard, V Wakelam, P Gratier, E Dartois, M Chabot, et al. Chemical constraints on the dynamical evolution of the cold core L694. Astronomy & Astrophysics - A&A, 2025, 698, pp.A278. 10.1051/0004-6361/202553796. hal-05153672

HAL Id: hal-05153672

https://hal.science/hal-05153672v1

Submitted on 9 Jul 2025

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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Chemical constraints on the dynamical evolution of the cold core L694*

A. Taillard<fnref n="1" />,**, V. Wakelam<fnref n="2" />, P. Gratier<fnref n="2" />, E. Dartois<fnref n="3" />, M. Chabot<fnref n="4" />, J. A. Noble<fnref n="5" />, and L. Chu<fnref n="6" />

Received 17 January 2025 / Accepted 24 April 2025

## ABSTRACT

In star-forming regions, molecular cloud history and dynamics set the trend in the chemical composition. Ice formation, in particular, is affected by the evolution of physical conditions, which can lead to different ice compositions within the same cloud. In cold cores with medium densities >10⁴ cm⁻<fnref n="3" />, low temperatures <15 K, and low UV radiation <G₀, most complex organic molecules are formed on dust grain surfaces and are released back into the gas phase through non-thermal mechanisms such as sputtering or heating by cosmic-rays, photodesorption, or chemical desorption. Studying both gas- and solid-phases can help observers to add constraints on the chemical and dynamical evolution of cold cores. We present a study of the cold core L694, observed with the IRAM 30m single-dish radio telescope. Observed species include CO (and its isotopologues) and CH₃OH, a key chemical species precursor of more complex organic molecules. We applied an inverted non-local thermal equilibrium radiative transfer code, previously used on observations of the pre-stellar core L429-C, in order to obtain gas-phase abundances by deriving the column densities of the detected species from the spectroscopic parameters of the targeted molecular transitions (intensity, line width), and from physical parameters derived from archival observations (H₂ volume density and gas temperature). This allowed us to probe the molecular abundances as a function of density and visual extinction. In parallel, we ran chemical models (both static and dynamic) to constrain the evolution time of the core by directly comparing the observations with the model outputs. We then compared the compositions of the cold cores L429-C and L694. The gas-phase abundances in L694 all exhibit a common depletion profile (with high variability in the depletion factor), as the core is identified to be in a more advanced (infalling) state compared to L429-C. The physical parameters of the two cores are, however, very similar, leading to close evolutionary timescales in our static models. The dynamical model fails to reproduce the CO gas-phase abundances at high density, predicting an evolutionary timescale that is too short compared to static models. A more detailed study on the parameter constraining the CO freeze-out could help to better constrain the timescale.

**Key words.** astrochemistry – ISM: abundances – ISM: clouds – ISM: molecules – ISM: individual objects: L694

## 1. Introduction

Complex organic molecules (COMs, composed of six atoms or more, as defined by Herbst & van Dishoeck 2009) begin to form in the early phases of star formation. These steps are crucial for setting the trend in the chemical composition for the protostellar environment. Cold cores are dependent on the physical conditions of their original molecular cloud. They are characterised by densities higher than 10³ cm⁻<fnref n="3" />, temperatures lower than 20 K, and visual extinctions higher than 3 mag (Bergin & Tafalla 2007). Since the temperature is low, the chemistry in the gas phase is limited, while a rich chemistry takes place at the surface of interstellar dust grains. Grains act as catalysts and play the role of the 'meeting point' for atoms, radicals, and other molecules: species at the surface diffuse and react with other molecules to form more complex products (Cuppen et al. 2017). In transluc-ent and dense clouds with Av > 1.5, the grains are covered by

a mantle where H₂O ice is the dominant component (Boogert et al. 2015). A shift occurs at densities >10⁵ cm⁻<fnref n="3" />, where volatile CO efficiently freezes out onto the grains (Bacmann et al. 2002; Qasim et al. 2018). The CO freeze-out leads to the formation of COMs on short timescales of a few tens of thousands of years (Cuppen et al. 2009), with the build-up of a CO-dominated layer over the H₂O-dominated mantle. Infrared observational studies have revealed cold core regions presenting high column densities of COMs, such as CH₃OH (Boogert et al. 2011; Chu et al. 2020; Perotti et al. 2020; McClure et al. 2023).

These regions are of particular interest, since they show a large variability in molecular composition despite sharing sim-ilar physical conditions. For instance, this variety can be found within a single source (e.g. in ices, Noble et al. 2017). This is the case in the Taurus Molecular Cloud (TMC), located at a dis-tance that varies between ~ 128 and 198 pc depending on the part of the cloud being studied (Onishi et al. 2002; Galli et al. 2019). Within the cloud, the widely targeted filament TMC-1 (Agúndez & Wakelam 2013, and reference therein) exhibits dif-ferent chemical behaviour at different locations. Over a wider group of targets, Fuente et al. (2019) showed that for the same

<footnote n="1">¹ Centro de Astrobiología (CAB), CSIC-INTA, Carretera de Ajalvir km4, 28850 Torrejón de Ardoz, Madrid, Spain</footnote>

<footnote n="2">² Laboratoire d'Astrophysique de Bordeaux (LAB), Univ. Bordeaux, CNRS, B18N, allée Geoffroy Saint-Hilaire, 33615 Pessac, France</footnote>

<footnote n="3">³ Institut des sciences Moléculaires d'Orsay, CNRS, Université Paris-Saclay, Bât 520, Rue André Rivière, 91405 Orsay, France</footnote>

<footnote n="4">⁴ Université Paris-Saclay, CNRS/IN2P3, IJCLab, 91405 Orsay, France</footnote>

<footnote n="5">⁵ Physique des Interactions Ioniques et Moléculaires, CNRS, Aix Marseille Univ., 13397 Marseille, France</footnote>

<footnote n="6">⁶ NSF NOIRLab, 950 N. Cherry Ave., Tucson, AZ 85719, USA</footnote>

<footnote>* This work is based on observations carried out under project number ID 111-21 with the 30m telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain).</footnote>

<footnote>** Corresponding author: angele.taillard@cab.inta-csic.es</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Left: H₂ column density map (N<sub>H</sub>₂ in cm⁻²) of L694, computed from Herschel maps at 250, 350, and 500 µm. The position of the continuum peak is represented by the blue cross. The red cross is at an off-position, used to show the different velocity profile in Fig. 2. Right: H₂ volume density (n<sub>H</sub>₂ in cm⁻³) obtained from the method presented in Marchal et al. (2019). In both maps, the temperature (in K) appears in white contours, and the spatial scale is represented by the horizontal white bar (10 000 AU).</figcaption></figure>

density of ~10⁵ cm⁻³ in various clouds, the CS gas-phase abundance varies by over more than one order of magnitude. Studying such differences through chemical modelling, which employs similar physical parameter inputs, can help convey the most probable hypothesis for each cloud. Yet, the chemistry within a model can be completely affected by a single parameter that may not be well constrained.

In this work, we studied the chemical composition of the cold core L694 using gas-phase observations obtained with the IRAM 30m telescope and through chemical modelling. By comparing observed and modelled abundances, we aimed to put constraints on the formation timescale of cold cores. This paper follows a methodology used in a previous study of the chemical composition of the cold core L429-C (Taillard et al. 2023). The source, the gas-phase observations, and the observational results are described in Sect. 2. In Sect. 3, we present and discuss the molecular abundances derived from observations. The chemical models – both static and dynamic – as well as the comparison of their results with the observations are given in Sect. 4. Lastly, in Sect. 5, we compare our results from these observations in L694 with those from our previous study on the cold core L429-C and discuss the time estimation derived from our chemical models.

**2. Source description and new observations**

**2.1. The LDN 694 cold core**

LDN 694 (sometimes referred to as L694-2, but hereafter L694) is a cold (T < 18 K), dense (n<sub>H</sub>₂ > 1 × 10⁶ cm⁻³), and isolated core in the Aquila Rift, located at a distance of approximately 205 pc (Kim et al. 2022). The core has been labelled as infalling, as molecular line profiles exhibit redshifted self-absorption – strong evidence of inwards motions (Lee et al. 2001) – as well as an indication of a velocity gradient, possibly due to rotation (Lee et al. 2004) and faster infall motions with supersonic speeds (Lee et al. 2007; Seo et al. 2013). Another study of the N₂H⁺ emission shows a centrally peaked dense core, with a line width increasing towards the centre, possibly due to infall (Crapsi et al. 2005). The N₂H⁺, CS, and DCO⁺ spectra shown in this above studies present blue-shifted components on the source position. The core is centrally condensed, with an estimated mass of 1 M⊙

(Lee et al. 2001), and is expected to form a protostellar object within 10⁴ years (Visser et al. 2002; Williams et al. 2006). No radio point source has been detected so far (Harvey et al. 2002), and the core presents a steep density profile that drops quickly (over 0.1 pc) as L694 exhibits an elongated, filament-like structure (Chu et al. 2020; Kim et al. 2022) (see Fig. 1). In Chu et al. (2020), the authors deduced a visual extinction up to ~27 mag (from the JHK photometric bands) towards the reddest of three detected stars, located at less than 0.1 pc from the centre of the core.

The core presents an advanced chemistry: it has a large column density of ortho-N₂D⁺ (Caselli et al. 2008, 3.2 × 10¹³ cm⁻²) and a high deuterium fractionation (N(N₂H⁺)/N(N₂D⁺) ~ 0.26) (Crapsi et al. 2005). Studies of C³⁴S by Kim et al. (2022) show that the signal is not centrally peaked. Moreover, the emission distribution suggests high depletion in the high-density region. The core also exhibits a CO depletion (with a depletion factor of 11 to 16) that is larger towards the nucleus (Wiström et al. 2016). A particularity of L694 is that the ice observations made by Chu et al. (2020) using the NASA infrared telescope facility (IRTF) show strong infrared absorption at a distance of 0.1 pc from the peak density towards at least one background star near the edge of the core. The authors measured a CH₃OH/H₂O ratio of 14.7% (larger than in L429-C for similar observations, Boogert et al. 2011). Moreover, only a small amount of CO is frozen onto the dust grain (≤15% with respect to CO<sub>gas</sub>), but ~30% of the CO ice is mixed with methanol – the highest such ratio measured in their core sample. This high ratio indicates that the mixture traces the densest parts of the cloud, where the CO depletion is the highest and the conversion from CO to CH₃OH occurs at a faster rate. L694 also presents a remarkably high CH₃OH/ice/CO<sub>tot</sub> ratio (0.064) with respect to the entire sample, which is a mixture of background stars and YSO observations within Barnard 59, LDN 483, LDN 694, and Barnard 335. In a follow-up study of five molecular cores (Chu & Hodapp 2021), the authors measured the total hydrogen volume densities to constrain ice formation processes within the same sample and derived Av maps using NIR photometry of a large population of background stars (i.e. L694 data used below in Appendix E). This study shows that the density required to form CH₃OH is only reached in less than 2% of the total core in the densest region.
<table><caption>**Table 1.** Detected lines and associated spectroscopic information.</caption><thead><tr><th>Molecule</th><th>Frequency (MHz)</th><th>Transition</th><th>E<sub>up</sub> (K)</th><th>g<sub>up</sub></th><th>A<sub>ij</sub> (s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>N<sub>2</sub>H<sup>+</sup></td><td>93173.9</td><td>(1–0)</td><td>4.47</td><td>15</td><td>3.62 × 10<sup>-5</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>96739.3</td><td>(2–1)</td><td>12.5</td><td>5</td><td>2.55 × 10<sup>-6</sup></td></tr><tr><td>CH<sub>3</sub>OH</td><td>96741.3</td><td>(2,0)–(1,0)</td><td>7</td><td>5</td><td>3.40 × 10<sup>-6</sup></td></tr><tr><td><sup>34</sup>SO</td><td>97715.3</td><td>(2,3)–(1,2)</td><td>9.1</td><td>7</td><td>1.07 × 10<sup>-5</sup></td></tr><tr><td>CS</td><td>97980.9</td><td>(2–1)</td><td>7.1</td><td>5</td><td>1.67 × 10<sup>-5</sup></td></tr><tr><td>SO</td><td>109252.2</td><td>(3,2)–(2–1)</td><td>21.1</td><td>5</td><td>1.08 × 10<sup>-5</sup></td></tr><tr><td>C<sup>18</sup>O</td><td>109782.1</td><td>(1–0)</td><td>5.27</td><td>3</td><td>6.26 × 10<sup>-8</sup></td></tr><tr><td>HNCO</td><td>109905.7</td><td>(5,0,5)–(4,0,4)</td><td>15.82</td><td>11</td><td>1.80 × 10<sup>-8</sup></td></tr><tr><td>NH<sub>2</sub>D</td><td>109782.1</td><td>(1–0)</td><td>5.27</td><td>3</td><td>6.26 × 10<sup>-8</sup></td></tr><tr><td><sup>13</sup>CO</td><td>110201.3</td><td>(1–0)</td><td>5.29</td><td>3</td><td>6.29 × 10<sup>-8</sup></td></tr><tr><td>C<sup>17</sup>O</td><td>112358.7</td><td>(1–0)</td><td>5.39</td><td>3</td><td>6.69 × 10<sup>-8</sup></td></tr></tbody></table>

## 2.2. New IRAM 30 m observations

During November 2021, the core was observed for 23 hours with IRAM 30m to provide short-spacing data. We centred the map at RA = 19h41m11s.63, Dec = -10°56′27.7″ (J2000). Our observation setup focused on four frequency bands corresponding to: 93–94.8 GHz, 96.2–98.1 GHz, 108.6–110.4 GHz, and 111.9–113.7 GHz. With a beam of ~28.5″, we obtained ~300″ × 300″ maps. Our previous study on a different cold core, which included interferometer data (NOEMA), showed no detectable small scale emission in L429-C (Taillard et al. 2023). We did not expect the beam dilution and the sub-beam source structure to have a significant effect on our results for L694. Each cube contained 70 velocity channels of 0.15 km s⁻<fnref n="1" />. On average, we observed noise sensitivity varying between 0.15 to 0.30 K depending on the molecule. The data reduction was performed using the Gildas package CLASS<fnref n="1" />. We used the Cologne Database for Molecular Spectroscopy<fnref n="2" /> (CDMS, Müller et al. 2001) and the Jet Propulsion Laboratory catalogue<fnref n="3" /> to conduct line identification. Moreover, we detected all the molecules (and their lines) listed in Table 1. As the setup differed slightly from our previous study, some molecules (i.e., N₂H⁺, NH₂D, and HNCO) were newly observed, while others (CCS, HC₃N, CN, and H₂S) were not. The molecules targeted give us constraints on the chemical evolution timescale (mainly CO and its isotopologues – CS and CH₃OH) and the dynamical state of the cloud (N₂H⁺, CO, and its isotopologues).

### 2.2.1. Detected lines and intensity maps

In the frequency band targeted, we detect a total of 11 lines (24, if we considered the hyperfine structure of N₂H⁺) coming from ten different molecules. We considered a line detected when its intensity was three times the mean noise level in the data cube. The following molecules we observed did not fulfil this condition: CCS (6–7) at 93 870 MHz, S¹⁸O (1–2) at 93 267 MHz, OCS (8–7) at 97301 MHz, and C³⁴S (2–1) at 96 412 MHz. The integrated intensity maps were built from the data cubes using the Spectral Cube python package (Ginsburg et al. 2019). Molecular emission can be divided into two categories: (1) emission widely observed across the filament (C¹⁸O, C³⁴SO, CH₃OH, and CS), and (2) emission centred on the continuum peak (HNCO,

<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Upper: C¹⁸O and CS spectra on the core position represented by the dark blue cross on the maps. Lower: spectra of the same molecules at the position marked by the red cross in Fig. 1. The grey dashed line is the centroid velocity.</figcaption></figure>

NH₂D, N₂H⁺, and SO). Integrated intensity maps can be found in Appendix A (Fig. A.1), with the position of the continuum peak marked by a blue cross. In the following subsection, we use the cardinal point convention with respect to right ascension to describe the positions on the map.

### 2.2.2. Kinematic analysis and infall

Most molecules present simple emission profiles that can be fitted by a single Gaussian. The channel velocity maps can be found in Appendix B (Figs. B.1). Figure 2 shows the C¹⁸O and CS spectra at two positions on the map: the dust continuum peak emission and the off-position in the eastern portion of the map (blue and red crosses, respectively, in Fig. 1). C¹⁸O presents a

<footnote n="1">¹ http://www.iram.fr/IRAMFR/GILDAS</footnote>

<footnote n="2">² https://cdms.astro.uni-koeln.de/classic/</footnote>

<footnote n="3">³ https://spec.jpl.nasa.gov</footnote>
Gaussian profile over most of the cloud, with stronger intensity along the filament. Two distinct regions exhibit a deformed Gaussian profile: on the lower border of the cloud, the red wing is slightly stronger than the blue one, while around the continuum peak, the Gaussian profile flattens, and the blue component becomes stronger (see Fig. 2). As discussed in the next section, this coincides with the depletion profile observed in abundance. The CS spectrum at the continuum peak exhibits a bimodal profile with a strong component (between 8.9 and 9.6 km s⁻¹) and a weaker component (between 9.6 and 10.4 km s⁻¹). The bimodal profile is evident throughout the entire filament, with the second component weakening until a high intensity peak is reached in the eastern part of the cloud (red cross in Fig. 1), where the CS profile becomes a single Gaussian. This velocity structure is observed in infalling cores such as L1544 (Tafalla et al. 1998). The authors hypothesized that this signature is due to inward motions towards the core. It was also observed towards a similar region, the TMC, by Lee et al. (2001). L694 was previously labelled as a 'strong infall candidate' in Lee et al. (2004), where a strong blue component and a less intense red component in the (3–2) CS transition at 146.96 GHz were observed at the position of the source. In both studies, the authors observed that N₂H⁺ emission is typically found in a significantly more compact region of the cloud compared to the more diffuse CS emission. In the region traced by N₂H⁺, CS exhibits a double-peaked profile, meaning that these two profiles show inward motions. In our case, the N₂H⁺ intensity is strongest at the continuum peak position (see Fig. A.1). CS, moreover, shows a double-peaked profile. All these elements combined confirm the infalling nature of L694.

## 2.3. Physical structure of the source and data at hand

We used two methods to investigate the physical parameters of the source based on observations, called method one and method two. The data used in method one to derive the physical structure of L694 were taken from the Herschel Science Archive<fnref n="4" /> (OBSIDS: 1342230846). L694 was observed with SPIRE at 250, 350, and 500 µm. Different studies estimated the core temperature with values ranging from 8 to 14 K (Evans II et al. 2001; Harvey et al. 2003; Williams et al. 2006; Seo et al. 2013; Chu & Hodapp 2021; Lin et al. 2023). Using the same pipeline as described in Taillard et al. (2023), we derived the dust temperature, the H₂ column density, and the H₂ volume density across the entire region (Fig. 1) observed with the IRAM 30m telescope. We obtained a dust temperature ranging from 9 K at the dust continuum peak to ~15 K at the core edges. In Fig. C.1 (left), we plot the volume density as a function of the Av obtained from the Herschel data.

Method one follows the standard approach of balancing radiative heating and cooling of grains, as used in numerous studies (Stutz et al. 2010; Launhardt et al. 2013; Lippok et al. 2016; Sadavoy et al. 2018). However, it only probes a certain grain size population, whereas larger grains are known to be cooler and smaller ones warmer (Herbst & Cuppen 2006). For simplicity, it is often assumed that grains have the same size, typically adopting the canonical grain radius of 0.1 µm, which represents the mass-weighted average over the grain size distributions (Mathis et al. 1977; Kruegel 2003)). Dust emission can be fitted by spectral energy distribution (SED), and one can deduce a temperature approximation from a modified blackbody model.

We used a second method ('method two') to check if the temperature derived with method one accurately traced the core temperature. We describe the process in Appendix C and the resulting abundances in Appendix E, using background star photometric measurements and the Av determination method explained in Chu & Hodapp (2021) for L694. We then converted the Av to N<sub>H₂</sub>. Based on the method described in Hocuk et al. (2017), we estimated T<sub>dust</sub> using the H₂ column density. We computed the temperature, Av, and H₂ column density, shown in the central panel of Fig. C.1. The temperature obtained is, overall, 3 K lower over the density grid, ranging from 7 to 12 K. Compared to method one, Av increases less steadily but covers the same range of values (2–30 mag). The density is slightly higher by a factor of 2 towards the edge of the cloud but is the same at the continuum dust peak position. The main issue with method two is that it does not probe the deepest region of the cloud. Therefore, we rely on the kernel density estimation to smooth the values around the dust continuum peak (see Chu & Hodapp 2021, for more details). Overall, this method gives us almost the same range of values for Av and the H₂ column density as in method one, but gives a slightly lower temperature. The physical conditions derived from these two methods are used to study the chemistry discussed in Sect. 3 and Appendix E.

# 3. Chemical composition of L694

## 3.1. Molecular abundance maps

In order to obtain molecular abundance (X<sub>mol</sub> = N<sub>mol</sub>/N<sub>H₂</sub>) maps from our gas-phase observations, we used the method described in Taillard et al. (2023), employing Herschel maps of the regions to derive the physical parameters (i.e. method one). The isotopic ratio we used for the conversion between C<sup>18</sup>O and CO is 557 (Wilson 1999). Each molecule was treated using ndRADEX (Taniguchi 2020), with collision rates extracted from the following original studies: NH<sub>2</sub>D from Müller et al. (2001), N<sub>2</sub>H<sup>+</sup> from Schöier et al. (2005), SO from Lique et al. (2005), CS from Lique et al. (2006), CO from Yang et al. (2010), CH<sub>3</sub>OH from Rabli & Flower (2010), and HNCO from Sahnoun et al. (2018).

The maps obtained using method one are presented in Figs. 3 and 4. As can be seen from the decreasing abundances around the blue cross, all the molecules present a depletion profile at the continuum peak.

The CO shows a strong depletion profile at the continuum peak, with its lowest abundance of 3.9 × 10⁻⁵ and its maximum of 1.9 × 10⁻<fnref n="4" />. The depletion between the two density extremes observed here is approximately a factor of 7. We can compute its depletion factor, f<sub>CO</sub>, by comparison to its canonical abundance measured in molecular clouds. The depletion factor is computed as f<sub>CO</sub> = f(X<sub>can</sub>/X<sub>12CO</sub>), where X<sub>can</sub> = 8.5 × 10⁻⁵ is the canonical abundance of <sup>12</sup>CO measured by Frering et al. (1982) in Taurus and ρ Oph. We obtain f<sub>CO</sub> = 2.23. This value is lower than that obtained by Crapsi et al. (2005), who used the same canonical abundance to find f<sub>CO</sub> = 11 at the core position. The difference between the results mainly arises from the two methods adopted: Crapsi et al. (2005) used Eq. (1) from Caselli et al. (2002b), which depends on the integrated intensity of C<sup>18</sup>O and the observed brightness at 1.2 mm. They note that the depletion factor is subject to multiple uncertainties, mainly due to the determination of the canonical abundance and the amount of the molecule actually probed in the densest region. CS exhibits optically thick emission lines from our RADEX computation, meaning that the abundance we determined is an upper-limit on the actual quantity. CS likewise exhibits a significant depletion

<footnote n="4">⁴ https://archives.esac.esa.int/</footnote>
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Gas-phase abundances of observed molecules with respect to H₂ on a logarithm scale (log(X<sub>mol</sub>)). The green cross marked on the CH₃OH map indicates the position of the background star towards which methanol ice was observed by Chu et al. (2020). The blue star represents the continuum peak position.</figcaption></figure>

gradient along the filament, where the abundance decreases to
$1.6 \times 10^{-10}$ at the continuum peak. The maximum value, located
west of the source, is $1.7 \times 10^{-7}$: three orders of magnitude
higher than the minimum value. For CH₃OH, the abundance
is highest east of the continuum peak (in a similar position to
SO and HNCO), with a value of $9.8 \times 10^{-9}$, and reaches a min-
imum of $5.2 \times 10^{-10}$ at the continuum peak. The observational
constraint on the ice abundance by Chu & Hodapp (2021) corre-
sponds to a volume density of $\sim 4.3 \times 10^4 \text{ cm}^{-3}$. Comparing this
methanol ice column density with our gas-phase one at the same

position, we find a gas-to-ice ratio of approximately 0.003%.
In L429-C, at a similar density of ~2.4 × 10⁴ cm⁻³, the ratio
was approximately 0.002%. N₂H⁺ exhibits hyperfine structure
with three distinct peaks. All three transitions trace the same
region, with critical densities at 12 K for the transitions span-
ning a narrow range from 2 × 10⁴ to 6.3 × 10⁵ cm⁻³. In this
case, we computed the column density by fitting only the mid-
dle component (with the rest frequency listed in Table 1), after
checking that all hyperfine components were optically thin and
traced the same spatial distribution. Maximum abundance can
<figure><img src="image_6.png" /><figcaption>Fig. 4. Maps similar to those presented in Fig. 3 for HNCO and NH₂D.</figcaption></figure>

be found at the filament border, with a value of 7.7 × 10⁻¹⁰, while minimum abundance, 5.5 × 10⁻¹¹, occurs at the continuum peak. As the literature (Caselli et al. 2002a; Lippok et al. 2013) finds, a correlation exists between the CO depletion profile and the presence of N₂H⁺ tracing the dense part of the core. N₂H⁺ is considered a late depleter (Bergin et al. 2002), and its low abundance at the peak position suggests an aged core. The late depletion of N-bearing species has been confirmed in recent studies (Redaelli et al. 2019; Caselli et al. 2022; Lin et al. 2020, 2023). In Pineda et al. (2022), the authors suggest that the depletion of another N-bearing species, NH₃, begins at densities higher than 2 × 10⁵ cm⁻³. In our case, N₂H⁺ depletion occurs at a lower density (~10⁴ cm⁻³). In Lin et al. (2020), the authors find an N₂H⁺ depletion of ~27 in the cold core L1512 at a central density of 10⁵ cm⁻³, while the depletion we observe in L694 is by a factor of ~15 at a density of the same order of magnitude. SO is detected in an even smaller area than N₂H⁺, showing maximum abundance to the east of the source, with a value of 4.8 × 10⁻⁸. Its minimum abundance is found at the dust continuum peak with a value of 9.7 × 10⁻¹⁰. HNCO and NH₂D are detected only in a small area immediately surrounding the core. These two molecules exhibit a tighter range of values between their extrema. HNCO abundance at the continuum peak is 1.4 × 10⁻¹⁰, with a maximum found east of this position – similar to SO – at 3.2 × 10⁻¹⁰. For NH₂D, the maximum abundance is 5.6 × 10⁻¹⁰, measured next to the continuum peak position, where the abundance is 1.5 × 10⁻¹⁰.

<table><caption>Table 2. Maximum abundance and depletion factor.</caption><thead><tr><th>Molecule</th><th>Maximum abundance</th><th>Depletion factor<sup>(*)</sup></th></tr></thead><tbody><tr><td>CO</td><td>1.9 × 10<sup>-4</sup></td><td>7</td></tr><tr><td>CH<sub>3</sub>OH</td><td>9.8 × 10<sup>-9</sup></td><td>19</td></tr><tr><td>CS</td><td>1.7 × 10<sup>-7</sup></td><td>1062</td></tr><tr><td>SO</td><td>4.8 × 10<sup>-8</sup></td><td>50</td></tr><tr><td>HNCO</td><td>1.4 × 10<sup>-10</sup></td><td>2</td></tr><tr><td>N<sub>2</sub>H<sup>+</sup></td><td>7.7 × 10<sup>-10</sup></td><td>14</td></tr><tr><td>N<sub>2</sub>D</td><td>5.6 × 10<sup>-10</sup></td><td>4</td></tr></tbody></table>

Notes. (*)Depletion factor is the maximum abundance over the lowest abundance found at the dust continuum peak.

## 3.2. Molecular abundance as a function of physical parameters

When comparing these molecular abundances (and their standard deviation) to the physical parameters in Figs. 5 and 6, all molecular abundances decrease with both volume density and visual extinction. We summarise maximum abundances and depletion factors (i.e. the maximum abundance over the abundance at the continuum peak position) in Table 2. CO shows a depletion gradient of more than a factor of 5 between volume densities of 2 × 10³ cm⁻³ and 2 × 10⁶ cm⁻³. Its abundance at a lower density seems to plateau where the abundance saturates to a value close to 3.1 × 10⁻⁴. CS has the strongest depletion profile along the filament, with an abundance varying from 1.7 × 10⁻⁷ to ~1.6 × 10⁻¹⁰, a factor >1000. CH₃OH shows a strong depletion gradient varying from 9.8 × 10⁻⁹ to 5.2 × 10⁻¹⁰ (by a factor of ~20) at densities increasing from 3 × 10³ to 10⁶ cm⁻³. SO shows a similar profile to CH₃OH, with strong depletion between 10⁴ and 10⁶ cm⁻³, as its abundance drops from 4.8 × 10⁻⁸ to 9.7 × 10⁻¹⁰, corresponding to a decrease by a factor of approximately 50. Despite their detection in relatively small regions of the map, HNCO, N₂H⁺, and NH₂D exhibit a depletion profile with factors of ~2, 14 and 4, respectively.

## 4. Comparison with chemical models

To compare with our observational results, we ran two sets of simulations using the Nautilus gas-grain model (see Sect. 4.1). As we explain in Sect. 4.2, we first adopted the fixed physical conditions observed in our region, as used in the radiative transfer analysis. Next, (see Sect. 4.3), we used time-dependent physical parameters as computed by a hydrodynamical model to follow the formation of the cores. For both sets of simulations, we compared the model results to the observed abundances over the density range where all the molecules were detected (from 3 × 10⁴ to 10⁶ cm⁻³).

### 4.1. The Nautilus gas-grain model

To run both sets of simulations, we used the gas-grain chemical code Nautilus (Ruaud et al. 2016; Wakelam et al. 2024). The code is available through a github repository<fnref n="5" />. Nautilus is a three-phase model that uses the rate-equation approach to compute, at each time step, the abundance of hundreds of chemical species in the solid (i.e. mantle and surface) and gas phase. The gas-phase reaction rates were computed using the latest

<footnote n="5">⁵ https://astrochem-tools.org/codes/</footnote>
<figure><img src="image_7.png" /><figcaption>**Fig. 5.** Left: abundance as a function of hydrogen volume density $n_{H_2}$ for CO, CH$_3$OH, CS, and SO. Right: abundance as a function of visual extinction, Av. Colour coding indicates the temperature. The mean abundance and standard deviation for each bin appear in green.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Fig. 6. Same as Fig. 5, for HNCO, N₂H⁺, and NH₂D.</figcaption></figure>

kida.uva.2024 chemical network (Wakelam et al. 2024). The surface reaction network was based on the rate equation method (as described in Hasegawa et al. 1992; Ruaud et al. 2016; Wakelam et al. 2024). The mobilities of the species at the surface of the grains and in the mantle were set as a fraction of the binding energies (0.8 in the mantle and 0.4 at the surface). As a result, diffusion was more efficient on the surface than in the mantle. As photodissociation (both by direct and secondary UV photons) can occur both at the surface and in the mantle. Mechanisms such as whole grain cosmic-ray heating, cosmic-ray sputtering, photo-desorption, and chemical desorption were included in addition to thermal desorption. All details on the processes included in the model are given in Ruaud et al. (2016) and Wakelam et al. (2024).

In our static and dynamical simulations, we used a set of standard parameters described below. We used the initial chemical composition listed in Table 1 of Wakelam et al. (2021b), except for the sulphur elemental abundance, which was set to

8 × 10⁻⁸ to better reproduce the S-bearing molecules, (following the low-metal abundances from Graedel et al. 1982). For the cosmic-ray ionisation rate (ζ), we used the prescription by Wakelam et al. (2021b) (i.e. Eqs. (6) and (7)), which depends on visual extinction. As in the observational analysis and as described in Sect. 2.3, we assumed that the gas temperature is equal to that determined using method one. As in previous studies (Taillard et al. 2023; Clément et al. 2023), we did not use this temperature for the dust, since these maps trace the outer region of the cloud where smaller and warmer grains reach temperatures higher by a few Kelvin. A slightly higher grain temperature can strongly affect chemical abundances, especially in the ices. We therefore adopted the parametric expression for dust temperature derived by Hocuk et al. (2017), in which the grain temperature varies as a function of Av. Following Clément et al. (2023), we added 1 K to the dust temperature obtained from this method.
<table><caption>**Table 3.** Physical parameters of the static models.</caption><thead><tr><th>Model</th><th>log(Density) (cm<sup>-3</sup>)</th><th>T<sub>gas</sub> (K)</th><th>Av (mag)</th></tr></thead><tbody><tr><td>1 & 9<sup>(*)</sup></td><td>4.5</td><td>11.3</td><td>8.3</td></tr><tr><td>2 & 10</td><td>4.7</td><td>10.8</td><td>11.1</td></tr><tr><td>3 & 11</td><td>4.9</td><td>10.4</td><td>14.8</td></tr><tr><td>4 & 12</td><td>5.1</td><td>10.0</td><td>19.1</td></tr><tr><td>5 & 13</td><td>5.3</td><td>9.8</td><td>23.3</td></tr><tr><td>6 & 14</td><td>5.6</td><td>9.6</td><td>26.7</td></tr><tr><td>7 & 15</td><td>5.8</td><td>9.5</td><td>28.7</td></tr><tr><td>8 & 16</td><td>6.0</td><td>9.4</td><td>29.8</td></tr></tbody></table>

**Notes.** (*) For each set of physical conditions (n(H₂), T, Av), we ran two models using either the H₂O-rich or the CO₂-rich CR sputtering yield, see text.

## 4.2. Static models

### 4.2.1. Physical parameters

Using the three-phase model described above, we first ran the model for a grid of physical parameters similar to those obtained with method one (and used to analyse the observed lines), considering only the physical regions where CH₃OH is detected.

We ran two sets of eight models. Each of the eight were sampled on the density range (i.e. from 3 × 10⁴ to 10⁶ cm⁻³) previously computed with method one. Visual extinction and gas temperature were computed from these different density samples (T<sub>gas</sub> decreasing from 11.3 to 9.4 K, and Av increasing from 8.3 to 29.8 mag). These physical parameters were fixed during the entire integration time. The sputtering yield was measured in experiments as a function of the predominant molecule in ices, showing that it decreases in H₂O-rich ices, while CO₂-rich ices exhibit the opposite behaviour (Dartois et al. 2018, 2020). This parameter is particularly important for the non-thermal desorption of certain molecules, especially methanol (Wakelam et al. 2021a; Taillard et al. 2023). For each set of eight models, we used either the H₂O-rich or CO₂-rich sputtering yield. In total, there are 16 unique models, which we summarise in Table 3. We then ran the models for each of the physical conditions, starting from the initial conditions described in Sect. 4.1, over a period of 10⁷ yr.

### 4.2.2. Results

To select the time that best reproduces our observations for each physical condition, we used the distance of disagreement – as described in Wakelam et al. (2006) – applied to three observed molecules (CS, CO, and CH₃OH), which are detected over the largest range of density in the cloud. The distance of disagreement defines the ‘best time’ as the chemical integration time at which the modelled abundances are the closest to the ones observed. The observed abundances used for comparison with the models are the average values within each density bin, where the observed density was binned according to the density parameters adopted in the models (see Table 3). The best time as a function of density is shown in Fig. 7, where we observe a general trend across all the sets of models: the time required to reproduce the observations decreases with density, from 1.3 × 10⁵ yr at 3.1 × 10⁴ cm⁻³ to 1.0 × 10⁴ yr at 1.0 × 10⁶ cm⁻³. The two sets of models yield similar results. As in Taillard et al. (2023), where we applied the same method to another core, the constraint is mainly given by CO, whose abundance was primarily set by

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

**Fig. 7.** Best time from the static models, obtained using the distance of disagreement based on three observed molecules (CO, CS, and CH₃OH), as a function of density (logarithmic scale). Brown circles correspond to L694, and blue circles to L429-C, as obtained from our previous study (Taillard et al. 2023).

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

**Fig. 8.** Methanol gas-to-ice ratio obtained for the two different sets of models as a function of density (in logarithm scale). The low and high CR sputtering yield appear in blue and orange, respectively. In black, we plot the IRTF observation from Chu et al. (2020) and its associated error.

the timescale required for depletion onto the grains (see also discussion Section 3 in Wakelam et al. 2021b), rather than by the cosmic-ray sputtering yield, which only affects abundance after depletion has occurred. In Appendix D, we discuss the level of agreement for the main molecules at the best times.

In Fig. 8, we compare the gas-to-ice ratio measured in L694 by Chu et al. (2020) with our best models. We find that the model with the low sputtering yield (shown in blue) produces a ratio closer to the observed value (0.008% at 5.1 × 10⁴ cm⁻³, as indicated by the black dot), while the models with the higher sputtering yield (shown in orange) differ from the observations by a factor of 10. At the density corresponding to the IRTF observations, we also compared our model results with the CO ice abundance (not shown in the figure). We find that all the models reproduce the observed CO ice abundance exactly, at approximately 1.2 × 10⁻⁴.

In Appendix E, we present the results of the static chemical models using the physical parameters derived with method two. We also compare these results with the observed gas-phase abundances.
### 4.3. Dynamical models

#### 4.3.1. Physical parameters

In previous studies, we showed that the timescale evolution between the diffuse medium and the dense core impacts the computed chemical composition relative to the static models described above (Ruaud et al. 2018; Wakelam et al. 2019; Clément et al. 2023). As suggested by the comparison between the static chemical models and the observations (see also Taillard et al. 2023), the densest parts of the cold core may experience a faster evolution of their density than the external parts. We therefore compared our observations to chemical simulations in which the physical parameters evolve with time, rather than remain fixed. The time dependent physical parameters were computed with the Smooth Particle Hydrodynamical (SPH) model presented in Bonnell et al. (2013). These simulations follow the gas dynamics in a galactic potential including spiral arms. In a $250 \times 250$ pc$^2$ region of the simulation, we identified 12 cold cores that formed at the end of the simulation (i.e. with a density larger than $10^5$ cm$^{-3}$). We selected all the particles in a region of 0.5 pc in a radius around the peak density and reconstructed the history of each of these particles in terms of density, gas temperature, and visual extinction. The number of particles per core varies between 80 and 350. A full description of these simulations can be found in Ruaud et al. (2018). These time-dependent physical conditions were then used as inputs for the Nautilus gas-grain model. The other parameters (i.e. $\zeta$, elemental abundances, and computation of dust temperature) are as described in Section 4.1. In the next section, we present one core as an example; the conclusions are consistent across all simulated cores.

#### 4.3.2. Model results

To compare the model results with the observations, we used the output of the chemical composition at the time of core formation. As previously shown (Ruaud et al. 2018; Clément et al. 2023), the various histories of the different particles produce a spread in the chemical abundances within each core. As such, for a sampling of the final densities, we computed the mean molecular abundances and standard deviation for each core. In Fig. 9, we show the model results for one selected core as a function of density and superimpose the observed abundances in L694 (this work) and L429-C (Taillard et al. 2023). We discuss the comparison between these two cold cores in the next section. In the simulations, the abundance of molecules (i.e. CO, CS, CH$_3$OH, HNCO, and N$_2$H$^+$) increase with density before decreasing. For CO, the models reproduce the observations up to a density of a $10^4$ cm$^{-3}$, but fail at high density. In fact, the modelled CO gas-phase abundances, computed with the H$_2$O-rich cosmic-ray sputtering yield, are lower by more than one order of magnitude at $10^5$ cm$^{-3}$ and above. Using the higher cosmic-ray sputtering yield appropriate for CO$_2$-rich ice, the predicted CO gas-phase abundance at high density increases, but remains significantly lower than that observed. For CS, the models (for both values of the sputtering yield) reproduce the decrease in abundance at densities above a $10^4$ cm$^{-3}$, but fail to reproduce the abundances at lower densities, where the modelled values are much lower than the observations. Both sets of models nevertheless reproduce N$_2$H$^+$ well. The constraint on the CR sputtering yield is particularly strong for HNCO, which is also well reproduced by the models using the H$_2$O-rich CR sputtering yield. However, since both HNCO and N$_2$H$^+$ are only observed at high density, the model predictions at low density cannot be tested. Methanol

is the only molecule for which the agreement with the dynamical simulations is poor across the parameter space. Further discussion on the static and dynamical model results is provided in Sect. 5.2.

## 5. Discussion

### 5.1. Comparing the L429-C and L694 cold cores

In our previous study (Taillard et al. 2023), we applied the same methodology to L429-C, a cold core located in the Aquila Rift, in a less advanced stage (not yet infalling). Comparing the physical parameters for both sources using method one, the visual extinction towards the central part of the core is approximately a factor of two higher in L429-C than in L694, with $\text{Av} \sim 60$ in L429-C and $\sim 30$ in L694. In L694, we also probe regions with lower densities ($10^3$ cm$^{-3}$) compared to L429-C ($5 \times 10^3$ cm$^{-3}$). The minimum temperature is similar in the two sources, reaching down to $\sim 9$ K in L694 and to $\sim 10$ K in L429-C. The observed abundances as a function of density in the two sources can be seen in Fig. 9. The CO gas-phase abundance is less depleted in L694 than in L429-C at both low ($5 \times 10^3$ cm$^{-3}$) and high density ($10^6$ cm$^{-3}$), and is similar between these densities. The CH$_3$OH abundance is higher in L694 for densities lower than $3.1 \times 10^5$ cm$^{-3}$, while it is similar in the two sources above. Unlike in L429-C, the methanol gas-phase abundance is not flat but decreases with density. The CS gas-phase abundance is similar in both sources. HNCO and N$_2$H$^+$ were observed in L429-C. We compare the models' best times in Fig. 7, which shows the results for the two cores: L429-C (in blue) and L694 (in red). For both pre-stellar cores, the high and low sputtering yield models overlap. The difference between the two sources is minor. Nevertheless, it appears that the timescale of evolution is a bit shorter in L429-C at low density. The external parts of this source might be chemically younger than those of L694. This would be consistent with the fact that this source is at a less evolved stage. We note that, despite the CO abundance difference observed at high density (by a factor of $\sim 2$), the model is unable to find a different evolutionary time.

### 5.2. Using molecular abundances to trace the dynamical evolution of the cores

With the static chemical models, we find that the densest parts of L694 (up to densities of $10^6$ cm$^{-3}$) should take approximately $10^4$ years to form from the diffuse medium to reproduce the observed molecular abundances. This constraint was obtained by comparing the observed and modelled molecular abundances, under the assumption that the gas is initially atomic. In the dynamical simulations, the timescale of the density increase is computed by the SPH simulations and thus depends on the physics and the kinematics of the regions. If we consider the CO gas-phase abundance, the dynamical model strongly underproduces the observed values at this density (Fig. 9). In these simulations, CO forms early once the density exceeds a few $10^3$ cm$^{-3}$. Figure 10 shows the results of one trajectory from the dynamical simulations of the selected cloud. In this example, the final drop in the CO gas-phase abundance occurs when the density becomes higher than a few $10^4$ cm$^{-3}$. It takes $5 \times 10^5$ years to go from $10^4$ to $2 \times 10^5$ cm$^{-3}$. This timescale is much longer than the 'best' time obtained from the static simulations. Carbon monoxide has a robust chemistry, and its gas-phase abundance is a good indicator of the depletion processes at play (see also the discussion in Sect. 3 in Wakelam et al. 2021b). Our results
<figure><img src="image_11.png" /><figcaption>**Fig. 9.** Mean abundances computed with the dynamical model (filled circles) as a function of density (see text). Standard deviation for the cold core is shown with vertical bars. Red and blue stars mark the observed abundances in L694 (this study) and L429-C (Taillard et al. 2023), respectively. Crosses and diamonds denote the abundances obtained at the best times with the static models, which only probe the density where the three molecules (CO, CS, and CH₃OH) used to determine the 'best time' are detected. Results obtained from using both cosmic-ray sputtering yields are shown for all the models.</figcaption></figure>

show that the time-dependent density evolution in the dynamical simulations we used is incompatible with the observed high abundance of CO, assuming that CO sticking is correctly reproduced in our model. Three parameters can impact CO depletion: density, time, and the sticking coefficient. Since we have observational constraints on density, time and the sticking coefficient

represent good avenues for model improvement. Apart from these physical parameters, both the chemistry and physics of CO can also directly impact its gas-phase abundance. For example, the efficiency of CO conversion into other species (e.g. COMs, Fuchs et al. 2009; Chuang et al. 2016; Simons et al. 2020; Molpeceres et al. 2024) will modify its solid-phase abundance.
<figure><img src="image_12.png" /><figcaption>**Fig. 10.** Abundance of CO in the gas phase (black points) and density (blue points) as a function of time computed for one core in the dynamical simulations (see text).</figcaption></figure>

Furthermore, its subsequent desorption via non-thermal mechanisms will depend on the efficiency of those processes (Mennella et al. 2004; Muñoz Caro et al. 2010; Fayolle et al. 2011; Dartois et al. 2019; Kruczkiewicz et al. 2024).

## 6. Conclusions

We observed the cold core L694 with the IRAM 30m telescope and detected ten molecules via 11 transitions (24 when accounting for the hyperfine structure of $N_2H^+$), including CH$_3$OH and CO isotopologues. We then derived molecular abundance maps (300'' × 300'') for all the molecules, allowing us to study how abundances depend on the physical conditions (primarily density) and to compare with a previous study of L429-C (Taillard et al. 2023). These observations were compared to predictions from the Nautilus chemical model to constrain the time evolution of the cloud. Our main results are:

- With the exception of CS (2–1) which exhibits two velocity components, emission spectra obtained from millimetre observations show simple velocity structure. Considering that the $N_2H^+$ emission is strongest at the densest part of the core and that the C$^{18}$O emission spectra widened in this area, we confirm the infalling nature of L694;

- All the molecules show strong depletion at the dust continuum peak, indicating a chemically advanced and dynamically active region. We find a CO depletion factor $f_{CO} = f(X_{can}/X_{12CO})$ of 2.23 at the dust continuum peak;

- We ran both static and dynamical Nautilus models to compare with our observations. Using the static models, we demonstrated that the timescale required to reproduce our gas-phase observations decreases with density, with the density itself varying by more than a factor of ~30 across the full range explored. Overall, the abundances are reproduced within a factor of 10, which we find satisfactory;

- We observed a methanol gas-to-ice ratio of 0.003% at a density of ~4.3 × 10$^4$ cm$^{-3}$. We need more data points in order to recover the increasing trend with the density observed in L429-C in our previous study. This value is reproduced within a factor of 3 by the best model (with low sputtering yield);

- The physical parameters derived for the two cold cores studied (L429-C and L694) using Herschel data show little difference. The gas phase, however, exhibits a more

advanced dynamic in L694, where all the molecules show depletion profiles. This is consistent with our static models, which predicted the need for shorter times to reproduce the observations in L429-C;

- The dynamical SPH models exhibit a range of chemical abundances that depend on the histories of individual particles. These models fail to reproduce the abundances of most of the molecules at low density and gas-phase CO abundance at high density. The depletion trends and abundances of CS and CH$_3$OH at high densities, however, are well reproduced;

- Using molecular abundances to trace the dynamical evolution of the core, we show that our dynamical models under-predict the CO gas-phase abundance at high density by at least one order of magnitude. In fact, the slow evolution of the density in these simulations results in a high depletion of CO, which is not supported by our observations. A deeper study on the parameters impacting CO depletion (sticking coefficient, density, and time) could lead to improvements.

As this source is the target of a JWST GTO proposal (PID 1187), this study was conducted in preparation for the chemistry that will be revealed by the telescope. A complete ice map of the source would enable a direct comparison with the gas-phase maps we obtained and would be significantly beneficial for studying the entire ice and/or gas interface.

*Acknowledgements.* We thank our referee, Olli Sipilä, for their useful comments that made the manuscript better. AT acknowledges funding from the European Research Council under the European Union's Horizon 2022 research and innovation program (grant agreement No. 101096293 SUL4LIFE). AT, VW, PG, ED, MC, and JN acknowledge the thematic Action "Physique et Chimie du Milieu Interstellaire" (PCMI) of INSU Programme National "Astro", with contributions from CNRS Physique & CNRS Chimie, CEA, and CNES.

## References

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

Bacmann, A., Lefloch, B., Ceccarelli, C., et al. 2002, A&A, 389, L6

Bergin, E. A., & Tafalla, M. 2007, ARA&A, 45, 339

Bergin, E. A., Alves, J., Huard, T., & Lada, C. J. 2002, ApJ, 570, L101

Bonnell, I. A., Dobbs, C. L., & Smith, R. J. 2013, MNRAS, 430, 1790

Boogert, A. C. A., Huard, T. L., Cook, A. M., et al. 2011, ApJ, 729, 92

Boogert, A., Gerakines, P., & Whittet, D. 2015, Annu. Rev. Astron. Astrophys., 53, 541

Bron, E., Daudon, C., Pety, J., et al. 2018, A&A, 610, A12

Caselli, P., Benson, P. J., Myers, P. C., & Tafalla, M. 2002a, ApJ, 572, 238

Caselli, P., Walmsley, C. M., Zucconi, A., et al. 2002b, ApJ, 565, 344

Caselli, P., Vastel, C., Ceccarelli, C., et al. 2008, A&A, 492, 703

Caselli, P., Pineda, J. E., Sipilä, O., et al. 2022, ApJ, 929, 13

Chu, L. E. U., & Hodapp, K. W. 2021, ApJ, 918, 2

Chu, L. E. U., Hodapp, K., & Boogert, A. 2020, ApJ, 904, 86

Chuang, K. J., Fedoseev, G., Ioppolo, S., van Dishoeck, E. F., & Linnartz, H. 2016, MNRAS, 455, 1702

Clément, A., Taillard, A., Wakelam, V., et al. 2023, A&A, 675, A165

Crapsi, A., Caselli, P., Walmsley, C. M., et al. 2005, ApJ, 619, 379

Cuppen, H. M., van Dishoeck, E. F., Herbst, E., & Tielens, A. G. G. M. 2009, A&A, 508, 275

Cuppen, H. M., Walsh, C., Lamberts, T., et al. 2017, Space Sci. Rev., 212, 1

Dartois, E., Chabot, M., Barkach, T. I., et al. 2018, A&A, 618, A173

Dartois, E., Chabot, M., Barkach, T. I., et al. 2019, A&A, 627, A55

Dartois, E., Chabot, M., Bacmann, A., et al. 2020, A&A, 634, A103

Evans II, N. J., Rawlings, J. M. C., Shirley, Y. L., & Mundy, L. G. 2001, ApJ, 557, 193

Fayolle, E. C., Bertin, M., Romanzin, C., et al. 2011, ApJ, 739, L36

Frerking, M. A., Langer, W. D., & Wilson, R. W. 1982, ApJ, 262, 590

Fuchs, G. W., Cuppen, H. M., Ioppolo, S., et al. 2009, A&A, 505, 629

Fuente, A., Navarro, D. G., Caselli, P., et al. 2019, A&A, 624, A105

Galli, P. A. B., Loinard, L., Bouy, H., et al. 2019, A&A, 630, A137

Ginsburg, A., Koch, E., Robitaille, T., et al. 2019, Zenodo

Graedel, T. E., Langer, W. D., & Frerking, M. A. 1982, ApJS, 48, 321

Harvey, D. W. A., Wilner, D. J., Di Francesco, J., et al. 2002, AJ, 123, 3325

Harvey, D. W. A., Wilner, D. J., Myers, P. C., & Tafalla, M. 2003, ApJ, 597, 424
Hasegawa, T. I., Herbst, E., & Leung, C. M. 1992, ApJS, 82, 167

Herbst, E., & Cuppen, H. M. 2006, PNAS, 103, 12257

Herbst, E., & van Dishoeck, E. F. 2009, ARA&A, 47, 427

Hocuk, S., Szűcs, L., Caselli, P., et al. 2017, A&A, 604, A58

Kim, S., Lee, C. W., Tafalla, M., et al. 2022, ApJ, 940, 112

Kruczkiewicz, F., Dulieu, F., Ivlev, A. V., et al. 2024, A&A, 686, A236

Kruegel, E. 2003, Phys. Today, 57, 53

Launhardt, R., Stutz, A. M., Schmiedeke, A., et al. 2013, A&A, 551, A98

Lee, C. W., Myers, P. C., & Tafalla, M. 2001, ApJS, 136, 703

Lee, C. W., Myers, P. C., & Plume, R. 2004, ApJS, 153, 523

Lee, S. H., Park, Y.-S., Sohn, J., Lee, C. W., & Lee, H. M. 2007, ApJ, 660, 1326

Lin, S.-J., Pagani, L., Lai, S.-P., Lefèvre, C., & Lique, F. 2020, A&A, 635, A188

Lin, Y., Spezzano, S., Pineda, J. E., et al. 2023, A&A, 680, A43

Lippok, N., Launhardt, R., Semenov, D., et al. 2013, A&A, 560, A41

Lippok, N., Launhardt, R., Henning, T., et al. 2016, A&A, 592, A61

Lique, F., Spielfiedel, A., Dubernet, M.-L., & Feautrier, N. 2005, J. Chem. Phys., 123, 134316

Lique, F., Spielfiedel, A., & Cernicharo, J. 2006, A&A, 451, 1125

Marchal, A., Miville-Deschênes, M.-A., Orieux, F., et al. 2019, A&A, 626, A101

Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425

McClure, M. K., Rocha, W. R. M., Pontoppidan, K. M., et al. 2023, Nat. Astron., 7, 431

Mennella, V., Palumbo, M. E., & Baratta, G. A. 2004, ApJ, 615, 1073

Molpeceres, G., Furuya, K., & Aikawa, Y. 2024, A&A, 688, A150

Muñoz Caro, G. M., Jiménez-Escobar, A., Martín-Gago, J. Á., et al. 2010, A&A, 522, A108

Müller, H. S. P., Thorwirth, S., Roth, D. A., & Winnewisser, G. 2001, A&A, 370, L49

Noble, J. A., Fraser, H. J., Pontoppidan, K. M., & Craigon, A. M. 2017, MNRAS, 467, 4753

Olofsson, S., & Olofsson, G. 2010, A&A, 522, A84

Onishi, T., Mizuno, A., Kawamura, A., Tachihara, K., & Fukui, Y. 2002, ApJ, 575, 950

Perotti, G., Rocha, W. R. M., Jørgensen, J. K., et al. 2020, A&A, 643, A48

Pineda, J. E., Harju, J., Caselli, P., et al. 2022, AJ, 163, 294

Predehl, P., & Schmitt, J. H. M. M. 1995, A&A, 293, 889

Qasim, D., Chuang, K.-J., Fedoseev, G., et al. 2018, A&A, 612, A83

Rabli, D., & Flower, D. R. 2010, MNRAS, 406, 95

Redaelli, E., Bizzocchi, L., Caselli, P., et al. 2019, A&A, 629, A15

Ruaud, M., Wakelam, V., & Hersant, F. 2016, MNRAS, 459, 3756

Ruaud, M., Wakelam, V., Gratier, P., & Bonnell, I. A. 2018, A&A, 611, A96

Ryter, C. E. 1996, Ap&SS, 236, 285

Sadavoy, S. I., Keto, E., Bourke, T. L., et al. 2018, ApJ, 852, 102

Sahnoun, E., Wiesenfeld, L., Hammami, K., & Jaidane, N. 2018, J. Phys. Chem. A, 122, 3004

Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F., & Black, J. H. 2005, A&A, 432, 369

Seo, Y. M., Hong, S. S., & Shirley, Y. L. 2013, ApJ, 769, 50

Simons, M. A. J., Lamberts, T., & Cuppen, H. M. 2020, A&A, 634, A52

Stutz, A., Launhardt, R., Linz, H., et al. 2010, A&A, 518, L87

Tafalla, M., Mardones, D., Myers, P. C., et al. 1998, ApJ, 504, 900

Taillard, A., Wakelam, V., Gratier, P., et al. 2023, A&A, 670, A141

Taniguchi, A. 2020, https://doi.org/10.5281/zenodo.4139707

Visser, A. E., Richer, J. S., & Chandler, C. J. 2002, AJ, 124, 2756

Wakelam, V., Herbst, E., & Selsis, F. 2006, A&A, 451, 551

Wakelam, V., Ruaud, M., Gratier, P., & Bonnell, I. A. 2019, MNRAS, 486, 4198

Wakelam, V., Dartois, E., Chabot, M., et al. 2021a, A&A, 652, A63

Wakelam, V., Gratier, P., Ruaud, M., et al. 2021b, A&A, 647, A172

Wakelam, V., Gratier, P., Loison, J. C., et al. 2024, A&A, 689, A63

Williams, J. P., Lee, C. W., & Myers, P. C. 2006, ApJ, 636, 952

Wilson, T. L. 1999, Rep. Progr. Phys., 62, 143

Wiström, E. S., Charnley, S. B., Cordiner, M. A., & Ceccarelli, C. 2016, ApJ, 830, 102

Yang, B., Stancil, P. C., Balakrishnan, N., & Forrey, R. C. 2010, ApJ, 718, 1062
**Appendix A: Integrated intensity maps**

We show the integrated intensity maps in Fig. A.1 obtained for the brightest transitions of each molecule observed. These maps were obtained by integrating the velocity channels where signal is found to be > 3 times the mean rms for each pixel.

<figure><img src="image_13.png" /><figcaption>Fig. A.1: Integrated intensity maps for each detected molecule. In the case of multiple transitions detected for one molecule, the brightest one is shown. The continuum peak is shown with a blue cross.</figcaption></figure>
**Appendix B: Velocity channel maps**

We present the velocity channel maps obtained for some molecules (C¹⁸O, SO, CS and CH₃OH), underlining the dynamical structure of the cloud, in Figs B.1.

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

Fig. B.1: Velocity channel maps of C¹⁸O (109782 MHz) (top, left), SO (99299 MHz) (top, right), CS (97980 MHz) (bottom, left) and CH₃OH (96741 MHz) (bottom, right), for the velocity channels from 8.85 to 10.05 km s⁻¹. The telescope beam is shown in hatched on the lower side of each panel. Contours (outer to inner) are 0.2, 0.5, 1, 1.5, 2.5, 3.5 K (T*).
## Appendix C: L694 physical parameters

In this section we compare how the physical parameters vary as a function of one another (Fig. C.1) for both sources (L429-C and L694). In the case of L694, both the temperatures obtained with method 1 (top) and method 2 (middle) are shown (see Sect. 2.3). The method 2 data used to derive the physical parameters is detailed in Sect. E. The temperatures determined with method 2 ranges from ~ 9.3 K on the core to ~ 15.5 K in the outer part of the cloud for L694 (top) and from ~ 11.7 to ~ 18 K for L429-C (bottom). Their densities are high at the continuum peak ( > 1.5 × 10⁶ cm⁻³) and lowers on the outer part of the cloud (~10³ cm⁻³). N_H₂ was obtained from the optical depth from Herschel with the same method as described in Appendix A in Taillard et al. (2023) and Av were computed from the conversion factor: N_H = 1.8 × 10²¹ × Av. The Av varies similarly, from 25 at the continuum peak to < 1 on the outer part of the cloud for L694. For L429-C, the Av is greater by a factor > 2, varying from ~ 9 to ~ 77 at high densities, which is a result of a higher N_H₂.

<figure><img src="image_15.png" /><figcaption>Fig. C.1: Av as a function of n<sub>H₂</sub> (cm<sup>-3</sup>) and in color-scale, the dust temperature (K) for, top to bottom, L694 computed from method 1, method 2 (Chu & Hodapp 2021), L429-C Taillard et al. (2023).</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Fig. D.1: Ratio between modelled (X<sub>mod</sub>) and observed abundances (X<sub>obs</sub>) in the gas phase for CO, CS and CH<sub>3</sub>OH as a function of density for each best time, in all four sets of models.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. D.2: Ratio between modelled (X<sub>mod</sub>) and observed abundances (X<sub>obs</sub>) in the gas phase for HNCO, SO and N<sub>2</sub>H<sup>+</sup> as a function of density for each best time, in all sets of models.</figcaption></figure>

## Appendix D: Goodness of fit for modelled molecules

In the Figs. D.1 and D.2, we compare the quality of the goodness of fit by looking at the ratio between models and observations in the gas phase. Overall, all molecules are reproduced within a factor of 10 (except for HNCO, within a factor of 20) and there is no difference in chemical composition when using two different CR ionisation yields. CO and CS very are well re-

produced in all models. For CH<sub>3</sub>OH, the models with a low sputtering are underproducing (by a factor 5) the molecule as density increases while in the high sputtering ones, the low density abundances are overproduced (by a factor > 6). It is however the models with low sputtering that are the closest to the observed methanol gas-to-ice ratio (See Sect. 4.2.2). SO is over-produced in all models, from a factor 3 to almost 10 at high density. N<sub>2</sub>H<sup>+</sup> is over-reproduced by a factor ~ 2 on the entire range of densities where it is detected. HNCO is the molecule with the worst X<sub>mod</sub>/X<sub>obs</sub> ratio, being over-produced by a factor ~ 30 on all the density range. In Wirström et al. (2016), the authors measured a O<sub>2</sub> abundance in L694 of < 1.6 × 10<sup>-7</sup> on the maximum of the dust continuum. With Nautilus, we predict an abundance of ~ 5 × 10<sup>-7</sup>, which is over-produced by a factor of 3.
## Appendix E: Deriving abundance maps using background stars observations (“method two”)

We used the background star data provided by Chu & Hodapp (2021) to construct physical parameters maps with a different method other than method 1. In this study, the authors have determined the extinction from background stars photometry (see section 4.1 of Chu & Hodapp 2021), and obtained a value of Av for each position of stars they observed. Using a similar smoothing method as they did, we managed to derive an almost identical extinction map as theirs (Fig. 1 in their study). The Av determined on the dust continuum peak has bigger uncertainties, as it is too dense to observe any background stars and we conveyed values by smoothing everything with a Gaussian kernel. The authors of the study also noted that in the region of highest extinction, the resolution under-samples the data, meaning that the Av determined in this area (shown as black pixels in their Figure 1) is interpolated from the three closest stars at most. Once we obtained the Av map, we derived the H₂ column density map by using the following equation (Predehl & Schmitt 1995; Ryter 1996; Olofsson & Olofsson 2010): $N_{H_2} = 1.8 \times 10^{21}/2$. The H₂ column density ranges from ~2 ×10²¹ to 2.6 ×10²² cm⁻². From $N_{H_2}$, we derive the H₂ volume density $n_{H_2}$ by using the same method as our previous work (derived from Bron et al. 2018). This gives us a volume density ranging from ~1.9 × 10³ to ~1.7 × 10⁷ cm⁻³. Lastly, to obtain the dust temperature, we used the same formula as described in section 4.1 from the parametrization derived in Hocuk et al. (2017). The temperature obtained range from ~7 K to ~12.5 K, which is the coldest of the temperatures presented in Figure C.1, plotted as the colour-map of the density as a function of Av. With these physical parameters, we computed the molecular column densities and then abundances for all the detected molecules. We found that the derived molecular abundances are similar (within a factor of 3) to the ones we obtained with method 1. We show in Figure E.1 the different abundances of CO, CS and CH₃OH obtained with the two methods.

We also ran the static chemical models described in section 4.2 using these new set of physical parameters and compared with the observed abundances derived using these parameters. The best time also varies linearly and are slightly smaller than the one we derived using method 1 ranging from 1 ×10⁶ to 7.7 ×10³ yr for all set of parameters. The goodness of fit for gas-phase species is slightly worse than with method 1 but all the predictions are still fitting within a factor 10. In Figure. E.2 we show the goodness of fit for CO, CS and CH₃OH. CO is reproduced within a factor 3 no matter the set of models or density. CS is reproduced within a factor of two at density < 10⁵ cm⁻³ and the agreement worsens as density increases. SO at low density is 10 times under-produced compared to observations but the ratio $X_{\text{mod}}/X_{\text{obs}}$ increasing with density and is overproduced at higher density within a factor 3. Methanol is under-produced (up to a factor 10) in low sputtering yield models at densities higher than 10⁵ cm⁻³. In high sputtering yield model, the fit is bad at densities lower than 10⁵ cm⁻³ and then is well reproduced at higher densities (within a factor 3). $N_2H^+$ is reproduced the same than what is shown in Fig. D.2. HNCO is better reproduced within a factor 10, favoured by the cold chemistry in these physical parameters as compared to the models presented in this study. CO ices and methanol ices are predicted the same and the CH₃OH gas-to-ice ratio is identical, both in behaviour and values in all sets of models. Overall, a few species are strongly impacted by low temperature chemistry (CS, SO) while others (CH₃OH, CO,

<figure><img src="image_18.png" /><figcaption>Fig. E.1: Comparison between the gas-phase abundances for C¹⁸O, CS and CH₃OH between method 1 (in orange) and method 2 (in blue).</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. E.2: Ratio between modelled ($X_{\text{mod}}$) and observed abundances ($X_{\text{obs}}$) in the gas phase for CO, CS and CH₃OH as a function of density for each best time, in all four sets of models for method 2.</figcaption></figure>

HNCO) have their behaviour changed compared to method 1 but are still well reproduced.