A high HDO/H2O ratio in the Class I protostar L1551 IRS5

Context. Water is a very abundant molecule in star-forming regions. Its deuterium fractionation provides an important tool for understanding its formation and evolution during the star and planet formation processes. While the HDO / H 2 O abundance ratio has been determined toward several young Class 0 protostars and comets


Introduction
Water is one of the most abundant molecules in the interstellar medium.Its role is significant in the star formation process as it cools down the warm regions and boosts the gravitational collapse.It has been detected in multiple environments, such as protostars (e.g., van Dishoeck et al. 2021), protoplanetary disks (e.g., Hogerheijde et al. 2011), comets (e.g., Mumma et al. 1986;Davies et al. 1997), and asteroids (e.g., Rivkin & Emery 2010).Deuteration is extremely sensitive to density and temperature conditions (Ceccarelli et al. 2014), and higher D/H ratios are expected in particular when water forms in cold and dense con-This work is based on observations carried out under project numbers W18AO and S16AE with the IRAM NOEMA Interferometer.IRAM is supported by INSU/CNRS (France), MPG (Germany), and IGN (Spain).
ditions (e.g., Coutens et al. 2012;Jensen et al. 2019).The comparison of the D/H abundance ratio of the terrestrial oceans (Vienna Standard Mean Ocean Water (VSMOW) ∼ 1.557 × 10 −4 , de Laeter et al. 2003) with comets and asteroids (a few 10 −4 ) suggests that these primitive objects could have delivered part of the water found on Earth (e.g., Hartogh et al. 2011;Altwegg et al. 2015;Trigo-Rodríguez et al. 2019).Chemical models of water deuteration in protoplanetary disks also show that the fractionation is not efficient enough to reproduce the terrestrial water D/H ratio without considering an inheritance from the parental molecular cloud (Cleeves et al. 2014).Studying the water deuteration in protostars can consequently help understand the origin of water in disks and planets (Furuya et al. 2016).
To reconstruct the full history of water during the star formation process and bridge the gap between the Class 0 protostellar stage and the formation of planets, comets, and asteroids, studies of more evolved protostars are needed.Even if the detection of water is usually difficult in such sources (Harsono et al. 2020), there are a few Class I objects in which water can be detected, for example after a strong increase in the accretion rate from the disk onto the protostar (e.g., FUor or EXor sources, Audard et al. 2014;Fischer et al. 2022).Such an event leads to a luminosity burst that warms up the disk and the envelope in which the icy grain mantles can sublimate (Cieza et al. 2016).However, the studies were limited to only two sources.A first study by Codella et al. (2016) showed the detection of HDO toward the EXor-like source SVS13-A with the IRAM 30 m telescope.The presence of high-energy transitions suggests that they come from the warm inner regions even though emission from largescale shocks cannot be excluded.From the nondetection of the H 2 18 O transition at 203.4 GHz, they estimated a lower limit of the HDO/H 2 O ratio of 10 −3 .More recently, a HDO/H 2 O ratio of (2.26 ± 0.63) × 10 −3 was obtained with ALMA in the disk of the FUor protostar V883 Ori (Tobin et al. 2023).The similarity of the water deuteration in isolated Class 0 protostars supports the inheritance of water from the envelope to the disk.
Here we report the interferometric determination of the HDO/H 2 O ratio in another Class I source, L1551 IRS5 (Chen et al. 1995).This protostellar binary, located in the Taurus molecular cloud (d ∼ 155 pc, Zucker et al. 2019), is classified as a FUor-like object by Connelley & Reipurth (2018) based on its near-infrared spectra (Mundt et al. 1985;Carr et al. 1987).The detection of complex organic molecules in its inner regions suggests high-temperature conditions that would be sufficient to thermally desorb the molecules frozen on grains, similarly to hot corinos (Bianchi et al. 2020;Cruz-Sáenz de Miera, in prep.).

Observations
Observations of two HDO transitions and one H 2 18 O transition were carried out toward L1551 IRS5 with the Northern Extended Millimeter Array (NOEMA) through the projects W18AO and S16AE (see Table 1).The HDO 3 1,2 -2 2,1 and H 2 18 O 3 1,3 -2 2,0 transitions at 225.9 GHz and 203.4 GHz, respectively, were observed simultaneously in the A configuration on 15 and 17 February 2019 and in the C configuration on 27 February 2019.High spectral resolution windows with channel spacing of 62.5 kHz (i.e., 0.08-0.09km s −1 ) were used to cover the two transitions.We resampled the data with a resolution dv of 0.5 km s −1 to increase the signal-to-noise ratio.The HDO 4 2,2 -4 2,3 transition at 143.7 GHz was covered with a poor resolution of 1.95 MHz (i.e., 4 km s −1 ) in other NOEMA observations carried out on 25 August and 31 August 2016 in the D configura-tion.All the data were calibrated and imaged using the CLIC and MAPPING packages of the GILDAS2 software.The list of calibrators for each observation is in Table A.1.The continuum was subtracted before the Fourier transform was applied to the line data.The data were imaged with the natural weighting mode.The synthesized beam sizes and rms obtained for each transition are listed in Table 1.The two components, N and S, of L1551 IRS5 are not fully resolved in our observations as they are only separated by 0.36 (i.e., ∼56 AU, Cruz-Sáenz de Miera et al. 2019).

Results
In a first step, we visualized the datacubes of the HDO and H 2 18 O lines with the Aladin 3 and CASSIS 4 softwares (Bonnarel et al.  2000; Glorian et al. 2021).We extracted the spectra at the positions of the northern source (α J2000 = 04 h 31 m 34 s .161,δ J2000 = +18 • 08 04.72 ) and the southern source (α J2000 = 04 h 31 m 34 s .165,δ J2000 = +18 • 08 04.359 , Takakuwa et al. 2020).To study variations in the line intensities with the position, the spectra were also extracted for a position between the two sources, at an equal distance from them (α J2000 = 04 h 31 m 34 s .163,δ J2000 = +18 • 08 04.550 ), that we call M. Because of the asymmetry in emission maps with respect to M, we include a position at a similar distance from N in the opposite direction (α J2000 = 04 h 31 m 34 s .159,δ J2000 = +18 • 08 04.890 ) called O. A comparison of these spectra is shown in Fig. 1, and their locations are indicated in Fig. 2.
Two velocity components are seen in the line profiles at 203.4 and 225.9 GHz, similarly to Bianchi et al. (2020) and Mercimek et al. (2022) for complex organic molecules.The component at ∼9 km s −1 is brighter at the N and O positions than at the M and S positions.The HDO component at ∼6 km s −1 is slightly brighter at the M position than the N and S positions, while it is less bright at the O position.The high-velocity peak corresponds to the N source or farther to the north, in agreement with Bianchi et al. (2020).The HDO component at 6 km s −1 seems to arise closer to the M position than the S position based on the comparison of the observed line intensities.The double peak created by these two velocity components could then possibly indicate rotation from a disk surrounding the N source.For the HDO line at 143.7 GHz, the spatial resolution is too low (∼2.0-2.5 ) compared to the binary separation (0.36 ) to observe spectral differences at the four positions (see Fig. 1).The two velocity components are not observed because of the low spectral resolution of the data (4 km s −1 ).
We used the CASSIS software to check possible blending with other molecules at velocities of both 6 km s −1 and 9 km s −1 .The spectroscopic data come from the CDMS and JPL databases (Pickett et al. 1998;Müller et al. 2005).No blending is observed for the component at ∼9 km s −1 either for HDO or H 2 18 O.However, the H 2 18 O component at ∼6 km s −1 is blended with a CH 3 OCH 3 3 3,1,1 -2 2,1,1 line emitted at ∼9 km s −1 .A faint line of CH 3 OCHO 6 6,0 -5 5,0 (emission at 9 km s −1 ) is observed in the wing of the HDO transition at 225.9 GHz.The HDO line at 143.7 GHz does not seem to be affected by blending.
Integrated emission maps of the HDO 225.9 GHz and H 2 18 O 203.4 GHz lines are shown in Fig. 2. The water emission is compact (less than 70 AU in diameter), as also seen for complex   organic molecules (Bianchi et al. 2020).As expected, the HDO component at 9 km s −1 (in red) peaks close to the N position, while the one at 6 km s −1 (in blue) peaks closer to the M position.A similar distribution is seen for the H 2 18 O component at 9 km s −1 .However, the one at 6 km s −1 peaks closer to the N position.This can be explained by the blending with the CH 3 OCH 3 line, which is emitted at a velocity of 9 km s −1 .It consequently means that the CH 3 OCH 3 emission contributes significantly to the flux observed at 6 km s −1 .The HDO and H 2 18 O emission sizes were estimated with circular Gaussian fitting of the lines in the (u,v)-plane.We found approximate source sizes of 0.35 and 0.45 for the components at 6 and 9 km s −1 , respectively (see Table A .3).

Molecule Transition Frequency
To determine the water deuteration for the two velocity components, we used the HDO 225.9 GHz and H 2 18 O 203.4 GHz lines.We focused our analysis on the spectra extracted at the position of the N source for the component at 9 km s −1 and the spectra extracted at the position M (between N and S) for the component at 6 km s −1 .The line fluxes of the two components were collected by fitting the HDO and H 2 18 O lines with double Gaussian profiles with the Levenberg-Marquardt fitter (see Fig. A.1).The results of the fits are presented in Table A.2. Adding a third Gaussian for the CH 3 OCHO peak does not affect the HDO flux within the uncertainties.The Gaussian fluxes also do not differ if we fix the same FWHM and v LSR for the HDO and H 2 18 O lines.Based on the Gaussian fits, we calculated the column densities of each isotopolog assuming local thermodynamic equilibrium (LTE), following the formalism of Goldsmith et al. (1999) and Mangum & Shirley (2015).Then we derived the HDO/H 2 O ratios assuming a 16 O/ 18 O ratio equal to 560 (Wilson & Rood 1994).We considered a range of excitation temperatures, T ex , between 100 and 300 K, which are characteristic of hot corinos (e.g., Jørgensen et al. 2018).Source sizes of 0.35 and 0.45 were assumed for M and N, respectively.The lines are optically thin (τ ≤ 0.23 for HDO and τ ≤ 0.15 for H 2 18 O) in both cases.C.1.The equivalent D/H ratio is indicated on the right axis.For L1551 IRS5 N (v LSR ∼ 9 km s −1 ), the ratio is plotted for an excitation temperature of 300 K. tion of the excitation temperature.The extreme values found for the HDO/H 2 O ratio are listed in Table 2.The ratio varies between (1.6 ± 0.6) × 10 −3 and (2.1 ± 0.8) × 10 −3 for the component at 9 km s −1 .Only a lower limit of >0.3 × 10 −3 can be derived for the component at 6 km s −1 given the blending of the H 2 18 O line with CH 3 OCH 3 .Using a different size does not significantly affect the HDO/H 2 O ratio (see Appendix B.1). Non-LTE analysis also leads to similar results (see Appendix B.2).
As the fluxes of the two components at 6 and 9 km s −1 are not disentangled in the HDO observation at 143.7 GHz, we could not use it to derive the HDO/H 2 O ratios.However, it can help us constrain the excitation temperatures of HDO.The analysis described in Appendix B.3 shows that a model with an excitation temperature for N equal to or lower than 220 K cannot reproduce the flux of the 143.7 GHz line within 20%.The HDO/H 2 O ratio derived for T ex = 300 K should consequently be considered more reliable than the value obtained for T ex = 100 K.

Discussion
Figure 3 shows a comparison of the HDO/H 2 O ratio of L1551 IRS5 with the values previously found in protostars and comets.For the component at 9 km s −1 , the HDO/H 2 O ratio is similar to those measured in isolated Class 0 protostars (Jensen et al. 2019) and in the disk of the Class I V883 Ori (Tobin et al. 2023), and significantly higher than the previously studied clustered sources.At 6 km s −1 , the lower limit does not provide sufficient constraints to discriminate between the two categories.
Water observed in the warm inner regions of protostellar objects results from the thermal desorption of the icy grain mantles.Deuteration takes place early in the star formation process when the temperature is low and the density is relatively high.As a result, the gas phase HDO/H 2 O ratio observed in the warm regions of Class 0 objects should reflect those of the icy grain mantles.Consequently, the observed deuterium fractionation ratio at the Class I stage should either be similar to the Class 0 stage or decrease if water is partially reprocessed in the hot corino.In the latter case the deuteration process would not be efficient enough to increase or maintain the overall deuteration level.From the HDO/H 2 O ratio found in the high-velocity component of L1551 IRS5 we can conclude that the HDO/H 2 O ratio does not decrease and is still quite high, similarly to what is found for V883 Ori (Tobin et al. 2023) and SVS13-A (Codella et al. 2016).

L17, page 4 of 11
A second conclusion is that the chemistry of L1551 IRS5 appears to be closer to isolated sources than to clustered sources.The similarity with isolated sources can be surprising at first sight, as L1551 IRS5 is part of the Taurus molecular cloud.However, it is located in the LDN 1551 group in the southern part of the cloud (e.g., Roccatagliata et al. 2020).This group shows a lower density of sources than other regions in Taurus (Gomez et al. 1993).The spatial density of sources in the Taurus molecular cloud is also much lower (only a few tens in ∼1 pc 3 , Gomez et al. 1993) than in Ophiuchus, which contains IRAS 16293−2422, and the NGC 1333 cloud in Perseus (a few 10 2 -10 3 , Bontemps et al. 2001;Jørgensen et al. 2006).As explained in Jensen et al. (2021), the higher water deuteration in isolated sources would be due to either a lower temperature or a longer prestellar phase compared to clustered sources.The environment around L1551 IRS5 could be colder than around IRAS 16293 and the NGC 1333 sources.The dense cores could also have more time to evolve before their collapse if L1551 IRS5 is on the edge of the cloud and if it is surrounded by fewer sources.As a consequence, the water deuteration in L1551 IRS5 would be closer to the Class 0 isolated sources.A similar explanation was proposed for the Class I protostar V883 Ori, which belongs to the Orion molecular cloud, but is in a relatively isolated area (Tobin et al. 2023).The lower values found in comets and in the protostars IRAS 16293−2422, NGC 1333 IRAS2A, IRAS4A, and IRAS4B suggest that our Solar System formed in a relatively dense clustered area.The similarity of the HDO/H 2 O ratio between the binary L1551 IRS5 and the singly Class I object V883 Ori also seems to indicate that the binarity does not influence the water deuteration.
In conclusion, this study highlights the importance of carrying out water deuteration measurements toward more sources, at different evolutionary stages, and in different environments in order to better understand the water formation and evolution during the formation process of solar-type stars.Investigations of the D 2 O/HDO ratio in Class I protostars would also be useful to check the similarity with Class 0 sources.

Appendix A: Additional tables and figures on observations
Table A.1 lists the calibrators used for the different observations.Gaussian fitting was performed on the spectra to extract the fluxes of the HDO and H 2 18 O lines (see Figure A.1).A summary of the fit results is presented in Table A.2.We also performed circular Gaussian fitting of the emission lines in the (u,v)plane to estimate the size of the sources.The results of the fits are summarized in Table A.3.Notes. (a) The two components at 6 and 9 km s −1 cannot be distinguished because of the poor spectral resolution of the data (1.95MHz).The total integrated flux is consequently only given for this transition.A similar flux is obtained for the N and M positions because of the low spatial resolution of the data (2.5 × 2.1 ).L17, page 7 of 11 L1551 IRS5 using the non-LTE radiative transfer code RADEX (van der Tak et al. 2007) with the Large Velocity Gradient (LVG) expanding sphere method.We varied the temperature between 100 and 300 K, similarly to the LTE analysis.For the H 2 density, we considered a range from 10 7 to 10 11 cm −3 .With a non-LTE analysis of three HDO lines, Coutens et al. (2014)

Fig. 1 .
Fig.1.Comparison of the spectra of the water isotopolog transitions.The green, blue, red, and orange solid lines represent the spectra extracted at the position S, M, N, and O, respectively.The dotted lines show velocities of 6.0 km s −1 (blue) and 9.5 km s −1 (red).Considering the synthetic beam sizes and the separation of the N and S components, the extracted spectra are not entirely independent from each other.For the HDO line at 143.7 GHz, the spatial resolution is too low to distinguish O from N and M from S. The x-axis scale is larger for this transition.
Figure B.2 shows the variation in the HDO/H 2 O ratio as a func-

Fig. 2 .
Fig. 2. Integrated emission maps of the HDO line at 225.9 GHz (top) and the H 2 18 O line at 203.4 GHz (bottom).The blue and red contours correspond to the integrated flux between 4.5 and 7.5 km s −1 and 7.5 and 12 km s −1 , respectively.The levels of contours start at 3σ with a step of 3σ for all except for the HDO contours in red (7.5-12 km s −1 ) which are spaced by 5σ.Dust continuum is indicated in grayscale.The beam sizes are shown in the bottom left corner.

Component 9 5 )Fig. 3 .
Fig. 3. Comparison of HDO/H 2 O ratios in comets and protostars.The values and references can be found in TableC.1.The equivalent D/H ratio is indicated on the right axis.For L1551 IRS5 N (v LSR ∼ 9 km s −1 ), the ratio is plotted for an excitation temperature of 300 K.
Fig. A.1.Gaussian line fits using the CASSIS software.The upper panel shows the spectra of the HDO transition at 225.9 GHz, while the lower panel is for the H 2 18 O transition at 203.4 GHz.The left plots are for the M position with the fit of the component at 6 km s −1 (blue), while the right plots are for the N position with the fit of the component at 9 km s −1 (red).The dotted lines indicate the central velocity of the respective Gaussian.The green line shows the corresponding double-Gaussian profiles (see Table A.2).
derived a H 2 density higher than 10 8 cm −3 in the warm inner regions of the low-mass protostar NGC1333 IRAS2A.The results are shown in Figure B.3.The HDO/H 2 O ratio derived with RADEX, (1.3 -1.8) × 10 −3 is in agreement with the LTE values within the uncertainties.B.3.Modeling of the HDO line at 143.7 GHz and constraints on the excitation temperatures The spectral resolution of the HDO transition at 143.7 GHz does not allow us to distinguish between the fluxes coming from the component at 6 km s −1 and the one at 9 km s −1 .However, we can use the total flux to constrain the excitation temperatures for the two components.Assuming the source sizes obtained with the fits in the (u,v)-plane, we calculated the column densities of HDO needed to reproduce the two components at 225.9 GHz for various excitation temperatures.Then we predicted the expected flux at 143.7 GHz at 6 and 9 km s −1 , and summed the two fluxes before comparing them to the observed flux (see Table A.2). Figure B.4 shows that a good agreement (< 20%) is only found if the excitation temperature for the component at 9 km s −1 is 240 K.The excitation temperature is not constrained for the component at 6 km s −1 .

FigFig
Fig. B.1.HDO/H 2 O ratios (× 10 −3 ) obtained for the component at 9 km s −1 (N) of L1551 IRS5 as a function of different excitation temperatures and sizes assuming optically thin LTE emission.
Fig. B.4. Percentage errors between the predicted and observed fluxes of the HDO line at 143.7 GHz as a function of excitation temperature for positions M and N. Source sizes of 0.35 and 0.45 are assumed for the M and N positions, respectively.

Table 2 .
Column densities and HDO/H 2 O ratio in L1551 IRS5 for different excitation temperatures.

Table A . 3 .
Results from the circular Gaussian fits performed in the (u,v)-plane.