H 2 mass-velocity relationship from 3D numerical simulations of jet-driven molecular outflows
Adriano Cerqueira, Bertrand Lefloch, Alejandro Esquivel, Pedro Rivera-Ortiz, Claudio Codella, Cecilia Ceccarelli, Linda Podio
AAstronomy & Astrophysics manuscript no. ms39269_last © ESO 2020December 15, 2020 H mass–velocity relationship from 3D numericalsimulations of jet-driven molecular outflows A.H. Cerqueira , B. Lefloch , A. Esquivel , P. R. Rivera-Ortiz , C. Codella , , C. Ceccarelli , and L. Podio LATO / DCET, Universidade Estadual de Santa Cruz, Rod. Jorge Amado, km 16, Ilhéus, BA, CEP 45662-900, Brazile-mail: [email protected] CNRS, IPAG, F-38000 Grenoble, FranceUniv. Grenoble Alpes, IPAG, BP 53, F-38041 Grenoble, Francee-mail: [email protected] Instituto de Ciencias Nucleares, Universidad Nacional Autónoma de México,Apartado Postal 70-543, 04510 Ciudad de México, México INAF, Osservatorio Astrofisico di Arcetri, Largo E. Fermi 5, 50125 Firenze, ItalyReceived: ______ ; accepted: ______
ABSTRACT
Context.
Previous numerical studies have shown that in protostellar outflows, the outflowing gas mass per unit velocity, or mass–velocity distribution m ( v ), can be well described by a broken power law ∝ v − γ . On the other hand, recent observations of a sample ofoutflows at various stages of evolution show that the CO intensity–velocity distribution, closely related to m ( v ), follows an exponentiallaw ∝ exp( − v / v ). Aims.
In the present work, we revisit the physical origin of the mass–velocity relationship m ( v ) in jet-driven protostellar outflows. Weinvestigate the respective contributions of the di ff erent regions of the outflow, from the swept-up ambient gas to the jet. Methods.
We performed 3D numerical simulations of a protostellar jet propagating into a molecular cloud using the hydrodynamicalcode Yguazú-a. The code takes into account the most abundant atomic and ionic species and was modified to include the H gasheating and cooling. Results.
We find that by excluding the jet contribution, m ( v ) is satisfyingly fitted with a single exponential law, with v well in therange of observational values. The jet contribution results in additional components in the mass–velocity relationship. This empiricalmass–velocity relationship is found to be valid locally in the outflow. The exponent v is almost constant in time and for a givenlevel of mixing between the ambient medium and the jet material. In general, v displays only a weak spatial dependence. A simplemodeling of the L1157 outflow successfully reproduces the various components of the observed CO intensity–velocity relationship.Our simulations indicate that these components trace the outflow cavity of swept-up gas and the material entrained along the jet,respectively. Conclusions.
The CO intensity–velocity exponential law is naturally explained by the jet-driven outflow model. The entrained mate-rial plays an important role in shaping the mass–velocity profile.
Key words.
Stars: formation – ISM: jets and outflows
1. Introduction
Outflows from young stellar objects (YSOs) can exhibit a greatvariety of morphological and physical characteristics. In theyoungest (10 yr) and deeply embedded Class 0 protostars (An-dré, Ward-Thompson & Barsony 2000), outflows are easilytraced using the CO molecule and their presence is ubiquitousin star forming regions, indicating that they are a common man-ifestation of both low- and high-mass star formation processes(Wu et al. 2004; Lee 2020). On the other hand, protostellar jetswere first associated with more evolved Class II objects, that is,optically revealed pre-main sequence objects that are still accret-ing (or classical T-Tauri stars). These jets are observed mainlythrough forbidden atomic emission lines, like [S II] and [N II],as well in H α (Reipurth & Bally 2001). In between these twolimiting cases, Class I protostars, with a typical age of 10 , mayshow evidence for both molecular outflows and protostellar jets at the same time (e.g., L1448 IRS 2 and IRS 3, see Bally, Devine& Alten 1997). Sometimes a fast and collimated molecular jetis also observed, as in Cep E-mm (Lefloch et al. 2015) or HH212, which are associated with a Class 0 source (Zinnecker, Mc-Caughrean & Rayner 1998; Lee et al. 2017).Nevertheless, the origin of the molecular outflows associ-ated with YSOs is still debated (see Lee 2020, for a recent re-view). They are believed to be the by-product of an interactionbetween a more collimated jet and / or wind, produced in or by thestar–disk interaction, and its surrounding medium (Bally 2016).As the wind and / or jet bow shock propagate into the ambientmedium, the gas of the excavated cavity walls advances andexcites a profusion of molecular emission lines. For low massYSOs in particular, this can take place via one of three mainmechanisms (see Arce et al. 2007, for a comprehensive review):(i) In wind-driven shell models , a wide-angle wind is supposedto accelerate the ambient medium gas. In this class of models, Article number, page 1 of 15 a r X i v : . [ a s t r o - ph . GA ] D ec & A proofs: manuscript no. ms39269_last both the wind and the surrounding medium are assumed to bestratified in density. (ii) In turbulent jet models, a jet subject todynamical Kelvin-Helmholtz instability can entrain gas througha growing turbulent layer, giving rise to an outflow. This mech-anism can also operate at the leading working surface. (iii) In
Jet bow-shock models, a collimated jet produces a leading bowshock that accelerates the ambient medium gas. Also, an inter-mittent jet may develop a set of internal working surfaces thatcan help in the process (Raga & Cabrit 1993).As emphasized in Arce et al. (2007), it is possible that morethan one mechanism is operating to produce a given molecularoutflow, or alternatively a given mechanism can dominate at dif-ferent epochs in the evolution of a given source. In any case, a pa-rameter that has been historically used to identify useful modelsis the slope of a power law that relates the mass of the outflowingmolecular gas with its velocity, or m ( v ) ∝ v − γ . Rigorously speak-ing, the mass–velocity relationship, sometimes called the massspectrum, is obtained by considering the mass in a given radialvelocity bin, meaning that the observed relationship is actually δ m ( v ) /δ v ∝ v − γ (see Arce & Goodman 2001, for a discussion).However, for the sake of simplicity, some authors refer to themass–velocity relationship as m ( v ) ∝ v − γ (Downes & Ray 1999;Downes & Cabrit 2003). The mass–velocity relationship givesus the mass of the outflow at a given radial velocity. We note thatwhat is actually observed is the intensity of a given emissionline, typically the CO (1-0) line profile, and that such an inten-sity correlates with the velocity as described above (Downes &Ray 1999; Arce & Goodman 2001). The intensity is then con-verted to mass (corrected or not by the opacity) to finally ob-tain the mass–velocity relationship. Molecular outflows seem todisplay a mass–velocity relationship that can be described by abroken power law (Bachiller & Tafalla 1999; Ridge & Moore2001; Arce & Goodman 2001), with shallower slopes ( γ < v <
10 km s − ) and steeper slopes ( γ >
3) atintermediate-to-high velocities ( v >
10 km s − ). In this way, nomatter the mechanism used to model a molecular outflow, themodel should account for the observed slopes.In the present paper, we focus on the molecular mass–velocity relationship for molecular outflows produced by a col-limated and supersonic jet using three-dimensional numericalsimulations. As the jet interacts with the ambient medium, jet-entrained gas and ambient gas swept up by the jet-driven bow-shock can in principle be disentangled. The appeal of such ascenario is two-fold: (i) there is increasing evidence that bothphenomena may coexist in Class 0 and Class I sources, as men-tioned in the previous paragraph, and (ii) molecular outflowsproduced by either jet entrainment or a jet-driven bow-shockcan e ff ectively end up in a power-law mass–velocity relationship(Chernin & Mason 1993; Zhang & Zheng 1997; Stahler & Palla2004). In the following section, we briefly compile some previ-ously important results obtained through numerical simulationsof jet-driven molecular outflows, and discuss some recent obser-vational findings that ultimately motivated the present work.
2. Previous numerical results and observedintensity–velocity relationships
Numerical simulations of molecular jets have been used exten-sively in the literature as an e ffi cient tool to investigate the kine-matical properties of jet-driven molecular outflows (see Downes& Ray 1999; Downes & Cabrit 2003; Rosen & Smith 2003,2004a; Smith & Rosen 2005, 2007). The mass–velocity andintensity–velocity relationships I CO ( v ) observed in low-J ( ≤
8) CO lines in molecular outflows have been studied by various au-thors as a possible test for discriminating between entrainmentmechanisms (see Downes & Cabrit 2003, hereafter, DC03).Previous works have described the CO intensity–velocitydistribution I CO ( v ) in outflows as a broken power law, I CO ( v ) ∝ v − γ with γ (cid:39) . v break ≈
10 – 30km s − and a steeper slope γ = ∼
20 km s − and of the temperaturedependence of the line emissivity (see DC03).Downes & Ray (1999) introduced the H molecule in theircalculations and found that the H mass–velocity distribution m ( v ) ∝ v − γ follows a similar relationship to the intensity–velocity relationship observed in the millimeter rotational linesof CO. Downes & Cabrit (2003) showed that the swept-upmolecular gas follows a mass–velocity relationship ( m H ( v ))similar to the intensity–velocity relationship I H ( v ) in the lowvelocity range ( v (cid:46)
30 km s − ). In contrast, in the high ve-locity range these latter authors found that I H ( v ) is shallowerthan m H ( v ), while I CO ( v ) is steeper than I H ( v ) but comparableto m H ( v ).Rosen & Smith (2003) focused on time-dependent jets, thatis, jets whose density varies as a function of time with respect tothe ambient medium. These latter authors found that the mass–velocity distribution is systematically shallower than the COintensity–velocity distribution. They also found that the indicesof the distributions are essentially unchanged when consideringatomic or molecular jet material.Keegan & Downes (2005) studied the temporal evolution ofthe power index γ, and their results are consistent with thoseof Downes & Cabrit (2003). Interestingly, Keegan and Downesfound that γ should increase slowly in time, attaining a limitingvalue after t ≈ I CO ( v ) ∝ exp( − v / v ), with v ∼ −
12 km s − , cameas a surprise (Lefloch et al. 2012). Further observational studiesconfirmed that these spectral signatures, with similar values of v , were detected in a plethora of molecular gas tracers, like CS(Gómez-Ruiz et al. 2015), HNCO and NH CHO (Mendoza etal. 2014), HC N (Mendoza et al. 2018), and HCO + (Podio et al.2014).The same analysis applied to the outflow sample of Bachiller& Tafalla (1999) observed in the CO J = I CO ( v ) ∝ exp( − v / v ), with values of v between 2and 12 km s − , well in the range of those determined in L1157-B1. Also, Lefloch et al. observed the trend that the more evolved,Class I outflows of the sample (Mon R2, L1551) display a shal-lower intensity–velocity distribution. Therefore, an exponentialrelation I CO ( v ) ∝ exp( − v / v ) is found to be a good approxima-tion of the observed intensity relation not only in L1157-B1 but Article number, page 2 of 15erqueira et al.: Jet-cavity kinematical relationship in several molecular outflows in general, with a reduced numberof free parameters compared to a broken power law. Mappingof the CO J = / surrounding protostarenvelope. Specifically, we present a new methodology based ondistinguishing the mixing level between the ambient mediumand the jet gas, which helped us to disentangle the distinct com-ponents that arise in the H mass–velocity relationship. The arti-cle is organized as follows. In §3.1 and §3.2 we provide details ofthe numerical setup and initial parameters for the simulations. In§3.4 we briefly compare our results with some previous numeri-cal studies of molecular jets. In §4 and 5 we present the results ofour numerical simulations and we discuss their implications forobservations of molecular outflows using L1157 as a reference.In §6 we present our conclusions.
3. Numerical simulations
The simulations presented here were performed using theYguazú-a code (Raga, Navarro-González & Villagrán-Muniz2000; Raga et al. 2002; Cerqueira et al. 2006). In its originalversion, the code was designed to solve hydrodynamic prob-lems with a chemical network for the following atomic and ionicspecies: HI, HII, HeI, HeII, HeIII, CII, CIII, CIV, NI, NII, NIII,OI, OII, OIII, OIV, SII and SIII. For the present work, we intro-duced the H molecule as a new species and added three disso-ciation reactions for molecular hydrogen:H + H → , (1)H + H → + H , (2)e + H → e + . (3)We used the collisional dissociation rates of molecular hy-drogen provided in Shapiro & Kang (1987) for these three reac-tions. We also calculated the cooling function considering boththe radiative and the dissociative processes. For the radiativecooling rate, we used the fit proposed by Lepp & Shull (1983),which considers both the rotational and vibrational cooling fromthe two reactions, H-H and H -H , in both high- and low-density regimes ( n < n cr ≈ cm − ). The dissociative cool-ing function was taken from Shapiro & Kang (1987). In Fig. 1we show the di ff erent cooling functions: atomic (blue line) andmolecular dissociative (green line) and radiative cooling (redline) . Fig. 1. Di ff erent contributions for the cooling: atomic emission lines(blue line), H dissociative cooling (green line), and H radiative cool-ing (red line). The cooling functions were calculated using the startingvalues (i.e., at t =
0) for the numerical densities, or n HI = − , n HII = .
08 cm − , n He = .
94 cm − and n H = .
51 cm − . Together with H , CO and H O have long been known toplay an important role in shocked gas cooling (Hollenbach &McKee 1979; Kaufman & Neufeld 1996; Flower & Pineau desForêts 2010). Detailed observational studies have confirmed thatline cooling from CO and H O can be as important as that fromH in protostellar outflows (see e.g. Nisini et al. 2010a; Busquetet al. 2014). Modeling of the structure of outflow shocks may besignificantly modified by the inclusion of additional terms suchas CO, H O, or even charged grains (see e.g. Flower & Pineaudes Forêts 2010), all of which are not taken into account here.In the present work, we have not included either the H O or theCO chemical networks or their related cooling terms. This wouldrepresent an e ff ort which is well beyond the state of the art of 2Dand 3D chemo-hydrodynamical codes such as WALKIMYA-2D(Castellanos-Ramírez et al. 2018, Rivera-Ortiz et al., in prep.).However, we note that Rosen & Smith (2004a) included equi-librium C and O chemistry in their numerical scheme in orderto calculate the CO, OH, and H O abundances, as well as toestimate the cooling expected from these molecules. They con-cluded that the mass–velocity relationship is always shallowerthan the intensity–velocity relationship, confirming previous re-sults based only on H (Downes & Cabrit 2003) . For that rea-son, the present work focuses on the entrained gas propertiesand we aim to revisit the H mass–velocity relationships, whoseproperties can be accessed following the present prescription.Our computational domain is a Cartesian 3D rectangular boxwith the following dimensions:( x , y , z ) = (2 , , × au , (4)and the jet propagates along the z − direction.The Yguazú-a is a multi-level binary adaptive grid code.Here we use a five-level grid which has ( x , y , z ) = (256, 256, In order to calculate each one of these curves we considered the ini-tial values for the numerical densities for the di ff erent species (atomic,ionic and molecular). The CO intensity–velocity is calculated implicitly in Downes &Cabrit (2003) using an analytical prescription and the local density, as-suming that CO density is 10 − of the H density.Article number, page 3 of 15 & A proofs: manuscript no. ms39269_last
Table 1.
Jet models
Model v j , A n a n j η T a T j τ P θ N tot τ e (km s − ) (cm − ) (cm − ) (K) (K) (years) ( ◦ ) (years)DR_SS 212 0 100 100 1 100 1000 - - - -DR 212 0.15 100 100 1 100 1000 - - 4 5, 10, 20 and 50DR_P 212 0.15 100 100 1 100 1000 200 6 4 5, 10, 20 and 50 Notes. v j , is the jet velocity in km s − ; A is the amplitude of variation in the jet velocity; n a and n j are the ambient medium and jet (numericalparticle) density; η = n j / n a is the jet to ambient medium density ratio; T a and T j are the ambient medium and jet temperatures (in K); τ P is theprecession period, N tot the number of di ff erent jet injection period, and τ e is the jet velocity variability period. ∆ x = ∆ y = ∆ z = .
13 au. The jet radius is ini-tially always given by R j =
391 au or ∼ ∆ x . The jet radiusis therefore compatible with those adopted in previous numeri-cal simulations (Downes & Ray 1999; Downes & Cabrit 2003)as well as with estimates for the HH jet radius (Reipurth et al.2000, 2002; Podio et al. 2006). Three cases were considered, which are summarized in Table 1: – model DR: an intermittent jet model for comparison withpreviously published simulations in the literature (DC03,Downes & Ray 1999); – model DR_SS: a steady state jet; – model DR_P: an intermittent, precessing jet model.With model DR_P, we aim to investigate the properties of theL1157 outflow, kinematical studies of which have revealed con-vincing evidence of precession (Gueth, Guilloteau & Bachiller1996; Podio et al. 2016). In this simulation (and for model DR)we assume that the jet velocity varies periodically with time, ac-cording to: v j = v j , · (cid:20) + A N tot (cid:88) i = sin (cid:18) πτ e , i · t (cid:19)(cid:21) , (5)where v j , is the jet velocity, A = ∆ v / v j , is the adopted amplitudevariation for the jet velocity variation, and τ e is the variabilityperiod ( t is the time). In our time-varying models, N tot = (cid:46) τ e , i (cid:46)
50 years (see Table 1). We note that although Eq.5 has been used here in an attempt to reproduce the model pre-sented in DC03, the idea that Herbig-Haro objects can be generi-cally explained by successive internal knots promoted by a sinu-soidal jet velocity variability is well established (e.g., Reipurth& Bally 2001). However, a detailed source modeling can requirea superposition of di ff erent sinusoidal terms, which have beendiscussed by Castellanos-Ramírez, Raga & Rodríguez-González(2018; see also Bally 2016), indicating that a multimode jet ve-locity variability may be important to explain the observed mor-phology and kinematics in some sources. For the precessing caseDR_P, we adopted a precessing angle of θ = ◦ and a precessingperiod of τ P =
200 years. The DR model has the same parame-ters as the model presented in DC03.In all models, we assume solar elemental abundances forboth ambient and jet material. The ratio n H / n H = n H = n HI + n HII ) is initially imposed for both the jetand the ambient medium (Downes & Ray 1999; Nisini et al.2010b). The helium fraction per hydrogen nuclei is assumed to be n He / [ n H + n H ] = . µ = . C V = .
25. The ionization fraction of hydrogen in the jet isinitially taken as f H = .
01 for T j = K in agreement withthe values inferred from atomic line observations of HH jets inthe optical (Podio et al. 2006).The ambient medium and jet parameters such as numericaldensity n , temperature T , and jet velocity are all given in Table1, along with the jet precessional and intermittence periods ofthe simulated models.Observationally, the jet temperature and density determina-tions span a wide range of values depending on the tracer used.Optical atomic line observations yield T j ∼ × − × K (Podio et al. 2006) while molecular line observations indicatelower values of about 10 − × K from near-infrared H rovibrational transitions (e.g., Caratti o Garatti et al. 2006), and T j ∼ n j ∼ − cm − , while(sub)millimeter line observations yield high values, n j (cid:38) cm − . This wide range of physical conditions reflects the intrin-sic complexity of the jet, which is often associated with inter-nal shocks that drive the formation of strong temperature anddensity gradients. Adopting single initial values for temperatureand density is most likely an oversimplified description of the jetphysical structure. We note however that the initial jet tempera-ture value adopted in the simulations are consistent with thoseobtained from jet molecular line observations (H , SiO, CO).The initial jet density in our simulations is 100 cm − (same as inDC03), which is lower than the values determined observation-ally. However, it is the jet-to-ambient density ratio which carriesthe most weight in modeling the dynamical evolution of the out-flow. This point was investigated in detail by Rosen & Smith(2004b). Based on their results, we do not expect significant dif-ferences in the simulations when adopting a higher density forthe jet, provided that the jet-to-ambient density ratio is kept con-stant. Our primary diagnostics are the mass–velocity relationship forboth the total mass m ( v ) and the molecular mass, m H ( v ), whichare computed for the whole computational domain or for a givenspatial region. In order to obtain the mass–velocity profiles, wefirst compute the column density N as the sum of particle densityalong the line of sight per velocity interval: N ( y (cid:48) , z (cid:48) ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) v = v CM = (cid:88) n (cid:48) ∆ x (cid:48) , (6) Article number, page 4 of 15erqueira et al.: Jet-cavity kinematical relationship
Fig. 2.
Sketch of the geometry of the flow with respect to the observer.The jet propagates along the z -axis, which is inclined by an angle i withrespect to the plane of the sky ( y (cid:48) − z (cid:48) ). In case of precession, the jetdraws a cone with a half angle of θ with respect to the z -axis. where v = v x cos i − v z sin i . (7)In Equation (6), ∆ x (cid:48) is the projection of the x − coordinate alongthe line of sight and n (cid:48) is the numerical particle density (total ormolecular) in the radial velocity range ( v − ∆ v / < v < ( v +∆ v / ∆ v = − . In Equation (7), i is theinclination angle with respect to the plane of the sky (see Fig. 2),meaning that v corresponds to the (observed) radial velocity.As the jet propagates, interaction with the ambient gas leadsto the formation of a mixed gas layer from jet and ambient ma-terial. In order to estimate the relative contribution of the di ff er-ent regions of the outflow, namely the swept-up ambient mediumand the mixed (jet plus ambient medium) material, we tagged thejet material with a passive scalar or tracer j . This passive scalaris set to j = j = j mix to indicate the levelof mixing considered, such that j ≤ j mix . Thus, material with j mix = j mix = . .
5, and that with j mix = We tested our numerical scheme against previous numerical sim-ulations published in the literature by running model DR (Table1). The parameters of model DR were chosen to mimic modelG in Downes & Ray (1999), as well as the “pulsed” model dis-cussed in Downes & Cabrit (2003).Briefly, the DR model has jet ( n j ) and ambient medium ( n a )numerical particle density n j = n a =
100 cm − , jet temperature T j = T a =
100 K, molec-ular hydrogen to hydrogen numerical density ratio n H / n H = , and a helium to hydrogen numerical density ratio of n He / n H = .
1. The only H dissociation process considered in this model is collisions with H atoms, and the rate coe ffi cient used here isgiven by: k D , H − H = × − exp( −
55 000 / T ) cm s − , (8)following Taylor & Raga (1995).Figure 3 shows the results of model DR at t =
400 yr. De-picted in this figure are the midplane ( x =
0) distribution of thetracer (panel a ), the temperature (panel b ), the velocity compo-nents along the z − axis (panel c ) and the y − axis (panel d ), thetotal particle density n (panel e ), and the molecular particle den-sity n H (panel f ). We plot a white contour line that separates theoriginal (quiescent and / or disturbed) ambient medium where thetracer j mix is zero (the region outside the contour line) from theinner jet where j mix = e ) and 3( f ) show the internal working surfaces at z = e and f ). Although densityenhancement is observed in the internal shocks (see Fig. 3) as ex-pected, neither the total density nor the molecular density attaintheir maximum values at the internal working surfaces. The re-gion where the total density is maximum is located at the head ofthe jet and inside the contour line, indicating that the jet and theambient medium have already been mixed. It is important to notethat while the molecular density is higher at the external edges ofthe laterally expanding trails of the leading bow shock, the totaldensity peaks are near to the jet axis and the apex of the leadingbow shock. This occurs because strong shocks at the jet head dis-sociate the molecular hydrogen, while the shocks are weaker atthe bow shock trails and the molecular hydrogen piles up withoutbeing dissociated. We already anticipate that the mass–velocityprofiles extracted from our simulations should present a high-velocity component ( v ∼
100 km s − ) if some mixing betweenthe jet and the ambient medium material is allowed. In this section, we compare the results of model DR (see Table1) with the results of Downes & Ray (1999). More precisely,we compare the results of the mass–velocity distribution at t =
300 yr and for an inclination angle of i = ◦ with respect to theobserver for a direct comparison with Fig. 3 in Downes & Ray(1999). In Fig. 4 we report the molecular mass–velocity profileobtained when integrating over a region defined by a minimumlevel of mixing between the jet and ambient medium materialfrom j mix = j mix = . m ( v ) ∝ v − γ in red. For the sake of clar-ity, we have made adjustments for two distinct radial-velocityintervals: 0 < v <
10 km s − and 10 < v <
100 km s − . Here-after, v is used to refer to the radial velocity.The mass–velocity relations in Fig. 4 are well described bya broken power law. At v ∼
10 km s − , the slope changes in allcases from j mix = j mix = .
0. As expected, the slopein the low-velocity range is always shallower than in the high-velocity range. However, there is an important and systematice ff ect in the slopes caused by the jet material removal from theintegration process that results in an overall steepening in themass–velocity relation.We can see in Fig. 4 that considering di ff erent levels of jetcontent in the computation of the mass–velocity relation inducesa similar e ff ect: in the high-velocity range ( v >
10 km s − ), the Article number, page 5 of 15 & A proofs: manuscript no. ms39269_last
Fig. 3.
Model DR at t =
400 yr. Distributions in the plane x = a ) j mix (tracer), ( b ) temperature, ( c ) v z , ( d ) v y , ( e ) total density n , and ( f )molecular gas density n H . The white line in each panel separates the region in the computational system filled by the ambient medium —wherethe tracer is equal to zero— from the mixing region and the jet itself. slopes are steeper for lower values of j mix . By contrast, in thelow-velocity range ( v ≤
10 km s − ) the slopes barely vary with j mix . We interpret these facts as a consequence of e ff ect of massaddition in a given velocity channel, when we go deeper into themixing layer. In the low-velocity range the profile is dominatedby the material of ambient origin, which can be either swept-upgas by the jet-driven bow-shock or entrained gas that barely in-teracted with the jet. The contribution of the jet and / or ambientinteracting material becomes increasingly apparent when con-sidering increasing j mix values in the high-velocity range.In their simulation run, Downes & Ray (1999) obtained amass–velocity distribution which was best fitted with a powerlaw of index γ = .
98 in the range 10 < v <
100 km s − (seetheir Fig. 3). In our simulation run with the same set of parame-ters with Yguazú-a, we obtain a similar mass–velocity relation-ship, which can be fitted by a power law (Fig. 4). However, wefind that the power-law index depends on the degree of mixingbetween the jet and the entrained ambient material, i.e., withthe value of j mix . For the swept-up (unmixed) molecular mate-rial ( j mix = γ = .
07, which is steeper than the valuefound by Downes & Ray (1999). However, if we take into ac-count the contribution of mixed ambient and jet material (e.g., j mix = .
6) to the mass–velocity relationship, we obtain a bestfit with a power law with a shallower index of γ = j mix = . j mix ≥ . Using a numerical code very similar to that of Downes & Ray(1999) and with similar initial conditions, DC03 further inves-tigated the mass–velocity and intensity–velocity relations in theCO J = S(1) v = molecular column density (leftpanels) and the mass–velocity relationships (right panels): themolecular m H ( v ) (black lines) and the total mass–velocity rela-tions m ( v ) (blue lines). For the sake of direct comparison with theresults presented by DC03 (see their Fig. 2), the molecular mass m H ( v ) and mass m ( v ) velocity relationships have been extractedfrom the DR model at t =
400 yr, considering an inclination an-gle of i = ◦ and di ff erent values of the jet and ambient gasmixing ratio j mix =
0, 0.3, 0.6, 1.0. In order to obtain the molec-
Article number, page 6 of 15erqueira et al.: Jet-cavity kinematical relationship
Fig. 4.
Model DR at t =
300 yr and an inclination angle of i = ◦ .Molecular mass–velocity relationship m H ( v ) (black) for the 0 km s − < v <
100 km s − radial velocity range and the best-fitting power law m ( v ) ∝ v − γ (red). The best-fitting index value γ is shown inside eachpanel for two intervals: 0 km s − < v <
10 km s − (left) and 10 km s − < v <
100 km s − (right). The j mix parameter is also indicated (from topto bottom, j mix =
0, 0.3, 0.6 and 1.0). The vertical axis displays thelogarithm of the molecular mass (in g). ular mass distribution of the outflow–jet system, we integratedover the full range of velocity, between + +
150 km s − ,and excluding the emission of the quiescent gas at rest velocity( v = ff erent values of j mix in Fig. 5 with respect to the swept-upgas (top panel) shows the presence of gas accelerated at v >
10 km s − already for j mix = . j mix = .
6, which results ina shallower mass–velocity distribution. This change of slope ismainly caused by the contribution of the massive gas clump thatdevelops near the apex of the bow shock inside the mixing layer.This e ff ect can be seen in Fig. 3(e).The total and molecular mass–velocity relationships are veryclose to each other in the swept-up gas (top panel in Fig. 5), Fig. 5.
Model DR at t =
400 yr and an inclination angle i = ◦ . Theresults are presented for four di ff erent values of the jet and ambientgas mixing degree j mix (from top to bottom): 0.0, 0.3, 0.6, 1.0. (left) Maps of H column density obtained by integration over the velocityrange, between + +
150 km s − . (right) Mass–velocity relationshipobtained for the molecular gas (black) and the total gas (blue). In thebottom panel ( j mix = . swept-up mass–velocity relationship obtained by DC03 (see their Fig. 2). Weshow the best-fitting power laws m ( v ) ∝ v − γ to the velocity intervals v ≤
10 km s − and 10 < v (cid:46)
20 km s − superposed in red dashed lines. Thepower-law index γ is given in each panel. The vertical arrows indicatethe position of bumps in the total mass–velocity profile (see text). which implies that the molecular gas dissociation can be ne-glected in the local gas acceleration (entrainment) process. Whentaking into account the jet–ambient gas mixing, the slope of thetotal mass–velocity relationship (blue in Fig. 5) starts to departfrom the molecular mass–velocity slope (black) for j mix = . ff erence increases with increasing velocity and in-creasing j mix (jet-ambient mixing ratio) values. The change ofslope between j mix = j mix = v ∼
50 km s − .We note that this velocity coincides with the projected veloc-ity of the high-density gas concentrated at the jet head. This re-gion can be seen in Fig. 3 (e) near z = . × au. The ve-locity component of this dense component along the jet axis is v z ≈
100 km s − , which corresponds to a projected (radial) ve-locity of 50 km s − . The increase of the mass at high velocities(in blue in Fig. 5) with increasing j mix values is essentially un- Article number, page 7 of 15 & A proofs: manuscript no. ms39269_last noticed in the molecular mass–velocity distribution (in black inFig. 5). This is consistent with e ffi cient H dissociation in theshocks at the jet head.In the case of full jet–ambient gas mixing ( j mix = v ≈
100 km s − (Fig. 5), which is the signature of the internalknots formed as a result of jet variability. This second bump ispresent in both the molecular (black lines) and the total mass(blue line) velocity distributions. Again, this is illustrated inFig. 3, where both total and molecular densities appear to be en-hanced behind internal shocks at a velocity consistent with theprojected velocities of the second bump. We emphasize that thissecond bump can only be seen if we consider the total jet mate-rial in the calculation ( j mix = . j mix .In the bottom right panel of Fig. 5 we have superimposedthe results of the Downes & Cabrit (2003) simulation with greenbullets. As we can see, the match between the swept-up H –velocity distribution of these latter authors and the results of oursimulations is very close. However, it should be noted that DC03claimed to have obtained the molecular mass–velocity relationfor the swept-up H , hence excluding the jet material from thecomputation (see their Fig. 2) . While DC03 claim that the jethas not been considered in their computation of the mass, wewere only able to reproduce their profile when taking into ac-count the jet contribution in the computation. The presence of ahigh-velocity bump in the H mass–velocity distribution is evi-dence that the jet has been considered in the integration proce-dure, as discussed above.To summarize, we carried out a detailed comparison of theYguazu-a code results with numerical simulations presented byDownes and Ray (1999) and DC03. We obtained excellent quan-titative agreement in both cases. This bolsters our confidence inthe use of this approach and in the results for the mass–velocitydistributions produced by our code.
4. Results
We first present the results of model DR_SS, which simulatesthe propagation of a nonprecessing jet under steady-state (SS)conditions. The parameters of the simulation are summarized inTable 1. Figure 6 shows the H column density spatial distri-bution, and the total and molecular mass–velocity distributionsobtained for di ff erent values of the jet-ambient gas mixing ratio j mix =
0, 0.3, 0.6, and 1.0 at t =
400 yr and for an inclinationangle of i = ◦ .The main di ff erence between model DR_SS and the simula-tion discussed above in Sect. 3 lies in the SS assumption, that is,the absence of intermittency in the mass-ejection process. Manysimilarities are therefore observed when comparing the mass–velocity relationships obtained in both simulations at a commonage of 400 yr, which are presented in Figs. 5-6. The main simi-larities can be summarized as follows: – The total mass and the molecular mass–velocity relations as-sociated with the swept-up gas ( j mix =
0) are very similar,and are separated by only a small and constant vertical o ff setover the whole velocity range ( (cid:46)
20 km s − ). – The gap between the total mass and the molecular mass be-comes more evident with increasing values of j mix and in-creasing velocity. For the case j mix = .
3, the curves end upat v (cid:39)
90 km s − showing a di ff erence in mass (total versusmolecular) of about two orders of magnitude. Fig. 6.
Model DR_SS at t =
400 yr and an inclination angle i = ◦ .The results are presented for four di ff erent values of the jet–ambientgas mixing degree j mix (from top to bottom): 0.0, 0.3, 0.6, 1.0. (left) Maps of H column density obtained by integration over the velocityrange between 1 and +
150 km s − . (right) Mass–velocity relationshipobtained for the molecular gas (black) and the total gas (blue). The best-fitting power laws for the three velocity intervals: v ≤
10 km s − (left;dashed red), 10 km s − < v <
30 km s − (middle; dashed green) and30 km s − < v <
60 km s − (right; dashed red) are drawn. The index γ for each one of these intervals are shown inside each panel. The verticalarrows indicate the position of bumps in the total mass–velocity profile(see text). – In the j mix = . v (cid:39)
50 km s − . The lack of detec-tion in the molecular mass–velocity distribution suggests itis mostly of atomic origin and that it traces the signature ofmaterial locally accumulated behind the jet head as a resultof shock compression. – A second bump of mainly atomic jet material is detected athigh velocity, namely v ∼
100 km s − (see panel j mix = . ff erences areseen between the two simulations when comparing the relativecontributions of the molecular and atomic material. Article number, page 8 of 15erqueira et al.: Jet-cavity kinematical relationship
The high-velocity bump detected at about 100 km s − is re-lated to the internal knots in the case of DR, whereas it traces thejet material in the SS model. Though mainly atomic, the bumpcontains a significant fraction of molecular material. The gapbetween the total (atomic + molecular) and the molecular mass iswider in the case of intermittent ejection (model DR; Fig. 5) thanin the SS regime (Fig. 6). Indeed, H dissociation is expected tooccur in the multiple internal shock knots produced in the pul-sating model; by comparison, in the SS model, only the leadingbow shock is expected to dissociate H . This is also well illus-trated by panels (e) and (f) in Fig. 3.The second bump detected at v ∼
100 km s − in the DRmodel (see Fig. 5) is fully developed. This indicates that thejet material reaches a higher velocity in the SS regime, ≥
100 km s − . As can be seen in Fig 6, in the SS regime, thejet reaches a velocity close to the maximum expected value, v = · sin30 ◦ =
106 km s − . On the contrary, in a pulsat-ing jet, the formation of internal shocks results in a decelerationof the jet material all along its axis.We determined the best-fitting power law to the molecularmass–velocity relationship for model DR_SS in the three veloc-ity intervals: v <
10 km s − , 10 km s − < v <
30 km s − and30 km s − < v <
60 km s − . The results are reported in Fig. 6,where we show the results of model DR_SS at t =
400 yr un-der an inclination i = ◦ with respect to the plane of the sky(see Fig. 2). The slopes of the molecular mass–velocity relation-ship decrease as j mix increases. In the low-velocity range, thebehavior of the relationship is very similar to that obtained for apulsating jet (model DR; Fig. 5). However, for the v >
10 km s − and j mix ≥ .
3, there is a region at intermediate radial velocities(10 km s − < v <
30 km s − ) with a moderate slope ( γ ∼ γ ∼ −
5) at high ve-locities ( v >
30 km s − ). Although a direct comparison with theslopes of the DR model (Figure 5) at intermediate-to-high veloc-ity is impossible because in that case the whole profile (from 10km s − < v <
60 km s − ) seems to be well described by a singlepower law, we can roughly estimate a mean γ value for DR_SSmodel using the two γ values obtained for v >
10 km s − , andwe obtain γ = .
92, 3.1, and 3.1 for j mix = .
3, 0.6, and 1.0, re-spectively. This suggests that the molecular mass–velocity pro-files at intermediate-to-high radial velocities for the SS modeltend to be shallower in comparison with those obtained for thepulsating model, indicating a high molecular mass content, a re-sult that is consistent with the scenario described in the previousparagraphs.
As many jets display hints of precession (e.g., Gueth, Guilloteau& Bachiller 1996, 1998; Podio et al. 2016; Santangelo et al.2015), we investigated the possible impact of this phenomenonon the mass–velocity relationship by running model DR_P (seeTable 1). We modeled precession with an angle θ = ◦ aroundthe jet axis and a period of τ =
200 years. The column densitymaps and the total and molecular mass–velocity distributions at t =
400 yr and i = ◦ , and di ff erent values of the jet–ambientgas mixing ratio j mix are presented in Fig. 7.The cavity created by the precessing jet appears broaderclose to the apex when compared with unprecesssing models,either intermittent (model DR; Fig. 5) or SS (model DR_SS;Fig. 6). Interestingly, both the molecular and total mass–velocitydistributions are very similar to those obtained in the case ofa pulsating, nonprecessing jet (Fig. 5) and a SS jet ( Fig. 6).The same two bumps at v ∼
50 km s − and v ∼
100 km s − are Fig. 7.
Model DR_P at t =
400 yr and an inclination angle i = ◦ .Results are presented for four di ff erent values of the jet–ambient gasmixing degree j mix (from top to bottom): 0.0, 0.3, 0.6, 1.0. (left) Mapsof H column density obtained by integration over the velocity rangebetween 1 and +
150 km s − . (right) Molecular mass–velocity relation-ship (black). The best-fitting power laws for the two velocity intervals: v ≤
10 km s − (left) and 10 < v (cid:46)
20 km s − are drawn in dashed red.The index γ is shown inside the panel. The exponential best fit is drawnin blue and the exponent v is given for each j mix value. also detected in the total mass–velocity distribution (not shownin Fig. 7). In summary, the mass–velocity relationship is not sig-nificantly altered by jet precession and is very similar to that ofnonprecessing, eventually pulsating jets. Above, we compare and analyze the results of our simulationsby modeling the mass–velocity distribution with a power law, m ( v ) ∝ v − γ . However, observational work on several molec-ular outflows by Lefloch et al. (2012) suggests that it couldbe possible to adopt another fitting, namely an exponential law m ( v ) ∝ exp( − v / v ). Based on our numerical simulations, we as-sessed the validity of this approach.We show the results of the fitting procedure for modelDR_SS (Fig. 8) and model DR_P (Fig. 9) at the di ff erent times t = , and, as in the preceding analyses, forfour values of the jet–ambient gas mixing ratio j mix = Article number, page 9 of 15 & A proofs: manuscript no. ms39269_last
Fig. 8.
Model DR_SS at t = i = ◦ . (Left) Distribution of moleculargas column density. (Right)
The mass–velocity relationship is depicted in black for each value of the j mix parameter, from j mix = . j mix = . v is given for each j mix value. i = ◦ . We firstconsider the H column density distribution and the molecularmass–velocity relationships resulting from the propagation of asteady-state jet (model DR_SS) . The value of the coe ffi cient v is indicated in each panel.Our first result is that it is indeed possible to obtain a verygood fit between the mass–velocity relationship and an expo-nential function m ( v ) ∝ exp( − v / v ) from early (400 yr) to late(2000 yr) computational times (see Figs. 8–9), although the shal-lower γ index in a power-law fitting translates into a higher v value. The value of v therefore increases with the jet–ambientgas mixing ratio. Our simulations for DR_SS (Fig 8) and DR_P(Fig. 9) show that higher values of v are found at earlier times(400 yr). Also, the best-fitting values of v for both modelsDR_SS (Fig. 8) and DR_P (Fig. 9) under the same inclinationangle are similar at the di ff erent times in the simulations, at 400,1200, and 2000 yr. In other words, regardless of the age of thesystem, once a minimum level of mixing occurs between the jetand the ambient medium gas ( j mix ≥ . v are almost unchanged. However, onecan observe a small decrease of v as time increases.Close inspection of the mass–velocity relationships revealsthat the exponential fitting (red curves) provides excellent solu-tions over the whole velocity range for j mix values of 0.3. Forhigher values of j mix , the high-velocity range of the distribution v ≥
10 km s − is still accurately fitted, with a higher value of v ,as can be seen in Fig. 9 (e.g., panels j mix = v (cid:38)
10 km s − in the case of swept-up H profiles (case j mix = .
0) and the presence of the jet as a bumpat v (cid:39)
40 km s − for j mix = .
0. Between these two extremeswe find great variability. We note that some residual emission isfound in the low-velocity (1 < v <
20 km s − ) range of the dis- tribution, which can be fitted by two exponential functions. Thispoint is addressed in more detail in the following paragraph. The results that we have discussed so far refer to the globalmass–velocity relationship computed over the whole computa-tional domain. However, the possibility of a local exponentialfitting was reported by Lefloch et al. (2012) thanks to CO multi-line observations of the shock position B1 in the southern lobe ofthe L1157 outflow. These latter authors noticed that similar spec-tral signatures in the CO J = t = i = ◦ and j mix = .
0. We selected four positions inthe outflow, labeled A1 to A4 and marked with black circles inthe map of molecular gas column density in Fig. 10. While po-sitions A1 and A2 are located inside the outflow cavity with A1close to the jet main axis, positions A3 and A4 are located atshocked positions at the interface between the outflow and theambient gas. We computed the molecular mass–velocity rela-tionships over circular areas of ∼ Article number, page 10 of 15erqueira et al.: Jet-cavity kinematical relationship
Fig. 9.
Model DR_P at t = i = ◦ . The results are presented for four di ff erent values ofthe jet–ambient gas mixing degree j mix (from top to bottom): 0.0, 0.3, 0.6, 1.0. (left) Maps of H column density obtained by integration over thevelocity range between 1 and +
150 km s − . (right) Mass–velocity relationship obtained for the molecular gas (black) and the fitted curve (red).The exponent v is given for each j mix value. fit m ( v ) ∝ exp( − v / v ). Hence, our simulations show that the ex-ponential shape of the mass–velocity distribution in the outflowis a rather general result. We note that distributions are somewhatirregular when considering only the ambient material ( j mix = v ∼ − . These values are also simi-lar to those of the global outflow mass–velocity distribution, asdisplayed in Fig. 9.A higher value of v of the order of 10 km s − is found atposition A4, close to the apex of the outflow cavity. Hence, itappears that if we leave aside the head of the outflow, the mass–velocity distributions display only modest variations across theoutflow cavity, and do not bear signatures of the ejection process(intermittency, precession).A closer look at the distributions shows that the mass–velocity distributions are better described by two componentstowards positions A1 to A3, with a change of slope (index) near v =
25 km s − . In the high-velocity range, a shallower distri-bution is observed, which is related to the jet-entrained mate-rial. In order to explore the sensitivity of the fitting parametersto the geometry and the age of the outflow, we extracted themass–velocity relationships at positions A1-A4 at three di ff erenttimes in the simulation of model DR_P, namely 400, 1200, and2000 yr, and for three values of the inclination angle with respectto the plane of the sky: i =
10, 30, and 60 ◦ . The fitting results ofthe mass–velocity relationships are summarized in Table 2.
5. Discussion
Lefloch et al. (2012) modeled the observed intensity–velocitydistribution I CO of the five outflow sources previously observedby Bachiller & Tafalla (1999): Mon R2, L1551, NGC2071,Orion A, and L1448. These latter authors showed that the ob- served I CO could be well fitted by an exponential function withvalues of v between 1.6 (Mon R2) and 12.5 (L1448) km s − .A quick inspection of the grid of models presented in Ta-ble 2 shows that these values fall well within the range of valuespredicted in our simulations, depending on the outflow age, theinclination angle, and the degree of jet–ambient gas mixing ra-tio. A value as low as 1.6 for Mon R2 is indeed easily accountedfor if the outflow propagates close to the plane of the sky, as pro-posed by DC03. For an inclination angle of 10 ◦ , and a time of1200-2000 yr, our modeling predicts low values v ∼ . − . j mix = v values of 6 . . − (according to Lefloch et al. 2012). These values of v are easily accounted for in the simulations at early ages (400-1000 yr) with an inclination angle of 60 ◦ . Solutions with a lowerinclination angle and a di ff erent degree of jet–ambient gas mix-ing ratio are possible. Detailed modeling of the sources is neces-sary to disentangle the impact of the di ff erent parameters. Inter-estingly, the spectra for L1448 presented in Fig. 4 in Lefloch etal. (2012) present a bump at v ∼
60 km s − , which we interpretas the signature of the driving jet.From comparison with the outflow sample of Bachiller& Tafalla (1999), a scenario emerges in which more-evolvedsources (like Mon R2) are better described by the swept-up gas(i.e, no entrainment; j mix = v ∼
10 km s − . In this section, we apply our numerical models to L1157 in or-der to better understand the origin of the CO intensity–velocitydistributions reported by Lefloch et al. (2012) in the southern
Article number, page 11 of 15 & A proofs: manuscript no. ms39269_last
Fig. 10.
Model DR_P at t = ◦ . (top) Map of H column density for j mix = .
0. The location of the aper-tures A1-A4 used to extract the mass–velocity relationships are drawnwith black circles. (bottom)
Mass–velocity relationship averaged overthe circular apertures A1 to A4. The best exponential fits are drawn indashed red. The value of v is given for each value of jet–ambient ma-terial mixing ratio j mix =
0, 0.3, 0.6, 1.0. outflow lobe of L1157. We note that our goal is not to provide adetailed modeling of the L1157 southern outflow lobe. For thisreason, we focus on the signatures associated with the B1 out-flow cavity.As mentioned in Section 3.2, the simulation parameters ofmodel DR_P were chosen to describe the behavior of a "typi-cal" precessing jet. For this reason, and taking into account thesimplicity of the underlying hypothesis of our model, we did notattempt to fine-tune the simulation parameters in order to obtainthe "best-fitting" model.
As mentioned in Section 2, several authors have investigated thedetails of the CO emission from the southern lobe of the L1157outflow. As first shown by Gueth, Guilloteau & Bachiller (1996),the precessing protostellar jet has shaped the southern lobe intotwo shells (outflow cavities) whose apexes are associated withthe molecular shock positions B1 and B2. Detailed modeling ofthe CO gas kinematics by Podio et al. (2016) showed that the jetprecesses on a cone inclined by 73 ◦ to the line of sight, with anopening angle of 8 ◦ on a period of 1640 yr. The modeling of the Table 2.
Model DR_P. Best exponential fitting parameters to the molec-ular mass–velocity distributions obtained towards positions A1–A4 inan aperture of 5000 au at t = i = ◦ , 30 ◦ , and 60 ◦ . j mix Age Aperture v ( km s − )(years) i = ◦ i = ◦ i = ◦ v i ± σ v i ( km s − ): 3.4 ± ± ± v i ± σ v i ( km s − ): 2.9 ± ± ± v i ± σ v i ( km s − ): 3.1 ± ± ± v i ± σ v i ( km s − ): 5.0 ± ± ± authors indicates that an angle of ≈ ◦ exists between the jetand the line of sight at the location of B1.The top-right panel of Fig. 11 shows the CO intensity–velocity distribution as observed in the J = Article number, page 12 of 15erqueira et al.: Jet-cavity kinematical relationship
Fig. 11. (Left)
Model DR_P at t = ◦ (towards the observer). (top) Map of H column density for j mix = .
9, in thevelocity intervals − [1 −
25 km s − ] (top) and − [30 −
45 km s − ] (bottom). The 5000 au aperture used to extract the mass–velocity distribution is drawnby a circle. (bottom) Molecular mass–velocity distribution extracted at position A4. We have superposed the best exponential fits m ( v ) ∝ exp( − v / v i )to the components associated with velocity intervals − [1 −
25 km s − ] (red) and − [30 −
45 km s − ] (blue). The exponent value v i is given for bothvelocity intervals. (Right) ASAI observations of L1157-B1. (top)
Intensity–velocity distribution obtained in the CO J = (cid:48)(cid:48) at the distance to the source). The line profile is fitted by a linear combination of two exponential functions g ∝ exp( − v / .
5) (blue), g ∝ exp( − v / .
4) (red) (from Lefloch et al. 2012). (bottom)
Associated mass–velocity distribution, adopting a standardCO-to-H abundance ratio and the excitation conditions derived by Lefloch et al. (2012) for the velocity components g and g . – The CO line profiles profiles are the sum of up to three com-ponents of specific excitation and velocity range, dubbed g , g , g , all of which can be modeled by an exponential lawwith a specific exponent v : 12.5, 4.4, 2.5, respectively. – The component of lowest excitation ( T ex =
23 K) and narrow-est velocity range, ( −
5; 0 km s − ), dubbed g , is detected overthe whole southern lobe and is the only component detectedtowards the southernmost, older cavity associated with theB2 shock. – The component of highest excitation ( T ex =
210 K) and high-est velocity range, [ − −
20] km s − , dubbed g , is detectedclose to the apex of the younger cavity associated with theB1 shock. – The component of intermediate excitation ( T ex =
64 K) andvelocity range, − < v < − − , dubbed g is detectedover the whole outflow cavity associated with B1. The excitation conditions of the CO gas in L1157 make it espe-cially easy to obtain the mass–velocity distribution from the COintensity–velocity distribution. This is because the CO J = g , g , g are indepen-dent of the velocity and the line emission is optically thin, butat velocities very close (a few km s − ) to ambient. Hence, the simple relation that exists at local thermodynamic equilibriumbetween N(CO) and the CO line flux (see e.g., Bachiller et al.1990) can be applied to the whole velocity range of emissionof each component. In practice, each component dominates overa specific velocity interval of the intensity–velocity distribution(see top right panel of Fig. 11). The mass–velocity distributionis therefore immediately obtained when considering the excita-tion conditions and the size of the main emitting gas componentas a function of velocity. The total mass–velocity relationshipis rigorously obtained by multiplying the relationship N(CO)(v)—which is derived from T b (CO)(v)— by the CO emission areaat each velocity interval. Despite the uncertainties in the over-lap region between components (e.g., near v = −
25 km s − ),the spectral slope of each component is found to be preservedin the derivation procedure from the CO intensity to the mass–velocity distribution. This is illustrated in the bottom-right panelof Fig. 11 in which we report the mass–velocity distribution to-wards L1157-B1. We note that we assume a standard abundanceratio CO / [H ] = − . Figure 11 presents the distributions of the molecular ma-terial in the velocity intervals [ − −
1] km s − (top) and Article number, page 13 of 15 & A proofs: manuscript no. ms39269_last [ − −
30] km s − , respectively. We also computed the molec-ular mass–velocity distribution measured in an aperture of 5000au close to the apex of the cavity (position A4). This is compara-ble to the beamwidth (HPBW) of the IRAM 30m telescope mainbeam at the frequency of the CO J = g and g over the B1 and B2 lobes,respectively. Our model DR_P is consistent with the interpreta-tion that g and g are associated with di ff erent ejection eventsresponsible for the formation of B1 and B2 cavities, respectively.The lower excitation conditions of g are consistent with an olderevent. The di ff erence of spectral slope ( v ) between the B1 andB2 cavities could be explained a priori by a higher jet inclina-tion to the line of sight in the direction of B2. However, thiscontradicts the kinematic modeling of the jet precession by Po-dio et al. (2016), which predicts an inclination angle to the lineof sight of about 25 ◦ at shock position B2, lower than towardsB1, and therefore favors a higher jet radial velocity than thatmeasured towards B1. The sensitivity of the millimeter CO linespectra available in the literature (see e.g., Lefloch et al. 2012)is not high enough to allow conclusions to be made about thepresence of the jet towards B2. On the other hand, our numericalsimulations show that the low-velocity ( v < − ) emissionactually arises from entrained ambient material in the outflowcavity walls, and has very little dependence on the driving jet.The jet that once created the B2 shock about 2500 yr ago (Po-dio et al. 2016) is now impacting the B1 cavity at the B1 and B0positions as a result of its precession. Hence, the emission fromB2 arises mainly from previously entrained, ambient material,which is now being slowed down. Inspection of Table 2 showsthat low values of v are also obtained in the swept-up ambientmolecular gas ( j mix = . As can be seen in Fig. 11, the mass–velocity distribution ex-tracted towards position A4, close to the apex of the cavityin the simulation, shows the presence of two distinct compo-nents associated with the velocity intervals [ − −
25 km s − ] and[ − −
45 km s − ], respectively. This situation is reminiscent ofthe CO (and the mass) intensity–velocity distribution observedtowards B1. Both components can be fitted by an exponentialfunction of exponent v (cid:39) . v (cid:39)
13, respectively. Thesevalues are in good agreement with those determined for compo-nents g and g towards the B1 position. We note that accordingto our modeling (see Table 2), these values are mainly sensitiveto the jet inclination angle to the line of sight. We note that theyweakly depend on the actual value of the jet–ambient materialmixing ratio, but the best agreement was obtained for j mix = −
30 km s − , as displayed in Fig. 11, isnot restricted to a few spots of shocked gas, such as for examplethe jet impact shock region at the apex of the outflow cavity. In-stead, it turns out that the high-velocity material ( | v | >
30 km s − )is tracing an elongated, collimated structure surrounding the jetthroughout the whole outflow cavity. This elongated structure issurrounded by the lower velocity material (0 < | v | <
25 km s − )of the outflow. In other words, the high-velocity material does We adopted the revised distance of 372 pc to L1157 (Zucker et al.2019). not trace only the jet shock impact region (the Mach disk).In our simulation, the amount of molecular material at the jethead strongly decreases as a result of molecular dissociation.Instead, the high-velocity component arises from material en-trained along the jet.
Observational evidence of the high-velocity component has beenreported in the millimeter rotational transitions of a few molecu-lar species, such as for example CO (Lefloch et al. 2012), HCO + (Podio et al. 2014), and SiO (Tafalla et al. 2010; Spezzano et al.2020). Unfortunately, the interferometric observations of L1157available in the literature focus mainly on the gas propagating atvelocities close to ambient, which is associated with the bow andthe outflow cavity walls. Therefore, the evidence for the high-velocity jet is still very scarce and unambiguous detection ofthe molecular material entrained along the jet is still missing.It is worth noting that Plateau de Bure observations of the SiO J = v < −
10 km s − . Single-dish observations of the SiO J =
6. Conclusions
Using the hydrodynamical code Yguazú-a, we performed 3Dnumerical simulations to revisit in a detailed manner the mass–velocity relationship in jet-driven molecular outflows. Great at-tention was paid to benchmark Yguazu-a against the hydrody-namical codes used by previous authors in the field (Downes &Ray 1999; Downes & Cabrit 2003). To do so, we modeled thepropagation of an intermittent jet adopting the same parametersas those of Downes & Ray (1999) and Downes & Cabrit (2003).We find excellent quantitative agreement between our simula-tions and those of these latter authors.Detailed comparison between our simulations and those ofDownes & Cabrit (2003) leads us to conclude that these latterauthors took into account the jet material contribution in the ob-tention of the mass–velocity distribution presented in their work.We find that the presence of a bump in the high-velocity range( v ∼
100 km s − ) is remarkable evidence of the presence of thejet and that all the previous works considered, to a greater orlesser extent, the presence of the jet in their mass–velocity pro-file computations.Overall, our simulations show that the mass–velocity distri-bution of the outflowing material can be successfully fitted byone exponential law m ( v ) ∝ exp( − v / v ). We systematically in-vestigated the signature of the mass–velocity distribution as afunction of time, depending on the jet inclination to the line ofsight and the degree of mixing between the jet and the ambi-ent material. We find that it may be necessary to introduce a Article number, page 14 of 15erqueira et al.: Jet-cavity kinematical relationship second component to fit the mass–velocity distribution in thehigh-velocity range. The spectral signature in the low-velocityrange is dominated by the contribution of material in the outflowcavity walls and is rather insensitive to the actual value of thejet–ambient gas mixing ratio. The synthetic mass–velocity dis-tributions from our simulations are good agreement with distri-butions derived from observations and we are able to reproducethe observational data when taking age and source geometry intoaccount.We verified that the profile of the mass–velocity distributioncomputed over a local area inside the outflow can still be wellfitted by an exponential function. The profiles and the v valuesare very similar over the outflow, but the distribution appearsmuch shallower at the apex of the outflow cavity as a result ofthe leading jet contribution.We performed a simple modeling of the L1157 southern out-flow cavity by simulating of a precessing jet with parameterssimilar to those reported by Podio et al. (2016). We were ableto reproduce the main features of the CO intensity–velocity dis-tributions observed in the southern outflow lobe of L1157 to asatisfactory degree. Our simulations suggest that the three com-ponents identified by Lefloch et al. (2012) are all related to theentrained gas and not the jet impact shock region, that is, theMach disk itself. High-angular-resolution observations with theNOEMA interferometer could easily test this conclusion. Acknowledgements.
We would like to thank to an anonymous referee, whosesuggestions has contributed to improve the presentation of the paper. A.H.Cerqueira would like to thanks the PROPP–UESC for partial funding (un-der project No. 073.6766.2019.0010667-91) as well as the Brazilian AgencyCAPES. B. Lefloch, P.R. Rivera-Ortiz and C. Ceccarelli acknowledge fundingfrom the European Research Council (ERC) under the European Union’s Hori-zon 2020 research and innovation programme, for the Project “The Dawn ofOrganic Chemistry” (DOC), grant agreement No 741002.
References
André, P., Ward-Thompson, D., & Barsony, M. 2000, in Protostars and PlanetsIV (Univ. of Arizona Press; eds Mannings, V., Boss, A.P., Russell, S. S.)Arce, H.G., & Goodman, A.A. 2001, ApJ, 551, L171Arce, H.G., Shepherd, D., Gueth, F. et al. 2007, in Protostars and Planets V,(Univ. of Arizona Press; eds. Reipurth, D. Jewitt, and K. Keil)Bachiller, R., Cernicharo, J., Martin-Pintado, J., et al. 1990, A&A, 231, 174Bachiller, R., & Tafalla, M. 1999, in
The Origin of Stars and Planetary Systems ,ed. C. J. Lada & N. D. Kylafis (Dordrecht: Kluwer Academic Publishers),227Bally, J. 2016, ARA&A, 54, 491Bally, J., Devine, D., & Alten, V. 1997, ApJ, 478, 603Busquet, G., Lefloch, B., Benedettini, M., et al. 2014, A&A, 561, 120Caratti o Garatti, A., Giannini, T., Nisini, B., et al. 2006, A&A, 449, 1077Castellanos-Ramírez, A., Raga, A.C., & Rodríguez-González, A. 2018, ApJ,867, 29Castellanos-Ramírez, A., Rodríguez-González, A., Rivera-Ortiz, P., et al. 2018,Rev. Mexicana Astron. Astrofis., 54, 409Cerqueira, A.H., Velázquez, P.F., Raga, A.C., et al. 2006, A&A, 448, 231Chernin, C.R., & Mason, C.R. 1993, ApJ, 414, 230Downes, T.P., & Cabrit, S. 2003, A&A, 403, 135 (DC03)Downes, T.P., & Ray, T.P. 1999, A&A, 345, 977Flower, D.R., & Pineau des Forêts, G., 2010, MNRAS, 406, 1745Frank, A., Ray, T.P., Cabrit, S., et al. 2014, in
Protostars and Planets VI , H.Beuther, R. Klessen, C.P. Dullemond and T. Henning (eds.), Univ. of ArizonaPress, Tucson, p.451-475Gómez-Ruiz, A.I., Codella, C., Lefloch, B., et al. 2015, MNRAS, 446, 3346Gueth, F., Guilloteau, S., & Bachiller, R., 1996, A&A, 307, 891Gueth, F., Guilloteau, S., & Bachiller, R., 1998, A&A, 333, 287Hollenbach, D., & McKee, C.F. 1979, ApJS, 41, 555Kaufman, M.J., & Neufeld, D.A. 1996, ApJ, 456, 611Keegan, R., & Downes, T.P. 2005, A&A, 437, 517Lee, C-.F., Ho, T.P., Li, Z-.Y., et al. 2017, Nature Astronomy, 1, 152Lee, C-.F. 2020, A&A Rev., 28, 1Lefloch, B., Cabrit, S., Busquet, G., et al. 2012, ApJ, 757, L25Lefloch, B., Gusdorf, A., Codella, C., et al. 2015, A&A, 581, L4 Lepp, S., & Shull, M. 1983, ApJ, 270, 578Mendoza, E., Lefloch, B., López-Sepulcre, A., et al. 2014, MNRAS, 445, 151Mendoza, E., Lefloch, B., Ceccarelli, C., et al. 2018, MNRAS, 475, 5501Nisini, B., Codella, C., Giannini, T., et al. 2007, A&A, 462, 163Nisini, B., Benedettini, M., Codella, C., et al. 2010a, 518, L120Nisini, B., Giannini, T., Neufeld, D.A., 2010b, ApJ, 724, 69Podio, L., Bacciotti, F., Nisini, B., et al. 2006, A&A, 456, 189Podio, L., Lefloch, B., Ceccarelli,C., et al. 2014, A&A, 565, 64Podio, L., Codella, C., Gueth, F., et al. 2016, A&A, 593, L4Raga, A.C, & Cabrit, S. 1993, A&A, 278, 267Raga, A. C., de Gouveia Dal Pino, E. M., Noriega-Crespo, A., et al. 2002, A&A,392, 267Raga, A.C., Navarro-González, R. & Villagrán-Muniz, M. 2000, Rev. MexicanaAstron. Astrofis., 36, 67Reipurth, B., & Bally, J. 2001, ARA&A, 39, 403Reipurth, B., Yu, K.C., Heathcote, S., et al. 2000, AJ, 120, 1449Reipurth, B., Heathcote, S., Morse, J., et al. 2002, AJ, 123, 362Ridge, N.A. & Moore, T.J.T. 2001, A&A, 378, 495Rosen, A., & Smith, M.D. 2003, MNRAS, 343, 181Rosen, A., & Smith, M.D. 2004a, A&A, 413, 593Rosen, A., & Smith, M.D. 2004b, MNRAS, 347, 1097Santangelo, G., Codella, C., Cabrit, S., et al. 2015, A&A, 584, 126Shapiro, P.R., & Kang, H. 1987, ApJ, 318, 32Smith, M.D., & Rosen, A. 2005, MNRAS, 357, 579Smith, M.D., & Rosen, A. 2007, MNRAS, 378, 691Spezzano, S., Codella, C., Podio, L., et al. 2020, arXiv:2006.05772Stahler, S.W., & Palla, F. 2004, in The Formation of Stars, Wiley-VCHTafalla, M., Santiago-García, J., Hacar, A., et al. 2010, A&A, 522, A91Taylor, S.D., & Raga, A.C. 1995, A&A, 296, 823Wu, Y., Wei, Y., Zhao, M. et al., 2004, A&A, 426, 503Zhang, Q., & Zheng, X. 1997, ApJ, 474, 719Zinnecker H., McCaughrean M.J., Rayner J.T. 1998, Nature 394, 862Zucker, C., Speagle, J.S., Schlafly, E.F., et al. 2019, ApJ, 879, 125, H.Beuther, R. Klessen, C.P. Dullemond and T. Henning (eds.), Univ. of ArizonaPress, Tucson, p.451-475Gómez-Ruiz, A.I., Codella, C., Lefloch, B., et al. 2015, MNRAS, 446, 3346Gueth, F., Guilloteau, S., & Bachiller, R., 1996, A&A, 307, 891Gueth, F., Guilloteau, S., & Bachiller, R., 1998, A&A, 333, 287Hollenbach, D., & McKee, C.F. 1979, ApJS, 41, 555Kaufman, M.J., & Neufeld, D.A. 1996, ApJ, 456, 611Keegan, R., & Downes, T.P. 2005, A&A, 437, 517Lee, C-.F., Ho, T.P., Li, Z-.Y., et al. 2017, Nature Astronomy, 1, 152Lee, C-.F. 2020, A&A Rev., 28, 1Lefloch, B., Cabrit, S., Busquet, G., et al. 2012, ApJ, 757, L25Lefloch, B., Gusdorf, A., Codella, C., et al. 2015, A&A, 581, L4 Lepp, S., & Shull, M. 1983, ApJ, 270, 578Mendoza, E., Lefloch, B., López-Sepulcre, A., et al. 2014, MNRAS, 445, 151Mendoza, E., Lefloch, B., Ceccarelli, C., et al. 2018, MNRAS, 475, 5501Nisini, B., Codella, C., Giannini, T., et al. 2007, A&A, 462, 163Nisini, B., Benedettini, M., Codella, C., et al. 2010a, 518, L120Nisini, B., Giannini, T., Neufeld, D.A., 2010b, ApJ, 724, 69Podio, L., Bacciotti, F., Nisini, B., et al. 2006, A&A, 456, 189Podio, L., Lefloch, B., Ceccarelli,C., et al. 2014, A&A, 565, 64Podio, L., Codella, C., Gueth, F., et al. 2016, A&A, 593, L4Raga, A.C, & Cabrit, S. 1993, A&A, 278, 267Raga, A. C., de Gouveia Dal Pino, E. M., Noriega-Crespo, A., et al. 2002, A&A,392, 267Raga, A.C., Navarro-González, R. & Villagrán-Muniz, M. 2000, Rev. MexicanaAstron. Astrofis., 36, 67Reipurth, B., & Bally, J. 2001, ARA&A, 39, 403Reipurth, B., Yu, K.C., Heathcote, S., et al. 2000, AJ, 120, 1449Reipurth, B., Heathcote, S., Morse, J., et al. 2002, AJ, 123, 362Ridge, N.A. & Moore, T.J.T. 2001, A&A, 378, 495Rosen, A., & Smith, M.D. 2003, MNRAS, 343, 181Rosen, A., & Smith, M.D. 2004a, A&A, 413, 593Rosen, A., & Smith, M.D. 2004b, MNRAS, 347, 1097Santangelo, G., Codella, C., Cabrit, S., et al. 2015, A&A, 584, 126Shapiro, P.R., & Kang, H. 1987, ApJ, 318, 32Smith, M.D., & Rosen, A. 2005, MNRAS, 357, 579Smith, M.D., & Rosen, A. 2007, MNRAS, 378, 691Spezzano, S., Codella, C., Podio, L., et al. 2020, arXiv:2006.05772Stahler, S.W., & Palla, F. 2004, in The Formation of Stars, Wiley-VCHTafalla, M., Santiago-García, J., Hacar, A., et al. 2010, A&A, 522, A91Taylor, S.D., & Raga, A.C. 1995, A&A, 296, 823Wu, Y., Wei, Y., Zhao, M. et al., 2004, A&A, 426, 503Zhang, Q., & Zheng, X. 1997, ApJ, 474, 719Zinnecker H., McCaughrean M.J., Rayner J.T. 1998, Nature 394, 862Zucker, C., Speagle, J.S., Schlafly, E.F., et al. 2019, ApJ, 879, 125