Context. IRAS 16293E is a rare case of a prestellar core being subjected to the effects of at least one outflow. Aims. We want to disentangle the actual structure of the core from the outflow impact and evaluate the evolutionary stage of the core. Methods. Prestellar cores being cold and depleted the best tracers of their central regions are the two isotopologues of the trihydrogen cation that are observable from the ground orthoHSUB2SUBDSUPSUP and paraDSUB2SUBHSUPSUP. We used the Atacama Pathfinder EXperiment APEX telescope to map the paraDSUB2SUBHSUPSUP emission in IRAS 16293E and collected James Clerk Maxwell Telescope JCMT archival data of orthoHSUB2SUBDSUPSUP. We compared their emission to that of other tracers including dust emission and analysed their abundance with the help of a 1D radiative transfer tool. The ratio of the abundances of orthoHSUB2SUBDSUPSUP to paraDSUB2SUBHSUPSUP can be used to estimate the stage of the chemical evolution of the core. Results. We have obtained the first complete map of paraDSUB2SUBHSUPSUP emission in a prestellar core. We compare it to a map of orthoHSUB2SUBDSUPSUP and show their partial anticorrelation. This reveals a strongly evolved core with a paraDSUB2SUBHSUPSUPorthoHSUB2SUBDSUPSUP abundance ratio towards the centre for which we obtain a conservative lower limit from 3.9 at 12 K to 8.3 at 8 K while the high extinction of the core is indicative of a central temperature below 10 K. This ratio is higher than predicted by the known chemical models found in the literature. ParaDSUB2SUBHSUPSUP and indirectly orthoHSUB2SUBDSUPSUP is the only species that reveals the true centre of this core while the emission of other molecular tracers and dust are biased by the temperature structure that results from the impact of the outflow. Conclusions. This study is an invitation to reconsider the analysis of previous observations of this source and possibly questions the validity of the deuteration chemical models or of the reaction and inelastic collisional rate coefficients of the HSUPSUPSUB3SUB ltinlineformula idFI1gtltalternativesgttexmath idtexeq1ltCDATAmathrmH3gttexmathmmlmath xmlnsmmlhttpwww.w3.org1998MathMathML displayinline idmmleq1mmlmsubsupmmlmrow classMJXTeXAtomORDmmlmi mathvariantnormalHmmlmimmlmrowmmlmn3mmlmnmmlmrow classMJXTeXAtomORDmmlmommlmommlmrowmmlmsubsupmmlmathltinlinegraphic xmlnsxlinkhttpwww.w3.org1999xlink idimgeq1 mimesubtypepng mimetypeimage xlinkhrefaa4735123eq1.pnggtltalternativesgtinlineformula isotopologue family. This could impact the deuteration clock predictions for all sources.
['astrochemistry', 'radiative transfer', 'ISM: abundances', 'ISM: molecules', 'ISM: structure', 'ISM: individual objects: IRAS 16293E', 'Astrophysics - Astrophysics of Galaxies']
First map of DSUB2SUBHSUPSUP emission revealing the true centre of a prestellar core Further insights into deuterium chemistry
{'First map of D 2 H + emission revealing the true centre of a prestellar core: further insights into deuterium chemistry': 'L. Pagani 1 , A. Belloche 2 , and B. Parise 2 *  \n- 1 LERMA & UMR8112 du CNRS, Observatoire de Paris, PSL Research University, CNRS, Sorbonne Universités, F- 75014 Paris, France\n- e-mail: laurent.pagani@obspm.fr\n- 2 Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany  \nReceived 04 / 07 / 2023; accepted ...', 'ABSTRACT': 'Context. IRAS 16293E is a rare case of a prestellar core being subjected to the e ff ects of at least one outflow.  \nAims. We want to disentangle the actual structure of the core from the outflow impact and evaluate the evolutionary stage of the core. Methods. Prestellar cores being cold and depleted, the best tracers of their central regions are the two isotopologues of trihydrogren cation which are observable from the ground, ortho-H2D + and para-D2H + . We used the Atacama Pathfinder EXperiment (APEX) telescope to map the para-D2H + emission in IRAS 16293E and collected James Clerk Maxwell Telescope (JCMT) archival data of ortho-H2D + . We compare their emission to that of other tracers, including dust emission, and analyse their abundance with the help of a 1D radiative transfer tool. The ratio of the abundances of ortho-H2D + and para-D2H + can be used to estimate the stage of the chemical evolution of the core.  \nResults. We have obtained the first complete map of para-D2H + emission in a prestellar core. We compare it to a map of orthoH2D + and show their partial anti-correlation. This reveals a strongly evolved core with a para-D2H + / ortho-H2D + abundance ratio towards the centre for which we obtain a conservative lower limit from 3.9 (at 12 K) up to 8.3 (at 8 K) while the high extinction of the core is indicative of a central temperature below 10 K. This ratio is higher than predicted by the known chemical models found in the literature. Para-D2H + (and indirectly ortho-H2D + ) is the only species that reveals the true centre of this core, while the emission of  \nother molecular tracers and dust are biased by the temperature structure that results from the impact of the outflow. This study invites to reconsider the analysis of previous observations of this source and possibly questions the validity isotopologue family. This  \nConclusions. of the deuteration chemical models or of the reaction and inelastic collisional rate coe ffi cients of the H + 3 could impact the deuteration clock predictions for all sources.  \nKey words. ISM: abundances - ISM: clouds - ISM: structure - Astrochemistry - Molecular processes - ISM: individual objects : IRAS 16293E', '1. Introduction': "IRAS 16293E is a well-known prestellar core (PSC) at a distance of 141 pc (Dzib et al. 2018) which is under the influence of at least one outflow coming from the multiple protostar system IRAS 16293-2422 (e.g. Wootten & Loren 1987, Castets et al. 2001, Lis et al. 2002, Stark et al. 2004). A second outflow, visible in the Spitzer / IRAC data and coming from the Young Stellar Object (YSO) WLY 2-69 was also reported by Pagani et al. (2016). This outflow points towards the south of the PSC, close to the HE2 position (Castets et al. 2001) but whether this represents a fortuitous alignment or reveals an actual impact is not yet established. IRAS 16293-2422 has attracted considerable attention as well, with the SMA (Jørgensen et al. 2011) and ALMA (program PILS, Jørgensen et al. 2016, and follow-up papers), and the first (and still unique) detections of para-H2D + and orthoD2H + in absorption in front of the Class 0 objects (Bruenken et al. 2014, Harju et al. 2017). IRAS 16293E has also gained new attention with the advent of ALMA observations combined with CSO maps (Lis et al. 2016) after the detection of ND3 (Roue ff et al. 2005), which confirm that most deuterated species such as  \nN2D + or DCO + peak at the outflow-PSC interface, with the exception of ND2H and ND3, the emission of which peak about 10 '' east of the DCO + emission peak. In a preliminary study we showed that Herschel continuum observations of the source with PACS and SPIRE reveal the temperature gradient induced by the outflow while ortho-H2D + emission traces the cold part of the core, away from the working surface, and approximately coincidental with the core seen in absorption in Spitzer IRAC midinfrared images (Pagani et al. 2016). The recent study by Kahle et al. (2023) of a large set of observations performed with the Atacama Pathfinder EXperiment (APEX) telescope towards both the YSOs and the PSC improves the census of molecular abundances and better di ff erentiates the various parts of this complex region, which confirms our understanding of its structure.  \nIt has been shown that the deuterium chemistry evolves in a rather monotonic way during the PSC formation phase and measuring the abundance of deuterated species could help to measure the time the cloud took to contract to the PSC core under scrutiny (Flower et al. 2006, Pagani et al. 2009, Pagani et al. 2013). Though it is easier to measure deuteration of species like NH3, HCO + , or N2H + , it has been advocated that measuring directly the abundance of the deuteration vector, namely the H + 3 isotopologues, would be more precise to characterize the evolu-  \ntion of the deuteration (Parise et al. 2011b, Harju et al. 2017). However single point measurements are often ambiguous, provide only line of sight averaged values, often decreasing the contrast between species. Only maps of the species allow us to fully understand the physics and chemistry of the PSCs.  \nIn this paper, we present the first map of para-D2H + , towards the prestellar core IRAS 16293E, and compare it with archival data of ortho-H2D + , dust emission and N2H + . Section 2 presents the observations. In Sect. 3, we analyse the emission of both deuterated isotopologues of H + 3 , and compare them to dust absorption and emission, and N2H + emission. We present a 1D modified Monte-Carlo model in Sect. 4 to quantify the abundance of both isotopologues, and discuss the various results in Sect. 5. Our conclusions are given in Sect. 6.", '2. Observations': "We observed para-D2H + (JK aKc : 110-101) at 691.6604434 GHz (Jusko et al. 2017) with the Swedish-ESO PI Instrument for APEX (SEPIA660) receiver on the APEX 12-m telescope equipped with Fast Fourier Transform (FFT) backends, each of 65536 channels of 61 kHz width. The observations were performed in several runs on ESO time (3 / 5 / 2021, 28-29 / 6 / 2021, project E-0105.C-0251A) and MPI time (2426 / 3 / 2022, 1 / 4 / 2022, 4 / 4 / 2022, project M-0109.F-9501C-2022). We also included observations with the 7 pixel CHAMP + camera taken on 7 / 9 / 2010 (one pointing, project M-085.F-00162010). The weather was exceptionally good (precipitable water vapour, PWV = 0.2-0.3 mm) during most of the SEPIA660 observations, and very good (PWV ≈ 0.5 mm) for the CHAMP + observations. Apart from the seven positions observed with CHAMP + and two positions observed with SEPIA660 in May 2021, which were observed in position switching mode, all observations were done in the OTF mode to cover a 50 '' × 50 '' region (June 2021) extended to 50 '' × 72 '' in 2022 as the preliminary results indicated an unexpected emission to the south. The various reference positions are summarized in Table A.1. The total observing time is about 25 hours. The angular resolution is ∼ 9 '' . The pointing was checked every hour, mainly on the nearby object RAFGL1922. The stability of the pointing is better than 3 '' . We used CLASS 1 to subtract first-order baselines computed on the velocity range -24 to + 31 km s -1 , and smoothed the data (spatially, using XY\\_MAP) to 14 '' and (spectrally) to 244 kHz ( ≈ 0.1 km s -1 ) to improve the sensitivity. The peak signalto-noise ratio (SNR) is 7 in brightness temperature and 12 in integrated intensity (in the range 2.81-4.18 km s -1 ). The rms noise level in integrated intensity is 6.7 mK km s -1 (T ∗ a scale) towards the centre of the map and reaches ≈ 11 mKkms -1 at the edges. We assumed a main beam e ffi ciency of η MB ≈ 0.48 ± 0.048 for CHAMP + observations 2 and η MB ≈ 0.46 ± 0.046 for SEPIA660 observations (Perez-Beaupuits, priv. comm.) which were applied to each set before merging them.  \nWe retrieved two sets of unpublished ortho-H2D + (JK aKc : 110-111) observations at 372.421340 GHz (Jusko et al. 2017) from the James Clerk Maxwell 15-m Telescope (JCMT) archive 3 , from projects M07AU29 (2007) and M09AU01 (2009). They were obtained with the 16 pixel Heterodyne Array Receiver Program (HARP) receiver with the Auto Correlation Spectral Imaging System (ACSIS) backend (61 kHz sampling). The atmospheric opacity at 225 GHz ranged from 0.052 to 0.034  \n(M07) and from 0.04 to 0.02 (M09). The M07 data consists of a 600-s long single pointing of 15 pixels (one corner detector was dead) towards the core centre. The M09 data (14 effective pixels) cover 232 '' × 113 '' , one square jiggle map centred on the core (1200 seconds integration) with an adjacent jiggle map to the East (300 seconds integration). The M09 central map was presented in our conference poster (Pagani et al. 2016). The original angular resolution of 13.2 '' was smoothed to 14 '' (after third-order baseline subtraction on a baseline extending from -21 to + 28 km s -1 , resampling and smoothing with CLASS / XY\\_MAP). The frequency sampling has been Hanningsmoothed to 122 kHz ( ≈ 0.1 km s -1 ). We assumed a main beam e ffi ciency of η MB ≈ 0.6 4 . The peak SNR is 11 in brightness temperature and 19 in integrated intensity (in the range 2.924.28 km s -1 ). The average integrated intensity rms noise is 23 mKkms -1 for the central map, twice as high for the eastern map of which we only use the edge to reach the emission limit of the core. The same observing run also provided us with N2H + ( J :4 3) observations at 372.6734633 GHz (Cazzoli et al. 2012), representative of many of the species observed towards the PSC (Lis et al. 2016, Kahle et al. 2023). Both D2H + and H2D + data sets are resampled on a 7 '' basis.  \nWe also used the JCMT SCUBA-II 850 µ m continuum observations obtained and published by Pattle et al. (2015). We refer the reader to that paper for a detailed description of the observations. We retrieved Spitzer IRAC data from the archive 5 . As far as we know, the PSC image from these observations has never been published apart from our conference poster (Pagani et al. 2016).", '3. Analysis': "Figure 1 shows the map of para-D2H + (left panel). This is the first time that a map of para-D2H + is presented. The doubly deuterated ion follows precisely the IRAC-1 absorption image seen at 3.6 µ m (central panel). Since D2H + is expected to peak in the central coldest part of clouds, and since the total absence of stars at 3.6 µ m reveals high extinction (typically AV ≥ 100 mag), these two maps do trace the centre of the IRAS16293E PSC, while the dust emission map peak (red cross), which is very sensitive to the temperature increase due to the outflow collision with the PSC edge, is displaced towards the impact point. Therefore the dust emission map does not reveal the true cold centre of the PSC. In the right panel, the upper part of the IRAC4 image of the PSC also seen in absorption is partially hidden (reduced contrast) by PAH emission in relation with another outflow mostly seen in the IRAC images. In the same panel, N2H + ( J :4 - 3) contours (in orange) also point towards the impact point, known as the DCO + peak (e.g. Lis et al. 2016).  \nThe left panel of Fig. 2 shows the ortho-H2D + emission superimposed on the IRAC-1 image. It is again clearly shifted to the East compared to the dust emission peak and N2H + ( J :4 - 3) line emission peak (marked by the red cross) as we already reported in Pagani et al. (2016) and has been confirmed by Kahle et al. (2023). It is remarkable that both dust and N2H + emissions peak at the same place, probably revealing the exact position of the outflow impact on the PSC. On the other hand, like for para-D2H + emission, dust absorption as seen by IRAC images, is correlated with ortho-H2D + emission. Dust in absorption is not sensitive to temperature e ff ects and this explains the  \nFig. 1. Continuum and line emission images of IRAS 16293E. Left: Colour map of para-D2H + integrated intensity emission with white contours (colour scale and white contours from 0.01 to 0.07 by 0.01 K km s -1 plus a contour at 0.075 K km s -1 ; contour spacing of 0.01 K km s -1 = 1.5 σ ) Centre: IRAC-1 image at 3.6 µ m. White contours trace para-D2H + as in the left panel. Red contours trace the SCUBA-II 850 µ m emission (20 to 160 by 20 [ ≈ 10 σ ] MJy sr -1 ). Right: IRAC-4 image at 8.0 µ m. White contours trace para-D2H + as in the left panel. Orange contours trace the N2H + ( J :4 - 3) integrated intensity (0.3 to 2.3 by 0.4 [ = 12 σ ] K km s -1 ). The yellow cross marks the centre of the D2H + emission ( α J2000: 16 h 32 m 30.51 s , δ J2000: -24 · 28 ' 53.7 '' , see also Fig. 2). The red cross marks the peak position of the dust and N2H + ( J :4 - 3) emissions. The black cross (left panel) marks the position of the observation of D2H + and H2D + by Vastel et al. (2004). The yellow arrow traces the direction of the outflow from the IRAS16293-2422 protostellar system and the red arrow the outflow from WLY 2-69. The grey and white disks represent the JCMT and APEX beams after smoothing to 14 '' .  \n<!-- image -->  \nFig. 2. Left: IRAC-1 image at 3.6 µ m. White contours trace orthoH2D + . The filled grey contour marks the strongest emission (contours from 0.10 to 0.34 by 0.06 [ = 2.6 σ ] K km s -1 ). Right: para-D2H + and ortho-H2D + maps. Same map of para-D2H + (black contours) as in Fig. 1 superimposed on the ortho-H2D + map (colour scale, plus white contours as in left panel). The yellow, black, and red crosses are as in Fig. 1. D2H + , H2D + , and N2H + spectra along the cut traced by the dashed line are displayed in Fig. 3  \n<!-- image -->  \ndecorrelation between dust absorption and emission. The dust absorption reveals the highest extinction part of the PSC and as expected this coincides with the emission peak of ortho-H2D + and para-D2H + . Another noticeable feature is the horse-shoe shape of the ortho-H2D + emission surrounding the PSC centre (the topmost grey-filled contour in Fig. 2, left side). This  \northo-H2D + horse-shoe shape clearly surrounds the emission peak of para-D2H + as shown in the right panel of Fig. 2. This indicates that the ortho-H2D + emission is weaker towards the centre of the PSC. Since it is only the top contour of ortho-H2D + emission which shows this feature, it means that only the inner part of the PSC shows a diminished emission and therefore a diminished ortho-H2D + abundance, compared to the outer layers. The intensity di ff erence between the ortho-H2D + peak emission at a right ascension o ff set of + 7 '' from the centre and the central position is 66 mK km s -1 , with a statistical noise of 29 and 21 mKkms -1 , respectively. The di ff erence is 2.3 σ , but if the ortho-H2D + abundance had been constant, the emission would have increased towards the centre, like for para-D2H + , instead of decreasing. Therefore the di ff erence with respect to a situation where the abundance would be constant is much higher than 2.3 σ and statistically significant. In Sect. 4, we will model the case of a constant ortho-H2D + abundance to quantify the expected di ff erence.  \nThe pointed observations of ortho-H2D + and para-D2H + that were previously reported by Vastel et al. (2004) towards the position marked by a black cross in Figs. 1 & 2 are found to be inconsistent with our more sensitive mapping observations towards the same position. For para-D2H + , they reported a peak temperature of 0.34 K ( ± 0.077 K ?) 6 and a linewidth of 0.29 ± 0.07 km s -1 , to be compared with our 83 ± 17 mK and a linewidth of 0.48 ± 0.07 km s -1 ; for ortho-H2D + , they report a peak temperature of 1.31 K ( ± 0.22 K ?) 6 and a linewidth of  \n0.36 ± 0.04 km s -1 , to be compared with our 0.31 ± 0.03 K and a linewidth of 0.61 ± 0.04 km s -1 . For a gas temperature in the range 16-25 K from former studies (Loinard et al. 2001, Castets et al. 2001, Lis et al. 2002, Stark et al. 2004, Lis et al. 2016) that applies to the DCO + emission peak region, the thermal linewidth alone is expected to be between 0.38 and 0.48 km s -1 for paraD2H + and between 0.43 and 0.54 km s -1 for ortho-H2D + , which are both inconsistent with the small linewidths reported by Vastel et al. (2004).", '4. Radiative transfer modelling of para-D 2 H + and ortho-H 2 D +': "Modelling the abundance of both H + 3 isotopologues is needed to quantify our findings. The shape of the PSC is not symmetrical but still relatively simple. There is however no easy way to derive its density and temperature structures, due to the complexity of its environment. Star counts or reddening measurements are impossible towards the centre of the core due to the lack of stars even in the sensitive Spitzer data (Fig. 1). 8 µ m continuum absorption measurements (Bacmann et al. 2000, Lefèvre et al. 2016) are hampered by the superposition of nebulosity emission from the IRAS 16293-2422A / B outflows and possibly some Polycyclic Aromatic Hydrocarbon (PAH) emission. All emissions, whether dust or molecular, are strongly biased by the temperature gradient induced by the outflow impact on the PSC, and the density profile is also possibly perturbed on that side, away from a simple Plummer-like profile (Whitworth & Ward-Thompson 2001). To disentangle density, temperature, and abundances for molecules, or dust properties for continuum emission, a thorough analysis in 3D of the core must be performed. This is beyond the scope of the present paper and will be addressed in the future. Here, we present a model combining three partially di ff erent one-dimensional (1D) spherical models.  \nThe first step is to build 1D gas density and temperature structures of the PSC. The only tracer from the edge to the centre of the PSC is dust, especially in emission around 300 GHz but its properties are not accurately known, and since the temperature and density profiles cannot be easily retrieved from the dust emission alone (Pagani et al. 2015), we have constructed the core structure from the following considerations:  \n- -except for the south-western extension, the emission of both para-D2H + and ortho-H2D + is relatively round and symmetrical. We consider that the bulk of the original core, before the outflow impacts its western side, is still present and close enough to a spherical object.\n- -The eastern part of the core is unperturbed and its H2 density can be described with a Plummer-like profile (Tafalla et al. 2002):  \nn(H2)(r) = 1 . 5 × 10 6 1 + ( r 3 × 10 16 ) 2 . 5 cm -3 . (1)  \nThe peak column density is 1.3 × 10 23 cm -2 , i.e. a visual extinction of AV ≈ 140 mag, compatible with the total absence of stars seen through the core even in the Spitzer MIR maps. It is also consistent with the dust emission at 850 µ mtowards the same spot (which is a weak constraint given the large uncertainty on the dust emissivity, Demyk et al. 2017a,b).  \nWe also apply this profile to the West side. Using a di ff erent density profile would locally change the abundance profile of the ions but would not impact the results for the central region of the core on first approximation.  \n- -In the centre of cold cores, temperatures are usually in the 8 - 12 K range (Caselli et al. 2008). We thus considered central temperatures of 8, 10, and 12 K and imposed a temperature of 12 K in the external layers on the eastern side (the exact value has no impact in the 11-16 K range).\n- -The N2H + ( J :4 - 3) emission peaks about 20 '' west of the PSC centre in a place that is warmed by the outflow impact. Former studies proposed temperatures between 16 and 25K (Loinard et al. 2001, Castets et al. 2001, Lis et al. 2002, Stark et al. 2004, Lis et al. 2016), and we display results for the median value of 20K but also checked the 16 and 25 K cases.\n- -The para-D2H + spectra from -14 '' to -28 '' (and possibly -35 '' ) in Fig. 3 are slightly stronger than their symmetrical ones eastward, and since the position of the peak of N2H + emission is in the middle of this range, we set the temperature increase at -20 '' .\n- -We divided the sphere in three parts: each time we ran a 1D spherical model but we extracted the emission from a fraction of the sphere only. We stitched the three extracted portions together and convolved the output to 14 '' . This is graphically explained in Appendix B. The spherical model is discretized into 20 layers of 1.14 × 10 16 cm thickness. This corresponds to 5.4 '' thickness per layer for a distance of 141 pc. The total size is 216 '' , which completely covers the core.  \nThough the warm part could be di ff erent in its density profile due to the possible compression from the outflow, if any (the outflow could be grazing the PSC and have a limited compression e ff ect only), it is not possible to assess the actual profile at this stage until we have built a 3D model to explore that possibility. However the warm part is far enough from the horse-shoe central region and has no influence on the modelling of the cloud centre, as our tests have shown. Similarly, the impact of the temperature between 16 and 25 K was found to be negligible.  \nThe modelled density, temperature, and isotopologue abundance profiles are displayed in Fig. 4. We adjusted the abundance, non-thermal linewidth and radial velocity to fit the widths and strengths of the ortho-H2D + and para-D2H + lines (see synthetic spectra in Fig.3). Linewidths show some inconsistencies as revealed by Gaussian fits to the observed spectra (for more details, see Appendix C). This pseudo-1D modelling is performed with a 1D Monte-Carlo code adapted from Bernes (1979)'s radiative transfer program. For the non-LTE modelling of the two ions, we used the inelastic collisional rate coe ffi cients from Hugo et al. (2009). The abundances of both ions relative to H2 are displayed in the central panel of Fig. 4 while their ratio is displayed in the bottom panel.  \nThe modelled lines are displayed in Fig. 3 and the same fit is traced in terms of integrated intensity in Fig. 5, where the centre of the ortho-H2D + horse-shoe shape is clearly visible. In Fig. 5, the observed intensities (full lines, with their individual error bars) are reproduced within 1 σ by the model (dashed line, same colour). Because the density keeps increasing towards the centre, the ortho-H2D + intensity drop needs a strong abundance decrease to explain the decorrelation between the two (unlike para-D2H + , the intensity of which keeps increasing). Since the para-D2H + abundance is found to be almost constant in the centre of the cloud (between 1.15 and 1.7 × 10 -10 , see below), we tested the case of a constant abundance for ortho-H2D + . The predicted intensity is shown in Fig. 5 with a dashed orange line. Towards the central position, the di ff erence reaches 7 σ and a uniform abundance of ortho-H2D + is clearly excluded.  \nFig. 3. From top to bottom: N2H + ( J :4 - 3), ortho-H2D + (JK aKc : 110-111) and para-D2H + (JK aKc : 110-101) spectra in black, extracted along the cut displayed in Fig. 2, with best-fit model of ortho-H2D + and para-D2H + overlaid in red. The three weakest para-D2H + spectra have been smoothed to 0.21 km s -1 . The velocity axis in each box runs from 1 to 6 km s -1 . Green dashed lines mark the 3 σ individual noise levels and the systemic velocity (3.6 km s -1 ).  \n<!-- image -->  \nFig. 4. Pseudo-1D model profiles of the core. Upper box: n(H2) Plummer-like profile in black. In red, temperature profile treated separately for the East (left) and the West (right) parts of the core, with three di ff erent minimum temperatures (10 K in plain line, 8 and 12 K in dashed lines). Middle box: ortho-H2D + (black) and para-D2H + (red) abundance profiles relative to H2 and computed separately for each side for the 10K model, error bars are 1 σ deviation from best model. Arrows mark upper limits. In particular, for ortho-H2D + in the centre of the core, the dotted line represents the abundance to stay just below 1 σ di ff erence with the observations, and the upper error bar represents a di ff erence of 2 σ with the observations. Lower box: para-D2H + / orthoH2D + ratio for the 10K model. Error bars are computed by taking the opposite error bars of both species when available, otherwise lower or upper limits are indicated by arrows. The long dashed and dash-dotted lines mark a ratio of 1 and 2, respectively, 2 being the upper limit of the models of Bovino et al. (2021) at 10 K.  \n<!-- image -->", '5. Discussion': "The H2D + / D2H + contrast reported in Sect. 3 is expected from chemical evolution of the deuteration of the H + 3 ion (Roberts et al. 2003, Walmsley et al. 2004) that indicates that H + 3 can be converted to H2D + , itself converted to D2H + and finally to D + 3 , the latter being potentially the most abundant ion in the end. The e ff ect is maximum at the peak density because the higher the density, the faster the chemical reactions, other parameters being equal. This is the first time that we directly witness the intermediate step of this transformation (H + 3 and D + 3 are not observable in dense clouds). This confirms that the original centre of the PSC is not coincident with the peak emission of species such as N2H + , N2D + , and DCO + , but 22 '' apart, which di ff ers from what is usually found in other objects (e.g., Tafalla et al. 2002, Pagani et al. 2007). This must result from the temperature enhancement produced by the molecular outflow hitting the PSC. This heating possibly explains the similar extent of H2D + and D2H + towards the working surface of the PSC, in the southwest direction. It may also explain the small secondary peak of para-D2H + as a result of the temperature increase compensating the drop of abundance, although further work is needed to reach a firm conclusion.  \nFig. 5. Observed intensities along the horizontal cut (full lines with their error bar, corresponding to the spectra from Fig. 3) and corresponding model intensities (dashed lines). For ortho-H2D + , we present an additional model (orange dashed line) for which the abundance is kept constant across the core. Para-D2H + (in red) observations and model have their own intensity axis, also in red, on the right side of the plot.  \n<!-- image -->  \nThe column densities of both species derived from the modelling in Sect. 4 are nearly equal towards the centre, N(orthoH2D + ) peaks at 2.1- 4.8 × 10 13 cm -2 (at 12 K - 8 K), and N(para-D2H + ) at 0.94 - 4.2 × 10 13 cm -2 . The maximum paraD2H + / ortho-H2D + column density ratio is 1.9 + 2 . 1 -1 . 2 , (at 10 K, ranging from 1.1 at 12 K to 3.1 at 8 K) 7 . However, the derived abundance profiles yield a stronger contrast in the inner shells as demonstrated in Fig. 5. For the two most inner shells (total radius of 10.8 '' ), the abundance of ortho-H2D + can be set to zero while the line of sight emerging intensity computed by the model is still too high (but less than 1 σ away from the observations). Since we asymptotically reach a di ff erence of + 0.8 σ for the central position by decreasing the abundance by almost two orders of magnitude, we selected the abundance that provides exactly 1 σ di ff erence as the reference and computed the upper limit of this abundance by setting the error to 2 σ . There is no lower limit since we cannot make the model weaker than the observations for this central spectrum. This is because the large ortho-H2D + abundance needed to reproduce the intensity of the line towards the walls of the horse-shoe cavity is also su ffi cient to reproduce the intensity of the line for the two sightlines across the horse-shoe hole, the same walls probably standing in front and in the back of the cavity, the observer having no privileged point of view (a situation similar to the case of L1506C, Pagani et al. 2010). In the centre of the cavity we find a very high para-D2H + / ortho-H2D + abundance ratio of ∼ 20. To check whether this ratio could be compatible with known models when taking into account the uncertainties of the line modelling, we computed a conservative lowest value of the uncertainty range of this ratio. We estimated an upper limit to the abundance of ortho-H2D + at 2 σ above the observations and found abundances of 1.8 × 10 -11 at 12 K, 3.05 × 10 -11 at 10 K, and 4.0 × 10 -11 at 8 K. We also computed a lower abundance for para-D2H + by departing from the best fit abundance by -1 σ (X[para-D2H + ] = 6.95 × 10 -11 , 1.3 × 10 -10 , and 3.3 × 10 -10 at 12, 10, and 8K respectively). The para-D2H + / ortho-H2D + ratio lowest limit ranges from 3.9 to 8.3 (at 12 K - 8 K).  \nThis para-D2H + / ortho-H2D + abundance ratio is high and beyond all predictions made by various models whether 0D (Pagani et al. 2009, Parise et al. 2011a), 1D (Pagani et al. 2013) or 3D (Bovino et al. 2019, 2021) which reach values around 1-2 in the densest parts of the cores, for a temperature of 10 K. Parise et al. (2011a) also found a too high ratio when analysing the HMM1 prestellar core in the LDN 1688 star forming region with a low temperature (7 K) but could reconcile the model and observations at 12 K. Here, the 12 K value is possibly marginally consistent with the Bovino models if we interpolate by eye between their 10 and 15 K cases (Bovino et al. 2021, their Fig. B.2). However, clouds with AV ≥ 100 mag. should be very cold in their centre, below 10 K following the models of Zucconi et al. (2001) and observations of Crapsi et al. (2007) and Pagani et al. (2007). These observations coupled to the present simple 1D modelling, if confirmed by a more accurate 3D modelling, either suggest that current chemical models are not yet close enough to reality to trace correctly the deuteration evolution of the tri-hydrogen ion in the cold medium or that the reaction rate coe ffi cients or inelastic collisional rate coe ffi cients of the H + 3 family with H2 and HD from Hugo et al. (2009), which were computed semiclassically, are not accurate enough and need to be improved. We found that increasing the collision rate coe ffi cients of paraD2H + with H2 by a factor 6 or increasing the ortho-H2D + + HD → para-D2H + + para-H2 reaction rate coe ffi cient by a factor 3 would change the para-D2H + / ortho-H2D + ratio by a factor 2. Other similar changes can also be considered.", '6. Conclusions': "We present in this article the first map of para-D2H + in the interstellar medium, obtained with the APEX telescope toward the prestellar core IRAS 16293E. The main conclusions of this work are:  \n- 1. Maps of para-D2H + , and to a lesser degree, of ortho-H2D + reveal the position of the core centre before it was compressed and heated by at least one outflow coming from nearby YSOs. The new reference position for the core centre is α J2000: 16 h 32 m 30.51 s , δ J2000: -24 · 28 ' 53.7 '' .\n- 2. We derived peak column densities of ∼ 1-3 × 10 13 cm -2 for both ions (at 10 K).\n- 3. We derived a conservative lower limit of the paraD2H + / ortho-H2D + abundance ratio in the centre of the PSC in the range 3.9-8.3 at 12 - 8 K while the core extinction (AV ≥ 100 mag) is indicative of a central temperature below 10 K. This ratio is higher than what published chemical models predict for temperatures below or around 10 K.\n- 4. Temperature e ff ects induced by the outflow impact at the hotspot make line emission of most species bright. Therefore an analysis limited to the hotspot could miss the real peak abundance of species, in particular ND3 and other deuterated species, the emission of which peaks away from the orthoH2D + and para-D2H + abundance peaks.  \nTo confirm the present model, a thorough modelling of the core to disentangle density and temperature variations induced by the outflow impact needs to take into account all available information (line emissions with several transitions per species, dust emission and absorption lower limit) in a 3D model, though the main di ffi culty of this model will be to constrain the line of sight distribution of all parameters.  \nThe high para-D2H + / ortho-H2D + ratio measured in IRAS 16293E could suggest, if confirmed, that the inelastic collisions  \nor reaction rate coe ffi cients of the H + 3 isotopologues need to be revised or that the chemical models need to be refined in order to improve our understanding of the deuteration process in star forming regions. This might have an impact on the deuteration chemical clock used to study many star forming regions.  \nAcknowledgements. Based on observations with the Atacama Pathfinder EXperiment (APEX) telescope. At the time of the first observations presented in this paper, APEX was a collaboration between the Max Planck Institute for Radio Astronomy, the European Southern Observatory, and the Onsala Space Observatory. This research has made use of NASA's Astrophysics Data System Abstract Service. It has also made use of observations from Spitzer Space Telescope and from the NASA / IPAC Infrared Science Archive, which are operated by the Jet Propulsion Laboratory (JPL) and the California Institute of Technology under contract with NASA. This research makes a large use of the CDS (Strasbourg, France) services, especially Aladin (Bonnarel et al. 2000), Simbad and Vizier Ochsenbein et al. (2000). We are indebted to C. de Breuck for the observations in the ESO time and to all the observers at APEX. We thank P.F. Goldsmith for useful comments and M. Tafalla for his helpful advices and editorial work.", 'References': 'Bacmann, A., André, P., Puget, J.-L., et al. 2000, A&A, 361, 555 Bernes, C. 1979, A&A, 73, 67  \nBonnarel, F., Fernique, P., Bienaymé, O., et al. 2000, A&AS, 143, 33 Bovino, S., Ferrada-Chamorro, S., Lupi, A., et al. 2019, ApJ, 887, 0 Bovino, S., Lupi, A., Giannetti, A., et al. 2021, A&A, 654, A34 Bruenken, S., Sipilä, O., Chambers, E. T., et al. 2014, Nature, 516, 219 Caselli, P., Vastel, C., Ceccarelli, C., et al. 2008, A&A, 492, 703  \n- Castets, A., Ceccarelli, C., Loinard, L., Caux, E., & Lefloch, B. 2001, A&A, 375, 40\n- Cazzoli, G., Cludi, L., Bu ff a, G., & Puzzarini, C. 2012, ApJ Supplement Series (ISSN 0067-0049), 203, 11\n- Crapsi, A., Caselli, P., Walmsley, M. C., & Tafalla, M. 2007, A&A, 470, 221 Demyk, K., Meny, C., Leroux, H., et al. 2017a, A&A, 606, A50\n- Demyk, K., Meny, C., Lu, X. H., et al. 2017b, A&A, 600, A123\n- Dzib, S. A., Ortiz-León, G. N., Hernández-Gómez, A., et al. 2018, A&A, 614, A20\n- Flower, D., Pineau Des Forêts, G., & Walmsley, C. M. 2006, A&A, 449, 621 Harju, J., Sipilä, O., Brünken, S., et al. 2017, ApJ, 840, 63\n- Hugo, E., Asvany, O., & Schlemmer, S. 2009, JCP, 130, 164302\n- Jusko, P., Töpfer, M., Müller, H. S., et al. 2017, Journal of Molecular Spectroscopy, 332, 33\n- Jørgensen, J. K., Bourke, T. L., Nguyen Luong, Q., & Takakuwa, S. 2011, A&A, 534, 100\n- Jørgensen, J. K., van der Wiel, M., Coutens, A., et al. 2016, A&A, 595, A117 Kahle, K. A., Hernández-Gómez, A., Wyrowski, F., & Menten, K. M. 2023, A&A, 673, A143\n- Lefèvre, C., Pagani, L., Min, M., Poteet, C., & Whittet, D. C. B. 2016, A&A, 585, L4\n- Lis, D. C., Gerin, M., Phillips, T. G., & Motte, F. 2002, ApJ, 569, 322\n- Lis, D. C., Wootten, H. A., Gerin, M., et al. 2016, ApJ, 827, 1\n- Loinard, L., Castets, A., Ceccarelli, C., Caux, E., & Tielens, A. G. G. M. 2001, ApJ, 552, L163\n- Ochsenbein, F., Bauer, P., & Marcout, J. 2000, A&AS, 143, 23\n- Pagani, L., Bacmann, A., Cabrit, S., & Vastel, C. 2007, A&A, 467, 179\n- Pagani, L., Lefèvre, C., Belloche, A., et al. 2016, EAS Publications Series, 7576, 201\n- Pagani, L., Lefèvre, C., Juvela, M., Pelkonen, V. M., & Schuller, F. 2015, A&A, 574, L5\n- Pagani, L., Lesa ff re, P., Jorfi, M., et al. 2013, A&A, 551, 38\n- Pagani, L., Ristorcelli, I., Boudet, N., et al. 2010, A&A, 512, A3 +\n- Pagani, L., Vastel, C., Hugo, E., et al. 2009, A&A, 494, 623\n- Parise, B., Belloche, A., Du, F., Güsten, R., & Menten, K. M. 2011a, A&A, 526, A31\n- Parise, B., Belloche, A., Du, F., Güsten, R., & Menten, K. M. 2011b, A&A, 528, C2\n- Pattle, K., Ward-Thompson, D., Kirk, J. M., et al. 2015, MNRAS, 450, 1094 Roberts, H., Herbst, E., & Millar, T. J. 2003, ApJ, 591, L41\n- Roue ff , E., Lis, D. C., van der Tak, F. F. S., Gerin, M., & Goldsmith, P. F. 2005, A&A, 438, 585\n- Stark, R., van der Tak, F. F. S., Sandell, G., et al. 2004, ApJ, 608, 341\n- Tafalla, M., Myers, P. C., Caselli, P., Walmsley, C. M., & Comito, C. 2002, ApJ, 569, 815\n- Vastel, C., Phillips, T. G., & Yoshida, H. 2004, ApJ, 606, L127\n- Walmsley, C. M., Flower, D. R., & Pineau des Forêts, G. 2004, A&A, 418, 1035 Whitworth, A. P. & Ward-Thompson, D. 2001, ApJ, 547, 317, aDS Bibcode: 2001ApJ...547..317W\n- Wootten, A. & Loren, R. B. 1987, ApJ, 317, 220  \nZucconi, A., Walmsley, C. M., & Galli, D. 2001, A&A, 376, 650', 'Appendix A: Reference positions for the APEX para-D 2 H + observations': 'Di ff erent reference positions were used for the various projects with APEX towards IRAS 16293E. We list them in Table A.1', 'Appendix B: Model building': "Fig. B.1. Steps to build an asymmetrical model from 1D spherical models. Based on a single density profile (eq. 1), we compute two models with di ff erent temperature profiles in the outer parts (either 12 or 20 K, top row, see Fig. 4). We apply two di ff erent abundance profiles (central row) for each ion for the East and West parts, the West part is computed for both temperature profiles. We extract the parts corresponding to each region of interest (sampled with a 1 '' resolution) and we assemble them in a cube (RA, Dec, Frequency, bottom row left). Finally, we convolve the cube with a 14 '' half-power beamwidth Gaussian beam to compare to the observations (bottom row right). We extract the spectra along the dash-dot line and compare them to the observations (see Fig. 3).  \n<!-- image -->  \nThe model is based on a single density profile (eq. 1) and we computed the emission on a 1 '' step basis (model diameter of 216 '' ) for various temperature and abundance profiles.  \n- 1. for positive o ff sets (East side) we used the pristine model and extracted the emission from the half sphere towards the East.\n- 2. for o ff sets between 0 and -15 '' , we used the same model but on the West side with a di ff erent abundance profile except in the central shell which is forcibly common to both sides, being unresolved, (the centre being common to the East and West sides, we technically use only the West side and therefore the East stops at 5 '' , and the West goes from -15 to + 5 '' )\n- 3. for the external heated part, selected beyond -15 '' from the centre, we used a 1D model with warm outer layers beyond a radius of 20 '' in the Western part.  \nFor the ortho-H2D + case shown in Fig. B.1, the strong dissymmetry is due to the fact that the emission is extended to the East  \nand not as much to the West. The ortho-H2D + abundance profiles (Fig. 4, central panel) reflect this di ff erence by decreasing beyond -30 '' on the Western side, while it keeps increasing beyond + 30 '' on the Eastern side. One should note that the warm layers beyond -20 '' contribute to positions closer to the cloud centre, because being a 1D onion model, they are present on the line of sight of these closer positions in front and behind. The line of sight crosses a cold layer inserted between two hot layers. This could be justified by the outflow interacting with the surface of the PSC rather than ploughing its way to the centre, since the ortho-H2D + and para-D2H + lines do not seem to be kinematically perturbed in this warmer region. This problem will be examined in detail in a future work. The model is not representative of the emission towards the North or the South of the PSC.", 'Appendix C: Linewidth and Gaussian fits of the D 2 H + APEX and H 2 D + JCMT spectra': "The lines of both ions are optically thin (peak opacity of 0.6 and 0.7 for para-D2H + and ortho-H2D + , resp.), and in first approximation consistent with a Gaussian shape. We have therefore adjusted a Gaussian to each of the spectra. These Gaussian fits are displayed in Fig. C.1. To adjust the Monte-Carlo model, we modified the isotopologue abundance profile across the core to minimize the integrated intensity di ff erence between the observations and the modelled line. Since the lines are optically thin, the exact shape of the line has little impact on its integrated intensity but of course, we take care that the linewidth be compatible with the observations. However, from the Gaussian fit, we found that the linewidth of the two isotopologues are not equal. At the edges, this is probably a noise limitation (e.g. the 21 '' para-D2H + spectrum), but in the centre, the para-D2H + spectra at o ff sets 0, -7 '' , and -14 '' are clearly larger than the ortho-H2D + ones, e.g., for the 0 position:  \n- -δ V(para-D2H + ) = 0.658 ± 0.074 km s -1\n- -δ V(ortho-H2D + ) = 0.524 ± 0.037 km s -1  \nH2D + being lighter than D2H + , one would expect the former to display lines slightly wider than the latter for similar thermal turbulence or macroscopic movements. The reverse could indicate that the core inside the horse-shoe, depleted in H2D + (see Sect. 3), has a signature of collapse of its own traced by D2H + alone, but we did not manage to reproduce this di ff erence with our model. The volume is probably too small to dominate the linewidth formation without perturbing the H2D + linewidth too. The other possibility is that this di ff erence be mostly due to noise (1.8 σ di ff erence) and we settled the model linewidth to be intermediate between these two values, making the D2H + modelled line slightly stronger than the Gaussian fit but narrower, and conversely for H2D + .  \nTable A.1. Reference positions for the APEX observationsNotes. ( a ) O ff sets are with respect to equatorial coordinates except if marked with (H) for horizontal wobbler switching  \nFig. C.1. JCMT spectra of ortho-H2D + J KK ' :110-111 (top row, velocity axis from 1 to 6 km s -1 ), and APEX spectra of para-D2H + J KK ' :110-101 (middle row, velocity axis from 1 to 6 km s -1 ) across the cloud (the position of the cut is indicated in Fig. 2). Gaussian fits are overlaid in red and green, respectively. The same Gaussian fits with same colour code are reproduced in the bottom row, after normalization to the peak temperature for the + 7 '' spectrum, for a better comparison of the ratio variation and linewidth di ff erences between the two isotopologues.  \n<!-- image -->"}