Turbulent mixing layers in supersonic protostellar outflows, with application to DG Tauri
Marc C. White, Geoffrey V. Bicknell, Ralph S. Sutherland, Raquel Salmeron, Peter J. McGregor
MMNRAS , 000–000 (0000) Preprint 6 September 2018 Compiled using MNRAS L A TEX style file v3.0
Turbulent mixing layers in supersonic protostellar outflows, withapplication to DG Tauri
M. C. White (cid:63) , G. V. Bicknell, R. S. Sutherland, R. Salmeron, P. J. McGregor † Research School of Astronomy & Astrophysics, The Australian National University, Canberra, ACT 2611, Australia
ABSTRACT
Turbulent entrainment processes may play an important role in the outflows from youngstellar objects at all stages of their evolution. In particular, lateral entrainment of ambientmaterial by high-velocity, well-collimated protostellar jets may be the cause of the multipleemission-line velocity components observed in the microjet-scale outflows driven by classicalT Tauri stars. Intermediate-velocity outflow components may be emitted by a turbulent, shock-excited mixing layer along the boundaries of the jet. We present a formalism for describingsuch a mixing layer based on Reynolds decomposition of quantities measuring fundamentalproperties of the gas. In this model, the molecular wind from large disc radii provides a con-tinual supply of material for entrainment. We calculate the total stress profile in the mixinglayer, which allows us to estimate the dissipation of turbulent energy, and hence the lumi-nosity of the layer. We utilize
MAPPINGS
IV shock models to determine the fraction of totalemission that occurs in [Fe II ] 1.644 µ m line emission in order to facilitate comparison toprevious observations of the young stellar object DG Tauri. Our model accurately estimatesthe luminosity and changes in mass outflow rate of the intermediate-velocity component ofthe DG Tau approaching outflow. Therefore, we propose that this component represents a tur-bulent mixing layer surrounding the well-collimated jet in this object. Finally, we compareand contrast our model to previous work in the field. Key words:
MHD – stars: individual: DG Tauri – stars: jets – stars: protostars – methods:analytical
Outflows are a near-universal component of young stellar objects(YSOs) throughout their evolution. They play a major role in starformation, and drive both the CO outflows seen in early-stage form-ing stars (e.g. Bachiller 1996; Reipurth & Bachiller 1997) and theHerbig-Haro flows emanating from more mature protostars (e.g.Reipurth & Bally 2001). These outflows are thought to be launchedeither from the protostellar surface (e.g. Ferreira et al. 2006), mag-netocentrifugally from magnetic reconnection points near the cir-cumstellar disc truncation radius (the X-wind and related models;Shu et al. 1994; Romanova et al. 2009) or from the disc surface at (cid:63) [email protected] † This work would not have been possible without the key involvement ofProfessor Peter McGregor. His leadership in the construction of the Near-infrared Integral Field Spectrograph on the Gemini North telescope, hiscommissioning observations of DG Tau, and his insightful contributionsto the subsequent research, were all vital to the success of this project. Petersadly passed away in March 2015. We will always remember him and aregrateful for his many contributions to astronomy, including his mentoringof students, as well as his kindness and good humour. He will be sorelymissed. larger radii (disc winds; Blandford & Payne 1982; Pudritz & Nor-man 1983). In fact, more than one launch mechanism may be inoperation (e.g. Larson 2003).The advent of the
Hubble Space Telescope and adaptive op-tics on ground-based telescopes has allowed the few hundred au ofprotostellar outflows closest to the protostar to be studied. The out-flows associated with optically-revealed T Tauri stars take the formof well-collimated ‘microjets’. The study of these microjets is im-portant because it is thought that they should not have interactedwith the wider interstellar medium so close to the star (althoughinteractions may occur if there is a remnant protostellar envelope,e.g. White et al. 2014b). If so, such observations may provide in-formation on the outflow before it significantly interacts with theambient medium.The small-scale outflows from YSOs typically show an onion-like kinematic structure in optical and near-infrared (NIR) forbid-den lines, with a well-collimated, high-velocity jet surrounded by aless-collimated, intermediate-velocity component (e.g. Hirth et al.1997; Woitas et al. 2002; Pyo et al. 2003; Coffey et al. 2008;Rodr´ıguez-Gonz´alez et al. 2012; Caratti o Garatti et al. 2013). Thenature of the high-velocity jets in many sources has been studiedextensively, including searches for signs of jet rotation c (cid:13) a r X i v : . [ a s t r o - ph . S R ] O c t M. C. White et al. (Bacciotti et al. 2000; Coffey et al. 2004, 2007; White et al. 2014a,hereafter Paper I), recollimation shocks (Paper I; G¨udel et al. 2005,2008; G¨unther et al. 2009; Bonito et al. 2011; Schneider et al. 2013)and studies of the propagation of shock-excited moving knots (Pa-per I; Burrows et al. 1996; Reipurth et al. 2002; Pyo et al. 2003;Agra-Amboage et al. 2011). The nature of the intermediate-velocityemitting material is still debated. Many authors attribute this emis-sion to the presence of an intermediate-velocity disc-wind outflowcomponent (e.g. Podio et al. 2011), which bridges the gap in launchradii between low-velocity ( (cid:46) km s − ) molecular winds (e.g.Takami et al. 2004; Beck et al. 2008; Agra-Amboage et al. 2014)and high-velocity ( > km s − ) jets. However, such an expla-nation does not provide a natural mechanism for the generationof forbidden optical- and NIR-line emission, which is attributed toshock excitation (e.g. Nisini et al. 2002). It has yet to be explainedhow a steady-state intermediate-velocity disc wind would undergoshock excitation with relatively uniform intensity along the observ-able length of the feature.It has been proposed that the intermediate-velocity forbidden-line emission components (IVCs) of small-scale protostellar out-flows result from the lateral entrainment of ambient material, or adisc wind, by the high-velocity jet (e.g. Pyo et al. 2003). This sug-gestion is based on the observation that in some objects, the spatialwidth of the intermediate-velocity component increases with dis-tance from the central star. Entrainment would cause the formationof a turbulent mixing layer between the supersonic jet and the ma-terial surrounding it, which would become shock-excited and emitin forbidden lines (e.g. Binette et al. 1999). Such a layer naturallygrows in thickness with distance along the jet, reproducing the ob-servations of low- to intermediate-velocity forbidden-line emissioncomponents in protostellar outflows (e.g. Cant´o & Raga 1991; Ragaet al. 1995).The entrainment explanation has fallen out of favour recentlyfor two reasons. First, jet simulations show that the jet pushes theambient medium aside as it is launched (L´opez-C´amara & Raga2010), preventing ambient material from interacting with the sidesof the jet. However, as suggested by Pyo et al. (2003) and Whiteet al. (2014a), the presence of a wide-angle molecular wind sur-rounding the central jet could supply a constant reservoir of mate-rial for entrainment by the jet. Secondly, hypersonic jets, such asprotostellar microjets, should not form lateral entrainment layersif they are regarded as high Mach number, purely hydrodynamicflows. The formation of turbulent mixing layers is driven by theaction of the Kelvin-Helmholtz (KH) instability at the jet-ambientmaterial interface, and the growth rate of the KH instability de-creases as the Mach number difference between the flows increases(Chernin et al. 1994; Trussoni 2008). However, protostellar jets areexpected to exhibit strong toroidal magnetic fields (e.g. Zanni et al.2007). The alignment of these fields with respect to the interfacebetween the jet and the surrounding material, and their perpendicu-larity to the flow, may destabilise the interface to the KH instability(Paper I; Miura & Pritchett 1982; Ray & Ershkovich 1983). There-fore, entrainment remains an open possibility in protostellar jets.For a more detailed discussion of this argument, see Paper I, § § §
3, we first com-pute a grid of shock models to determine the ratio between the ob-servable [Fe II ] line emission of protostellar jet mixing layers, andthe mixing layer bolometric luminosity estimated by our model.We then directly compare our model to the [Fe II ] IVC of the ap-proaching outflow from the YSO DG Tauri, and find that it is inexcellent agreement with observations. § § We construct an analytical, semi-empirical model of a two-dimensional turbulent entrainment layer in order to interpret the [Fe II ] 1.644 µ m IVC line emission observed in DG Tau. The modeldescribes the turbulent mixing layer that forms between a high-velocity jet and a low-velocity wider-angle wind, and depends onlyupon directly observable quantities, removing the requirement tospecify an ‘entrainment efficiency’ parameter (e.g. Cant´o & Raga1991; Raga et al. 1995). We use the observed spreading rate of thelayer to calculate the dissipation of turbulent energy in the entrain-ment layer, and its resulting luminosity. The first step in this processis the calculation of total turbulent stress, t xy , in the mixing layer.The model setup is shown in Fig. 1. A high-velocity jet withdensity ρ jet propagates at velocity v jet away from the star-disc sys-tem. The jet is surrounded by a wider-angle molecular wind, withdensity ρ w and velocity v w (cid:28) v jet . We henceforth refer to this windas the ‘ambient wind’. A turbulent mixing layer forms at the inter-face between the two flows, as a result of the KH instability. Weapproximate this interface with a two-dimensional model. In thismodel, the x -axis is defined as the unperturbed jet-wind boundary.This is the streamwise direction; the transverse coordinate is y . Themixing layer width, h ( x ) , increases monotonically with distancefrom the central star. We define the depth the mixing layer expandsinto the jet as h ( x ) , and the depth it penetrates the ambient windas h ( x ) , where h ( x ) < . Hence, h ( x ) = h ( x ) − h ( x ) .An averaging prescription is used to describe the mean flow. A two-dimensional model is an adequate representation of the jet edge,given that the mixing layer width is not significantly greater than the jetradius (see below). MNRAS , 000–000 (0000) urbulent mixing layers in protostellar outflows w i d e - a n g l e , l o w - v e l o c i t y w i nd circumstellar disk h i g h - v e l o c i t y , w e ll - c o lli m a t e d j e t mixing layer xy h ( x ) h ( x ) h ( x )jetwindmixing layer ce n t r a l s t a r v jet v w v jet Figure 1.
A representation of the model setup used throughout this paper. Ahigh-velocity, well-collimated jet is launched from a protostar-circumstellardisc system. This jet is surrounded by a wider-angle disc wind (top panel).The interface between the jet and the wind is approximated by a two-dimensional turbulent shear layer (bottom panel). The x -axis of the modelis placed where the jet-wind interface would lie in the absence of the mixinglayer, and is parallel to the direction of the jet. Dashed arrows show the flowdirection of the components. Model components are not to scale. We adopt the mass-weighted statistical averaging prescription ofFavre (1969), which was introduced to the study of astrophysicalflows by Bicknell (1984); see also Kuncic & Bicknell (2004). Alltime-varying quantities are decomposed into an average componentand a fluctuating component. Quantities such as pressure, p , den-sity, ρ , and magnetic field, B , are expressed in terms of mean (bar)and fluctuating (primed) components, such that the time-average ofthe fluctuating component (angle brackets) is zero: p = ¯ p + p (cid:48) where (cid:104) p (cid:48) (cid:105) = 0 ; (1) ρ = ¯ ρ + ρ (cid:48) where (cid:104) ρ (cid:48) (cid:105) = 0 ; and (2) B i = ¯ B i + B (cid:48) i where (cid:104) B (cid:48) i (cid:105) = 0 , (3)where subscripts i and j represent generalized coordinates. As pre-scribed by Favre (1969), the velocity, v i , is mass-weighted, and isexpressed as v i = ˜ v i + v (cid:48) i , where (cid:104) ρv (cid:48) i (cid:105) = 0 . (4)This approach has two advantages. First, mass is conserved in themean flow (Favre 1969). Secondly, it prevents the generation ofan excessive number of terms when the dynamical equations arestatistically averaged; e.g. the mean value of the momentum flux issimply expressed as (cid:104) ρv i v j (cid:105) = ¯ ρ ˜ v i ˜ v j + (cid:104) ρv (cid:48) i v (cid:48) j (cid:105) . (5) This approach is common in fluid dynamics, and has been used inthe theory of compressible turbulent jets and accretion discs (seeBicknell 1984; Kuncic & Bicknell 2004). Consider a compressible magnetized fluid with density ρ , velocity v , pressure p , magnetic field B , in a gravitational potential field φ G .Averaging the mass continuity and momentum conservation equa-tions of magnetohydrodynamics ( § A1) yields, for a quasi-steadystate system, ∂ (¯ ρ ˜ v x ) ∂x + ∂ (¯ ρ ˜ v y ) ∂y = 0 , and (6) ∂ (¯ ρ ˜ v i ˜ v j ) ∂x j = − ¯ ρ ∂φ G ∂x i − ∂ ¯ p∂x i + ∂ ( t Rij + t Bij ) ∂x j . (7)The magnetic stress tensor is defined as t Bij = (cid:104) B (cid:48) i B (cid:48) j (cid:105) π − δ ij (cid:104) B (cid:48) (cid:105) π , (8)assuming that the magnetic field is dominated by its turbulent com-ponent, so that ¯ B i = 0 . We define the Reynolds stress tensor as t R ij = −(cid:104) ρv (cid:48) i v (cid:48) j (cid:105) (9)(Kuncic & Bicknell 2004). The total turbulent stress is t ij = t R ij + t Bij . (10)The aim of our calculation is to estimate the mass entrain-ment rate and bolometric luminosity of the mixing layer, basedon directly observed parameters. We consider the i = x momen-tum equation, that is, the equation governing the streamwise evolu-tion of momentum resulting from the lateral transfer of momentumwithin the mixing layer. We neglect the streamwise pressure, mag-netic and gravitational gradients, so that ∂ (¯ ρ ˜ v x ) ∂x + ∂ (¯ ρ ˜ v x ˜ v y ) ∂y ≡ ¯ ρ ˜ v x ∂ ˜ v x ∂x + ¯ ρ ˜ v y ∂ ˜ v x ∂y = ∂t xy ∂y , (11)where t xy = t Rxy + t Bxy as per equation (10).We analyse the orders of magnitude of the terms in equations(6) and (11) in Appendix A2. We denote the characteristic advec-tive length scale of the mixing layer in the x -direction as L . In thisanalysis, we show that the following order of magnitude relation-ships exist in our model: ˜ v y ∼ hL ˜ v x , (12) v (cid:48) ∼ (cid:18) hL (cid:19) / ˜ v x , and (13) h (cid:48) ∼ hL . (14)This analysis also justifies the neglecting of the Reynolds stressterm governing the streamwise transfer of momentum, ∂t xx /∂x ,in equation (11).In order to achieve our goal of estimating the rate of dissipa-tion of turbulence, we define a pseudo self-similar variable in thetransverse direction, ξ ( x, y ) = yh ( x ) (15) ⇒ ξ = y ( x ) h ( x ) , and ξ ( x ) = y ( x ) h ( x ) . (16) MNRAS , 000–000 (0000)
M. C. White et al.
We note that ξ − ξ = 1 . It may also be shown that total (thermalplus turbulent) pressure balance is maintained across the mixinglayer ( § A2): p tot = ¯ p (cid:124)(cid:123)(cid:122)(cid:125) p thermal + (cid:104) ρv (cid:48) y (cid:105) (cid:124) (cid:123)(cid:122) (cid:125) p turbulent = const. (17)This allows us to specify the quantity W ( § A2, equation A8), andstreamwise velocity, ˜ v x , in the mixing layer as a function of trans-verse position. As the simplest approximation, we prescribe linearprofiles across the layer: ˜ v x = U ( ξ ) v jet , and (18) W = S ( ξ )( W jet − W w ) + W w , where (19) U ( ξ ) = S ( ξ ) = ξ − ξ = ξ − ξ + 1 . (20) We assume pressure equilibrium across the mixing layer. Let thejet-to-ambient wind density ratio be η = ρ jet /ρ w . The density isgiven by ¯ ρ = µmk W ¯ p tot , (21)where k is the Boltzmann constant, m is the atomic mass unit, and µ is the molecular weight of the gas. It follows that the mixing layerdensity profile is given by ¯ ρ ( η, ξ ) = ρ jet η + (1 − η ) S ( ξ ) . (22)It is important to note that the density and temperature profiles areintrinsically linked due to our assumption of constant pressure; it isnot possible for one profile to vary across the mixing layer unlessthe other profile also varies.We now calculate the transverse velocity and turbulent stressprofiles within the mixing layer. We transform equations (6) and(11) into the ( x, ξ ) coordinate system. The equation of continuity,equation (6), becomes ∂ (¯ ρ ˜ v y ) ∂ξ = h (cid:48) ( x ) ξ ∂ ¯ ρ ˜ v x ∂ξ . (23)We can calculate the tranverse velocity profile from this equation,after substituting in the density (equation 22) and streamwise ve-locity (equation 18) profiles. This gives: ˜ v y ( η, x, ξ ) = v jet h (cid:48) ( x ) [ η + (1 − η ) S ( ξ )] × (cid:90) ξξ ξ (cid:48) dd ξ (cid:48) (cid:20) U ( ξ (cid:48) ) η + (1 − η ) S ( ξ (cid:48) ) (cid:21) d ξ (cid:48) . (24)The integral factor in equation (24) occurs in many of the subse-quent expressions, and we define it as D ( η, ξ ) . This factor has aclosed-form solution ( § A3) after imposing the boundary condition ˜ v y ( ξ ) = 0 , so the transverse velocity varies smoothly from v y = 0 in the jet into the mixing layer. This is a reasonable boundary con-dition, because the supersonic jet will approach the mixing layerboundary so quickly that it will not be substantially deflected byturbulence prior to impacting the mixing layer.Following transformation into the ( x, ξ ) coordinate system, In the case of DG Tau, the presence of a stationary recollimation shockin the jet channel (Paper I) indicates that the jet is in pressure equilibriumwith its environs downstream of this shock. − − − Jet-to-wind density ratio, η − . − . . . D i m e n s i on l e ss c oo r d i n a t e ξ jetambient wind ξ ξ Figure 2.
Position of the mixing layer boundaries with the jet ( ξ ) and withthe ambient wind ( ξ ) as a function of the jet-to-wind density ratio, η . Theposition of the boundary between the jet and wind in the absence of a mixinglayer is shown by the dashed line. the equation of downstream momentum conservation, equation(11), can be rearranged to provide an equation for the turbulentstress: ∂t xy ∂ξ = ρ jet v j h (cid:48) ( x ) (cid:18) − ( ξ − ξ + 1) ξη + (1 − η )( ξ − ξ + 1) + D ( η, ξ ) (cid:19) .(25)Integration of equation (25) gives: t xy ( η, x, ξ ) = ρ jet v jet h (cid:48) ( x ) F ( η, ξ ) , (26)where the function F ( η, ξ ) is given in § A3. We set t xy ( ξ ) = t xy ( ξ ) = 0 since we expect the turbulence to be confined primar-ily to the region ξ < ξ < ξ . The condition t xy ( ξ ) = 0 is usedto compute the form of F ( η, ξ ) ( § A3); the condition t xy ( ξ ) = 0 ,and hence F ( ξ , η ) = 0 , allows us to calculate the position of thejet-mixing layer boundary, ξ , in ( x, ξ ) -space as a function of onlythe jet-to-wind density ratio η : ξ ( η ) = 2 η log( η ) + (4 − η ) η − η − for η (cid:54) = 1 , and = 13 for η = 1 . (27)The position of the mixing layer boundaries as a function of η isshown in Fig. 2. In the limit of a significantly underdense jet ( η → ), the mixing layer penetrates the jet and wind evenly. In the limitof a significantly overdense jet ( η → ∞ ), the mixing layer almostexclusively penetrates the ambient wind.Knowledge of the position of the mixing layer boundaries al-lows the forms of the transverse velocity and turbulent stress pro-files to be directly calculated as a function of position within themixing layer, ξ , and the jet-to-wind density ratio, η . The forms ofthese expressions are algebraically complex, and are given in § A4.The mixing layer density, transverse velocity and turbulent stressprofiles are shown in Fig. 3.As seen in Fig. 3(b), ambient wind material is pulled upwardsinto the mixing layer at the wind-mixing layer boundary ( ξ ) withan effective entrainment velocity, v ent = ˜ v y ( η, x, ξ ) (28) = v jet h (cid:48) ( x ) η ( η − η log( η ) − η − for η (cid:54) = 1 , and = v jet h (cid:48) ( x ) 16 for η = 1 . (29) MNRAS000
Position of the mixing layer boundaries with the jet ( ξ ) and withthe ambient wind ( ξ ) as a function of the jet-to-wind density ratio, η . Theposition of the boundary between the jet and wind in the absence of a mixinglayer is shown by the dashed line. the equation of downstream momentum conservation, equation(11), can be rearranged to provide an equation for the turbulentstress: ∂t xy ∂ξ = ρ jet v j h (cid:48) ( x ) (cid:18) − ( ξ − ξ + 1) ξη + (1 − η )( ξ − ξ + 1) + D ( η, ξ ) (cid:19) .(25)Integration of equation (25) gives: t xy ( η, x, ξ ) = ρ jet v jet h (cid:48) ( x ) F ( η, ξ ) , (26)where the function F ( η, ξ ) is given in § A3. We set t xy ( ξ ) = t xy ( ξ ) = 0 since we expect the turbulence to be confined primar-ily to the region ξ < ξ < ξ . The condition t xy ( ξ ) = 0 is usedto compute the form of F ( η, ξ ) ( § A3); the condition t xy ( ξ ) = 0 ,and hence F ( ξ , η ) = 0 , allows us to calculate the position of thejet-mixing layer boundary, ξ , in ( x, ξ ) -space as a function of onlythe jet-to-wind density ratio η : ξ ( η ) = 2 η log( η ) + (4 − η ) η − η − for η (cid:54) = 1 , and = 13 for η = 1 . (27)The position of the mixing layer boundaries as a function of η isshown in Fig. 2. In the limit of a significantly underdense jet ( η → ), the mixing layer penetrates the jet and wind evenly. In the limitof a significantly overdense jet ( η → ∞ ), the mixing layer almostexclusively penetrates the ambient wind.Knowledge of the position of the mixing layer boundaries al-lows the forms of the transverse velocity and turbulent stress pro-files to be directly calculated as a function of position within themixing layer, ξ , and the jet-to-wind density ratio, η . The forms ofthese expressions are algebraically complex, and are given in § A4.The mixing layer density, transverse velocity and turbulent stressprofiles are shown in Fig. 3.As seen in Fig. 3(b), ambient wind material is pulled upwardsinto the mixing layer at the wind-mixing layer boundary ( ξ ) withan effective entrainment velocity, v ent = ˜ v y ( η, x, ξ ) (28) = v jet h (cid:48) ( x ) η ( η − η log( η ) − η − for η (cid:54) = 1 , and = v jet h (cid:48) ( x ) 16 for η = 1 . (29) MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows − N o r m a li ze dd e n s it y , ¯ ρ / ρ j e t (a) − . − . . . . . . N o r m a li ze d t r a n s v e r s e v e l o c it y , ˜ v y / ( v j e t h ( x )) (b) Density ratio, η − ξ ambientboundary ξ jetboundaryDimensionless coordinate ξ − . . . . . . . . N o r m a li ze d t u r bu l e n t s t r e ss , t x y / ( ρ j e t v j e t h ( x )) (c) Figure 3.
Normalized (a) density, (b) transverse velocity and (c) turbulentstress profiles for the mixing layer described by this model. The profiles areplotted between the mixing layer boundaries ξ and ξ . The numerical val-ues that these boundaries take are different for different jet-to-wind densityratios η , as per equation (27) and Fig. 2. This is equivalent to the entrainment velocity specified in the mod-els of Cant´o & Raga (1991) and Raga et al. (1995). However, it oc-curs naturally as a result of the boundary conditions of the problem,rather than being specified by an experimentally-determined ‘en-trainment efficiency’ parameter. We compare the entrainment ve-locity of our model to the earlier work in §§ v y / ( v jet h (cid:48) ( x )) , is plotted as a func-tion of η in Fig. 4(a). We define the mixing layer mass flux to be ˙ M ( x ) ≡ (cid:90) ξ ξ ρ ( ξ ) v jet U ( ξ ) h ( x ) d ξ . (30)The contribution to the mass flux from intercepted jet material issimply given by ˙ M jet ( x ) = ρ jet v jet h ( x ) ξ ( η ) . (31)(cf. Raga et al. 1995). Therefore, the entrained mass flux is ˙ M ent ( x ) ≡ (cid:90) ξ ξ ρ ( ξ ) v jet U ( ξ ) h ( x ) d ξ − ρ jet v jet h ( x ) ξ ( η ) . (32) . . . . . . N o r m a li ze d e n t r a i n m e n t v e l o c it y , v e n t / ( v j e t h ( x )) (a) . . . . . . N o r m a li ze d l o ca l r a t e - o f- c h a ng e o f m i x - i ng l a y e r m a ss fl ux , ˙ M / ( ρ j e t v j e t h ( x )) (b)TotalJetWind − − − Jet-to-wind density ratio, η . . . . . G ( η ) (c) Figure 4.
Mixing layer parameters which only depend on the jet-to-ambient-wind density ratio, η . (a) Normalized entrainment velocity fromthe ambient wind, from equation (29). (b) Normalized rate-of-change of themixing layer mass flux, ˙ M (cid:48) . The contribution to the rate-of-change of themixing layer mass flux from jet interception and wind entrainment ( ˙ M (cid:48) ent )is shown by the dashed and dot-dashed curves, respectively. (c) Dimension-less function G ( η ) , for the determination of the rate of turbulent energyproduction per unit area, from equation (39). The mass entrainment rate from the ambient wind is simply thederivative of equation (32) with respect to x . It can be shown ( § A5)that ∂ ˙ M ent ∂x ≡ ˙ M (cid:48) ent = ρ w v ent , (33)as expected, since ambient wind material is being drawn into themixing layer with velocity v ent ( § ∂ ˙ M/∂x ≡ ˙ M (cid:48) , is shown in Fig. 4(b). The contribution of wind entrainment tothe mass flux of the mixing layer is greatest for an underdense jet(Fig. 5). The ultimate aim of this model is to determine the rate of tur-bulent energy production and subsequent dissipation and radia-tion in the mixing layer. The rate of turbulent energy production, ˙ E turb = 2 t xy s xy , where s xy is the shear in the mixing layer (Kun-cic & Bicknell 2004). The mean shear may be calculated from s xy = 12 (cid:18) ∂ ˜ v x ∂y + ∂ ˜ v y ∂x (cid:19) ≈ ∂ ˜ v x ∂y (34)since the average transverse velocity, ˜ v y , varies slowly with respectto x . MNRAS , 000–000 (0000)
M. C. White et al. − − − Jet-to-ambient wind density ratio, η . . . . . . P e r ce n t a g ec on t r i bu ti on t o r a t e - o f- c h a ng e o f m i x i ng l a y e r m a ss fl ux Jet interceptionWind entrainment
Figure 5.
Percentage contribution to the rate-of-change of mixing layermass flux from jet interception (dashed curve) and ambient wind entrain-ment ( ˙ M (cid:48) ent , dot-dashed curve) as a function of jet-to-ambient-wind densityratio η . The rate of turbulent energy production in the mixing layeris directly comparable to the observed luminosity of the mixinglayer, assuming that the cooling time of the gas is short, so that theturbulent energy produced is radiated efficiently. This is a reason-able assumption in DG Tau; based on a gas temperature of Kand using the cooling function of Sutherland & Dopita (1993), wedetermine a cooling time for the IVC of 2.5 yr, which results ina cooling length ∼ au ≈ . arcsec projected distance for aflow speed of km s − . This cooling length is short comparedto the length of the mixing layer, which is (cid:38) au ( § ˙ E turb ( x, ξ ) d V = t xy ( x, ξ ) ∂ ˜ v x ∂y = t xy ( x, ξ ) v jet h ( x ) (35) = ρ jet v jet h (cid:48) ( x ) h ( x ) F ( η, ξ ) . (36)We integrate over y to form an expression for the turbulent energy Kuncic & Bicknell (2004) identified the various terms in the internal en-ergy equation relating to the advection of internal plus turbulent energy,the generation of turbulent energy, and its emission via radiation. Since thecooling length in the mixing layer is much less than the advective length,this assumption effectively approximates the energy equation by equatingthe production of turbulent energy to the emission of radiation. This is sim-ilar to the approach taken in classical accretion disc theory (e.g. Shakura &Sunyaev 1973). produced per unit area:d ˙ E turb d A ( x ) = (cid:90) y y d ˙ E turb ( x, ξ ) d V d y (37) = ρ jet v jet h (cid:48) ( x ) G ( η ) , where (38) G ( η ) = 2 η + 3 η − η log( η ) − η + 112( η − for η (cid:54) = 1 , and = 124 for η = 1 . (39)For the purpose of comparing this model to observations of three-dimensional protostellar jet mixing layers, this is the turbulent en-ergy produced per unit circumference per unit length. Therefore,total turbulent energy production in the layer is calculated multi-plying equation (38) by πR mix ( x ) , where R mix ( x ) is the mixinglayer radius, and then integrating over the observed mixing layerlength, L : ˙ E tot = (cid:90) L πR mix ( x ) ρ jet v jet h (cid:48) ( x ) G ( η ) d x (40) = 2 πR mix Lρ jet v jet h (cid:48) ( x ) G ( η ) , (41)assuming that R mix and h (cid:48) ( x ) are independent of x . II ] 1.644 µ m Shock Modelling Our model provides estimates for the bolometric luminosity of aprotostellar outflow mixing layer ( § II ] line luminosities for a given total luminosity. Tothis end, a grid of shock models capable of heating their post-shockgas to between × K and × K were computed usingthe
MAPPINGS
IV code, version 4.0.1 (Sutherland & Dopita 1993;Allen et al. 2008; Nicholls et al. 2012, 2013; Dopita et al. 2013),covering both solar abundances (Asplund et al. 2009) and iron-depleted abundances (Jenkins 2009, 2013). We use this grid as arepresentative model of partially-ionized gas being heated and sub-sequently cooling; we are not concerned with the shock structureitself. Therefore, we chose pre-shock gas parameters, summarisedin Appendix C1, that yield densities and temperatures in the post-shock region that are comparable to those expected in the mixinglayer. These pre-shock parameters are not intended to be strictlyrepresentative of protostellar outflows, nor do we imply that theIVC emission is generated in a single flat-planar shock structure,which would not be a good approximation to the shocks occurringin a turbulent mixing layer.The shocks were driven into a pre-shock medium with a hy-drogen number density of cm − , which approximately matchesthe conditions in the DG Tau jet ( § MAPPINGS
IV collisional ioniza-tion excitation model. Selecting a value of plasma β , that is, the ra-tio of thermal to magnetic pressure in the pre-shock material, and apost-shock temperature fixes the shock velocity. For the post-shocktemperatures given above, the resulting shock velocity is between and km s − for β = 1 , in agreement with observations (Pa-per I). This computed shock velocity was then used to calculate theproperties of the post-shock, cooling gas. The final [Fe II ] 1.644 MNRAS000
IV collisional ioniza-tion excitation model. Selecting a value of plasma β , that is, the ra-tio of thermal to magnetic pressure in the pre-shock material, and apost-shock temperature fixes the shock velocity. For the post-shocktemperatures given above, the resulting shock velocity is between and km s − for β = 1 , in agreement with observations (Pa-per I). This computed shock velocity was then used to calculate theproperties of the post-shock, cooling gas. The final [Fe II ] 1.644 MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows − − − Ionization fraction, χ . . . [ F e II ] . m l u m i no s it y fr ac ti on , L . / L t o t × − Pre-shock plasma β × K)2.0 3.0 4.0 5.0 6.0Post-shock temperature ( × K)2.0 3.0 4.0 5.0 6.0
Figure 6.
Contribution of [Fe II ] 1.644 µ m line emission to total shockemission from a grid of MAPPINGS
IV models, where the iron abundance isnot depleted from solar abundance values and the pre-shock number density n H = 10 . The percentage of shock emission which is radiated as [Fe II ]1.644 µ m line emission is shown as a function of the pre-shock ionizationfraction, χ . Models have been computed for a range of post-shock tempera-tures and pre-shock values of plasma β as indicated. Points are joined usinga one-dimensional Akima spline interpolation (Akima 1970). The verticaldotted line denotes the known DG Tau jet ionization, χ = 0 . . µ m line emission from the shock is then expressed as a fractionof the total of all the line emission plus the two-photon emission,which can be a large contributor to total emission in these heated,partially-ionized models.The results from this model grid are shown in Fig. 6. Acrossa range of pre-shock ionizations and plasma β (a convenient proxyfor magnetic field strength), the fraction of shock luminosity emit-ted in the [Fe II ] 1.644 µ m line is ∼ − . For the purposes ofour model, this order-of-magnitude will suffice for comparing withobservations.The [Fe II ] total emission affects the model structure, so thatthe [Fe II ] 1.644 µ m luminosity does not simply scale with thegas phase abundance of iron, but nearly so (see Appendix C2 fordetails). The cooling fraction in [Fe II ] 1.644 µ m with undepletedgas is (cid:46) × − (Fig. 6), and between × − and × − forthe depleted models (Fig. 7). It is therefore reasonable that the [Fe II ] 1.644 µ m line emission as a fraction of total cooling emissionlies between − and − for a range of iron depletion factors.We parametrize the ratio between bolometric luminosity, L tot , and[Fe II ] 1.644 µ m luminosity, L . , as L . = (cid:18) − Fe depletion factor (cid:19) × L tot . (42)
We summarize the parameters of the DG Tau jet required as inputsto our model here. The main parameters are the jet density, ρ jet Note that zero depletion corresponds to a depletion factor of 1. − − − Ionization fraction, χ [ F e II ] . m l u m i no s it y fr ac ti on , L . / L t o t × − Pre-shock plasma β × K)2.0 3.0 4.0 5.0 6.0Post-shock temperature ( × K)2.0 3.0 4.0 5.0 6.0
Figure 7.
As for Fig. 6, but with the pre-shock gas having an iron depletionfactor of 100. ( § v jet ( § h (cid:48) ( x ) ( § L ( § η ( § MAPPINGS
IV models, we also need to know the iron depletion factor in theDG Tau IVC ( § The electron density may be determined from NIR observationsthrough the ratio of the [Fe II ] emission lines at . µ m and . µ m. Pesenti et al. (2003) computed a relationship betweenthis line ratio and electron density for a 16-level model of an Fe + atom. The similar BE99 technique makes use of the ratio of the[S II] emission lines at ˚A and ˚A in the optical regime(Bacciotti & Eisl¨offel 1999; Maurri et al. 2014). These techniqueshave been applied to the DG Tau jet (Paper 1; Bacciotti et al. 2000;Agra-Amboage et al. 2011; Maurri et al. 2014), yielding a typi-cal electron density n e ∼ cm − . Maurri et al. (2014) reportedhigher electron densities, up to cm − , close to the central star.However, we do not observe significant IVC emission at this posi-tion (Paper I, fig. 2.6 therein), so that we use the lower jet densitycorresponding to the region where we observe a mixing layer.Our model requires the jet mass density as an input. We con-vert the measured electron density into physical density as in PaperI, § n H = n e /χ e , where χ e is the ionization fraction of the gas. Although theionization fraction of the jet appears to vary with position (Maurriet al. 2014), an average ionization fraction of χ e = 0 . ± . is areasonable approximation (Bacciotti et al. 2000). This yields a hy-drogen number density n H = 3 . × cm − . The mass density ρ = 1 . mn H for a gas consisting of 90 per cent hydrogen and 10per cent helium, where m is the atomic mass unit. This calculationleads to a jet mass density ∼ − g cm − , which we use as afiducial value for our model. MNRAS , 000–000 (0000)
M. C. White et al. . . . . . . . . . . D i a m e t e r ( ) ± ± D i a m e t e r( AU )
50 100 150 200 250 300
Deprojected distance from DG Tau (AU) . . . . . . . x ( )0 . . . . . M i x i ng l a y e r w i d t h ( ) (b) M i x i ng l a y e r w i d t h , h ( x ) ( AU ) Figure 8.
Growth rates of the diameter of the approaching DG Tau outflowcomponents. (a) Diameter of the approaching jet (circles) and IVC (trian-gles) as a function of distance from the central star. Jet diameters were mea-sured using cross-jet Gaussian fits to the high-velocity component inten-sity of [Fe II ] 1.644 µ m, and are approximately deconvolved from the PSF(Paper I). IVC diameters were measured from direct inspection of imagesand cross-outflow intensity profiles of the [Fe II ] 1.644 µ m intermediate-velocity component. Linear fits to the growth of both components over theregion . – . arcsec from the central star are shown as solid lines; thelabels indicate the slope of the line (i.e. the growth rate of the componentdiameter). (b) Inferred mixing layer widths over the region . – . arcsecfrom the central star, as per equation (43). The velocity of the DG Tau jet varies with time (Paper I; Bac-ciotti et al. 2002; Pyo et al. 2003; Agra-Amboage et al. 2011),and this is the likely cause of the observed moving shock-excitedknots (e.g. Raga et al. 1990). The jet velocity is typically mea-sured from the high-velocity peak of line emission (Pyo et al. 2003;Agra-Amboage et al. 2011). More recently, we used multicompo-nent Gaussian fitting, coupled with a statistical F -test, to rigorouslyseparate the two [Fe II ] 1.644 µ m line-emission components in theapproaching DG Tau outflow (Paper I). These fits show that thehigh-velocity component of the outflow has a range of velocitiesfrom – km s − in the 2005 observing epoch. Therefore, weadopt an average jet velocity of km s − for use in our model. The deprojected length of the observed mixing layer, . arcsec ≈ au, can be directly measured from our data (Paper I). Theseparation of the line emission from the two approaching outflowcomponents allows for the calculation of the mixing layer growthrate. The diameter of the HVC (jet), D jet , was determined by fittinga Gaussian to the [Fe II ] 1.644 µ m line emission in the cross-jet di-rection, and approximately deconvolving the width of this Gaussianfrom the PSF via the formula D jet = FWHM obs − FWHM PSF (Pa-per I, § D IVC , was measuredfrom both an image of that component (Paper I, fig. 2.6d therein)and cross-outflow [Fe II ] 1.644 µ m IVC profiles at each down-stream position. IVC diameters could not be reliably determined beyond ∼ arcsec from the central star, due to incomplete linefitting coverage in this region, although conservative lower limitscould be inferred. The component diameters are shown as a func-tion of downstream position in Fig. 8(a). The inferred mixing layerwidth is simply the difference between the observed radii , r jet and r IVC , of the jet and IVC, h ( x ) = r IVC − r jet = D IVC − D jet . (43)The growth rate of the mixing layer is then h (cid:48) ( x ) = ( D (cid:48) IVC − D (cid:48) jet ) / . (44)We determine the growth rate of the mixing layer as fol-lows. We construct linear fits to the lateral growth of both thejet and IVC in the approaching DG Tau outflow in the 2005 ob-serving epoch over the region . – . arcsec from the central star(Fig. 8a). These fits give growth rates of D (cid:48) IVC = 0 . ± . , and D (cid:48) jet = 0 . ± . . From equation (44), the measured growth ratesimply a mixing layer growth rate of . ± . .The inferred mixing layer width as a function of distance fromthe central star is shown in Fig. 8(b). This is a noisier profile thanthe individual jet diameters; therefore, it is preferable to determine h (cid:48) ( x ) from equation (44). The density of the jet is well-defined (see § au along the outflow axis following deprojection, and au across the outflow direction, resulting in a total wind open-ing angle of ◦ . By considering the K -band extinction towardsDG Tau, and the ratio of H µ m emitting mass tototal H mass, they determined a minimum total wind mass in thisregion of . × − M (cid:12) . This corresponds to a minimum averagewind H number density of × cm − , assuming a filling factorof 1 and a conical geometry.We now make approximations about the flow geometry of thiswind in order to determine its density in the entrainment regionof the outflow. Consider a distance . arcsec ≈ au from thecentral star, which is halfway along the observed mixing layer. Ifthe wind undergoes no further collimation beyond what is observedin H emission, and maintains a conical geometry, it will have atotal width of au at this position. Assuming a constant windmass-loss rate, ˙ M , and wind velocity, v , the wind density, ρ , isinversely proportional to the wind radius, R , squared: ˙ M = ρπR v = const. ⇒ ρ ρ = R R . (45)Therefore, an increase in wind radius of a factor . would meana decrease in wind density of a factor of ∼ , resulting in an H number density by ∼ × cm − at au from the centralstar. Assuming a gas composition of 90 per cent hydrogen and 10per cent helium by number density, this results in a mass densityof . × − g cm − , and a jet-to-ambient-wind density ratio of12.2. The fits were made using the deprojected distance from the central starand the physical diameter of each outflow component, thereby accountingfor the projection of the DG Tau outflows to the line-of-sight ( ◦ ; Eisl¨offel& Mundt 1998). MNRAS000
Growth rates of the diameter of the approaching DG Tau outflowcomponents. (a) Diameter of the approaching jet (circles) and IVC (trian-gles) as a function of distance from the central star. Jet diameters were mea-sured using cross-jet Gaussian fits to the high-velocity component inten-sity of [Fe II ] 1.644 µ m, and are approximately deconvolved from the PSF(Paper I). IVC diameters were measured from direct inspection of imagesand cross-outflow intensity profiles of the [Fe II ] 1.644 µ m intermediate-velocity component. Linear fits to the growth of both components over theregion . – . arcsec from the central star are shown as solid lines; thelabels indicate the slope of the line (i.e. the growth rate of the componentdiameter). (b) Inferred mixing layer widths over the region . – . arcsecfrom the central star, as per equation (43). The velocity of the DG Tau jet varies with time (Paper I; Bac-ciotti et al. 2002; Pyo et al. 2003; Agra-Amboage et al. 2011),and this is the likely cause of the observed moving shock-excitedknots (e.g. Raga et al. 1990). The jet velocity is typically mea-sured from the high-velocity peak of line emission (Pyo et al. 2003;Agra-Amboage et al. 2011). More recently, we used multicompo-nent Gaussian fitting, coupled with a statistical F -test, to rigorouslyseparate the two [Fe II ] 1.644 µ m line-emission components in theapproaching DG Tau outflow (Paper I). These fits show that thehigh-velocity component of the outflow has a range of velocitiesfrom – km s − in the 2005 observing epoch. Therefore, weadopt an average jet velocity of km s − for use in our model. The deprojected length of the observed mixing layer, . arcsec ≈ au, can be directly measured from our data (Paper I). Theseparation of the line emission from the two approaching outflowcomponents allows for the calculation of the mixing layer growthrate. The diameter of the HVC (jet), D jet , was determined by fittinga Gaussian to the [Fe II ] 1.644 µ m line emission in the cross-jet di-rection, and approximately deconvolving the width of this Gaussianfrom the PSF via the formula D jet = FWHM obs − FWHM PSF (Pa-per I, § D IVC , was measuredfrom both an image of that component (Paper I, fig. 2.6d therein)and cross-outflow [Fe II ] 1.644 µ m IVC profiles at each down-stream position. IVC diameters could not be reliably determined beyond ∼ arcsec from the central star, due to incomplete linefitting coverage in this region, although conservative lower limitscould be inferred. The component diameters are shown as a func-tion of downstream position in Fig. 8(a). The inferred mixing layerwidth is simply the difference between the observed radii , r jet and r IVC , of the jet and IVC, h ( x ) = r IVC − r jet = D IVC − D jet . (43)The growth rate of the mixing layer is then h (cid:48) ( x ) = ( D (cid:48) IVC − D (cid:48) jet ) / . (44)We determine the growth rate of the mixing layer as fol-lows. We construct linear fits to the lateral growth of both thejet and IVC in the approaching DG Tau outflow in the 2005 ob-serving epoch over the region . – . arcsec from the central star(Fig. 8a). These fits give growth rates of D (cid:48) IVC = 0 . ± . , and D (cid:48) jet = 0 . ± . . From equation (44), the measured growth ratesimply a mixing layer growth rate of . ± . .The inferred mixing layer width as a function of distance fromthe central star is shown in Fig. 8(b). This is a noisier profile thanthe individual jet diameters; therefore, it is preferable to determine h (cid:48) ( x ) from equation (44). The density of the jet is well-defined (see § au along the outflow axis following deprojection, and au across the outflow direction, resulting in a total wind open-ing angle of ◦ . By considering the K -band extinction towardsDG Tau, and the ratio of H µ m emitting mass tototal H mass, they determined a minimum total wind mass in thisregion of . × − M (cid:12) . This corresponds to a minimum averagewind H number density of × cm − , assuming a filling factorof 1 and a conical geometry.We now make approximations about the flow geometry of thiswind in order to determine its density in the entrainment regionof the outflow. Consider a distance . arcsec ≈ au from thecentral star, which is halfway along the observed mixing layer. Ifthe wind undergoes no further collimation beyond what is observedin H emission, and maintains a conical geometry, it will have atotal width of au at this position. Assuming a constant windmass-loss rate, ˙ M , and wind velocity, v , the wind density, ρ , isinversely proportional to the wind radius, R , squared: ˙ M = ρπR v = const. ⇒ ρ ρ = R R . (45)Therefore, an increase in wind radius of a factor . would meana decrease in wind density of a factor of ∼ , resulting in an H number density by ∼ × cm − at au from the centralstar. Assuming a gas composition of 90 per cent hydrogen and 10per cent helium by number density, this results in a mass densityof . × − g cm − , and a jet-to-ambient-wind density ratio of12.2. The fits were made using the deprojected distance from the central starand the physical diameter of each outflow component, thereby accountingfor the projection of the DG Tau outflows to the line-of-sight ( ◦ ; Eisl¨offel& Mundt 1998). MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows Takami et al. (2004) notes that their estimates of H mass anddensity in the wider-angle wind are lower limits, given that cold gasmay be present in the outflow, and the filling factor of the wind maybe less than unity. Hence, our estimate of the jet-to-ambient winddensity ratio represents an upper limit to possible values for this pa-rameter. Therefore, whilst we consider it likely that η lies between1 and 10, we have investigated a parameter range of . (cid:54) η (cid:54) for completeness. The iron depletion in the approaching outflow components fromDG Tau was measured by Agra-Amboage et al. (2011) via compar-ison of [Fe II ] 1.644 µ m flux to the [O I] 6300 ˚A fluxes reportedby Lavalley-Fouquet et al. (2000). Through comparison with shockwave models (Hartigan et al. 2004), they determined that the irondepletion factor in the approaching DG Tau outflow is ∼ – in gasfaster than − km s − (the jet), and ∼ – for gas at speedsbelow − km s − (the IVC).Agra-Amboage et al. (2011) noted that their measurements aretentative, given that they were required to compare [Fe II ] and [OI] line fluxes obtained ∼ yr apart. Furthermore, it is unlikely thatthe DG Tau jet would exhibit any iron depletion, as we would ex-pect dust grains to be destroyed by passage through the strong rec-ollimation shock at the base of the approaching outflow (Paper I).However, it is reasonable that the slower, wider-angle outflow com-ponents would exhibit iron depletion, as they are launched fromwider disc radii and may be less shock-processed (Agra-Amboageet al. 2011). Indeed, higher depletion at lower flow velocities hasbeen observed in other YSOs (e.g. calcium in the HH 111 outflow;Podio et al. 2009). The jet iron depletion measurement of Agra-Amboage et al. may be contaminated by IVC emission, given thoseauthors made a simple velocity cut to separate outflow components,rather than using line component fitting (e.g. Paper I; Lavalley et al.1997). Therefore, we take a range of iron depletion factors, – ,for the DG Tau approaching IVC. We now compare the estimates from our model to our previousobservations of the approaching DG Tau outflow intermediate-velocity component (Paper I). These estimates are based on theoutflow parameters for DG Tau detailed above ( § The estimated mixing layer luminosity for DG Tau from our modelis shown in Fig. 9. We compute the total mixing layer bolometricluminosity as per § L mix,tot = 2 πR mix Lρ jet v jet h (cid:48) ( x ) G ( η ) . (46)We estimate the mixing layer radius, R mix , to be ∼ au (Fig. 8).The mixing layer luminosity is estimated using a jet velocity of km s − ( § − g cm − ( § h (cid:48) ( x ) =0 . ± . ( § §§ II ] 1.644 µ m luminosityfrom the approaching DG Tau outflow as follows. We consider ev-ery spaxel covering the approaching outflow that was successfully − Jet-to-ambient-wind density ratio, η E s ti m a t e d m i x i ng l a y e r [ F e II] . m l u m i no s it y ( e r g s − ) zero310100 DG Tau [Fe II]1.644 µm IVCluminosity, . × erg s − Figure 9.
Estimates of the DG Tau mixing layer [Fe II ] 1.644 µ m luminos-ity, from equations (42) and (46). Luminosities are calculated for a rangeof iron depletion factors and jet-to-ambient-wind density ratios, assuming ajet velocity of km s − and a jet density of . × − g cm − . Solidcurves show the estimated luminosity for h (cid:48) ( x ) = 0 . ; the surroundinggreyed regions indicate the range of luminosities at a given iron depletionfactor for . (cid:54) h (cid:48) ( x ) (cid:54) . ( § ∼ – ; (cid:46) η (cid:46) ). The thick dashed line shows the observed [Fe II ] 1.644 µ mluminosity of the DG Tau approaching IVC, . × erg s − . fit with two [Fe II ] 1.644 µ m emission-line components (Paper I,fig. 2.6 therein). We then calculate the flux from the fitted IVC com-ponent in each spaxel, and sum across the entire outflow to producea total IVC [Fe II ] 1.644 µ m luminosity of . × erg s − , as-suming a distance to DG Tau of pc (Elias 1978).Our model estimates for the luminosity of a mixing layer inthe approaching DG Tau outflow is in good agreement with ourobservations of the approaching IVC. For (cid:46) η (cid:46) ( § ∼ – , our model estimates a mixinglayer luminosity of (1 . – . × erg s − . This is a good levelof agreement between model and observations, and constitutes astrong indicator that the luminosity of this region of the outflows isdriven by turbulent dissipation. In our model, mass enters the mixing layer from both the jet via in-terception, and from the ambient wind via entrainment. From equa-tion (30), the rate at which material enters the mixing layer is ∂ ˙ M∂x ≡ ˙ M (cid:48) = ρ jet v jet h (cid:48) ( x ) (cid:18) − η + η log( η ) + 1( η − (cid:19) . (47)Multiplying by πR mix , where R mix ≈ au as per § L = 270 au ( § MNRAS , 000–000 (0000) M. C. White et al. sitions (cf. the calculation of the total turbulent energy productionin the mixing layer in §§ (cid:54) η (cid:54) , § (3 . – . × − M (cid:12) yr − for the observed mixing layer(Fig. 10a). By comparison, the mass-loss rate of the DG Tau jet is ∼ × − M (cid:12) yr − (Paper I; Agra-Amboage et al. 2011) from[Fe II ] emission-line ratios; the total mass-loss rate of all ionizedoutflow components (jet plus IVC) is (1 – × − M (cid:12) yr − fromthe VLA data of Lynch et al. (2013). The mass-loss rate of themolecular wind is lower, (cid:38) . × − M (cid:12) yr − (Takami et al.2004). However, in the overdense-jet regime, mass interceptionfrom the jet is the main contributor to the mass within the mixinglayer (Fig. 5). Therefore, we conclude that the total mass gain intothe mixing layer estimated by our model is less than the combinedmass-loss rates of the DG Tau jet and molecular wind, as requiredfor consistency.Recently, Maurri et al. (2014) performed an analysis of theDG Tau approaching outflow using the BE99 technique for deter-mining physical flow parameters from optical line ratios (Bacciotti& Eisl¨offel 1999). They found that, over the first . arcsec ofthe approaching outflow, the mass outflow rate of the jet (identi-fied as the high-velocity interval, or HVI, in their paper) decreasedby ∼ . dex. Over the same region, the mass outflow rate of themedium-velocity interval (MVI, which is comparable to our IVC)increased. This is what would be observed if the IVC/MVI repre-sents a turbulent mixing layer which is primarily gaining materialfrom the central jet/HVI.We compare the observations of Maurri et al. (2014) to ourmodel estimates for the rate-of-change of mixing layer mass fluxper unit length, and find them to be consistent. We performed a lin-ear fit to the MVI mass-loss rates of Maurri et al. (2014, fig. 15therein), and determine an increase in MVI mass-loss rate per unitlength of × − M (cid:12) au − yr − . For the observed parametersof the DG Tau outflow ( § (1 . – . × − M (cid:12) au − yr − , with the lower rate-of-change corresponding to a more overdense jet (Fig. 10a). The es-timates from our model strongly suggest that the approaching IVCof the DG outflows is consistent with being the signature of a tur-bulent mixing layer around the central jet.
In the jet entrainment models of Cant´o & Raga (1991) and Ragaet al. (1995), material is injected into the mixing layer from theambient medium/wind with a prescribed entrainment velocity. Thisvelocity is expressed as a fraction of the sound speed in the ambientwind, c w , as it was argued that the ambient wind would be incapableof supplying material at a velocity greater than the sound speed(Cant´o & Raga 1991, however, see § For (cid:54) η (cid:54) , mass entrainment from the ambient wind contributes20–33 per cent of the rate-of-change of mixing layer mass flux (Fig. 5), sothe mass entrainment per unit length from the ambient wind is (0 . – . × − M (cid:12) au − yr − . . . . . . . . . R a t e - o f- c h a ng e o f m i x i ng l a y e r m a ss fl uxp e r un it a r ea ( − g c m − s − ) (a) TotalJetWind T o t a l r a t e - o f- c h a ng e o f m i x i ng l a y e r m a ss fl uxp e r un itl e ng t h ( − M (cid:12) AU − y r − ) − Jet-to-ambient-wind density ratio, η E n t r a i n m e n t v e l o c it y ( k m s − ) (b) Figure 10.
Theoretical estimates for the DG Tau mixing layer. (a) Esti-mated rate-of-change of the mixing layer mass flux (solid line), from equa-tion (47); contributions from the jet (dashed curve) and ambient wind (dot-dashed curve) are also shown. Rate-of-change of mixing layer mass flux perunit length (right-hand axis) is calculated assuming a mixing layer radius of25 au ( § km s − , a jet density of − g cm − , and a mixing layer growth rate of . . Greyed regionsin both panels show the estimated parameters for a range of mixing layergrowth rates, . (cid:54) h (cid:48) ( x ) (cid:54) . . as the entrainment efficiency, (cid:15) (cid:54) , where the entrainment velocityis written as v ent = (cid:15)c w .The entrainment velocity of our model is given by equation(29), and is shown as a dimensionless function of η in Fig. 4(a).For jet-to-ambient-wind density ratios − (cid:54) η (cid:54) , the di-mensionless entrainment velocity, v ent / ( v jet h (cid:48) ( x )) varies between ∼ . and ∼ . . For a jet velocity of km s − and h (cid:48) ( x ) = 0 . – . , we estimate a range of entrainment veloci-ties, . (cid:54) v ent (cid:54) . km s − . For our inferred values of η (cid:46) and h (cid:48) ( x ) = 0 . , equation (29) estimates an entrainment velocity (cid:46) km s − (Fig. 10b).Assuming LTE, the H µ m emission observed in theapproaching DG Tau outflow has a temperature of × K (Becket al. 2008; Agra-Amboage et al. 2014). As the ambient wind is notdirectly observable in the region where entrainment is occurring,we may assume that the wind has cooled somewhat. Therefore, wetake an indicative temperature of K, which leads to a soundspeed of . km s − in the wind for a molecular gas with meanmolecular weight . m H .Direct comparison with the range of entrainment velocitiespredicted above implies a range of entrainment efficiencies be-tween 0.13 and 3.05 for the full range of possible values forthe mixing layer growth rate and jet-to-ambient wind density ra-tio. Adopting the best-fit value for the mixing layer growth rate, h (cid:48) ( x ) = 0 . , and assuming that the jet is likely to be overdenseby up to a factor of 10 ( § § MNRAS000
Theoretical estimates for the DG Tau mixing layer. (a) Esti-mated rate-of-change of the mixing layer mass flux (solid line), from equa-tion (47); contributions from the jet (dashed curve) and ambient wind (dot-dashed curve) are also shown. Rate-of-change of mixing layer mass flux perunit length (right-hand axis) is calculated assuming a mixing layer radius of25 au ( § km s − , a jet density of − g cm − , and a mixing layer growth rate of . . Greyed regionsin both panels show the estimated parameters for a range of mixing layergrowth rates, . (cid:54) h (cid:48) ( x ) (cid:54) . . as the entrainment efficiency, (cid:15) (cid:54) , where the entrainment velocityis written as v ent = (cid:15)c w .The entrainment velocity of our model is given by equation(29), and is shown as a dimensionless function of η in Fig. 4(a).For jet-to-ambient-wind density ratios − (cid:54) η (cid:54) , the di-mensionless entrainment velocity, v ent / ( v jet h (cid:48) ( x )) varies between ∼ . and ∼ . . For a jet velocity of km s − and h (cid:48) ( x ) = 0 . – . , we estimate a range of entrainment veloci-ties, . (cid:54) v ent (cid:54) . km s − . For our inferred values of η (cid:46) and h (cid:48) ( x ) = 0 . , equation (29) estimates an entrainment velocity (cid:46) km s − (Fig. 10b).Assuming LTE, the H µ m emission observed in theapproaching DG Tau outflow has a temperature of × K (Becket al. 2008; Agra-Amboage et al. 2014). As the ambient wind is notdirectly observable in the region where entrainment is occurring,we may assume that the wind has cooled somewhat. Therefore, wetake an indicative temperature of K, which leads to a soundspeed of . km s − in the wind for a molecular gas with meanmolecular weight . m H .Direct comparison with the range of entrainment velocitiespredicted above implies a range of entrainment efficiencies be-tween 0.13 and 3.05 for the full range of possible values forthe mixing layer growth rate and jet-to-ambient wind density ra-tio. Adopting the best-fit value for the mixing layer growth rate, h (cid:48) ( x ) = 0 . , and assuming that the jet is likely to be overdenseby up to a factor of 10 ( § § MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows Cant´o & Raga (1991) and Raga et al. (1995) utilized laboratoryexperiments (Birch & Eggers 1972) to estimate the entrainment ef-ficiency, (cid:15) , of protostellar jet mixing layers. Raga et al. proposethat (cid:15) ∼ . . Furthermore, both Cant´o & Raga (1991) and Ragaet al. (1995) claimed that (cid:15) (cid:54) , because the ambient wind shouldnot supply material at greater than the sound speed. However, ourmodel, based only on observable parameters of the protostellar out-flows, implies an entrainment efficiency (cid:46) (cid:15) (cid:46) . , in contradic-tion to earlier models. We argue below that our estimated entrain-ment efficiency is physically viable.A detailed analysis of the laboratory experiments of Birch &Eggers (1972) is beyond the scope of this paper. However, we makeseveral important points. First, the experiments of Birch & Eggersconcern adiabatic mixing layers. However, both observations (e.g.Bacciotti et al. 2002) and analytical estimates of the mixing layercooling length ( § § v ent . This ve- locity is supersonic with respect to the ambient wind sound speed,and probably transonic with respect to the shear layer sound speed,so shocks will result; however, these shocks will be standard tur-bulence internal to the mixing layer, and are therefore consistentwith our general approach. There will be no global blunt body-typeshock along the shear layer boundary.In reality, the transition from non-turbulent flow outside themixing layer to fully turbulent flow within would not be as abruptas we have modelled here. In particular, the sudden increase in y -velocity at the shear layer-ambient wind boundary is likely an arte-fact from our adoption of a linear shear layer temperature profile;the transition from low-density entraining molecular gas to high-density jet gas is likely to be more gradual. However, even if thereare comparable turbulent velocities within the molecular gas beingdrawn into the layer as in the gas well within the mixing layer, therate of dissipation per unit volume, ∼ ρv (cid:48) /l t , where ρ is the den-sity, v (cid:48) is the turbulent velocity, and l t is the turbulent eddy scalesize, is lower in the molecular gas than in the mixing layer, becauseof the lower density of the former. Jets that undergo lateral entrainment will eventually become com-pletely turbulent, as the inner boundary of the mixing layer expandsinto the jet and reaches the symmetry axis (e.g. Bicknell 1984;Dash et al. 1985). This does not appear to occur in the DG Tau jetwithin . arcsec ∼ au of the central star, as is evidenced bythe low-velocity-dispersion core of the approaching high-velocity[Fe II ] 1.644 µ m component (Paper I, fig. 2.6c therein). It is there-fore relevant to determine if our model predicts the DG Tau jetshould remain laminar within the NIFS field. Whilst a fully three-dimensional, axisymmetric model is formally required to make thiscalculation (Cant´o & Raga 1991), our model provides a useful pre-liminary exploration.The jet will become totally turbulent once the jet-mixing layerboundary, y , reaches the symmetry axis of the jet. The jet-mixinglayer boundary position is given by y ( x ) = h ( x ) ξ ( η ) = h (cid:48) ( x ) ξ ( η ) x , (48)assuming h ( x ) is linear in x . The downstream distance at whichthe jet becomes completely turbulent, x turb , is then simply x turb = r jet h (cid:48) ( x ) ξ ( η ) → x turb r jet = (cid:0) h (cid:48) ( x ) ξ ( η ) (cid:1) − . (49)For DG Tau, h (cid:48) ( x ) = 0 . ( § ξ ( η = 10) ∼ . ( §§ r jet,max ∼ au(Fig. 8a). The distance from the central star where the DG Taujet becomes completely turbulent is then ∼ au ≈ . arc-sec along the outflow axis, accounting for projection effects. Thisis well beyond the extent of the NIFS field, in agreement with ourearlier observations (Paper I). We have constructed a model of the turbulent lateral entrainment ofambient material by a supersonic, collimated jet ( § MNRAS , 000–000 (0000) M. C. White et al. turbulent mixing layer between the jet and the surrounding molec-ular wind, via calculation of the total turbulent stress within thelayer. This allows theoretical estimates of, e.g. the luminosity andentrainment rate of the mixing layer to be formed.We computed estimates for the bulk properties of the [Fe II ]1.644 µ m IVC observed in the approaching outflow from the YSODG Tauri ( § MAPPINGS
IV code, to facilitate comparison between the observed[Fe II ] luminosity of the component, and the estimated bolomet-ric luminosity from our model. Our model accurately estimates theluminosity and rate-of-change of mass flux of the DG Tau IVC,leading us to conclude that the IVC does indeed represent a turbu-lent mixing layer between the DG Tau high-velocity jet, and wider-angle disc wind.We compared our work with previous models of turbulent en-trainment by jets, specifically those of Cant´o & Raga (1991) andRaga et al. (1995). Unlike the previous models, our adoption ofan alternative semi-empirical approach means our model is not de-pendent upon an ‘entrainment efficiency’ parameter, which must beestimated from laboratory experiments. We argued that the require-ment for subsonic ‘entrainment velocities’ from the ambient windis not necessary in the context of our model. We also estimatedthe extent of laminar jet flow in DG Tau ( § ACKNOWLEDGEMENTS
Based on observations obtained at the Gemini Observatory, whichis operated by the Association of Universities for Research in As-tronomy, Inc., under a cooperative agreement with the NSF on be-half of the Gemini partnership: the National Science Foundation(United States), the National Research Council (Canada), CONI-CYT (Chile), the Australian Research Council (Australia), Min-ist´erio da Ciˆencia, Tecnologia e Inovac¸˜ao (Brazil) and Ministeriode Ciencia, Tecnolog´ıa e Innovaci´on Productiva (Argentina).We are extremely grateful for the support of the NIFS teamsat the Australian National University, Auspace, and Gemini Obser-vatory for their tireless efforts during the instrument integration,commissioning and system verification: Jan Van Harmleen, PeterYoung, Mark Jarnyk (deceased), Nick Porecki, Richard Gronke,Robert Boss, Brian Walls, Gelys Trancho, Inseok Song, ChrisCarter, Peter Groskowski, Tatiana Paz, John White, and JamesPatao. Thanks must also go to Tracy Beck for her support duringthe commissioning and observations.MW acknowledges the generous travel support fromAcademia Sinica to attend the conference Star Formation ThroughSpectroimaging at High Angular Resolution in July 2011, whichprovided useful information for this study. This work was sup-ported by the Australian Research Council through DiscoveryProject Grant DP120101792 (R. Salmeron).
REFERENCES
Agra-Amboage V., Dougados C., Cabrit S., Reunanen J., 2011, A&A, 532,A59Agra-Amboage V., Cabrit S., Dougados C., Kristensen L. E., Ibgui L., Re-unanen J., 2014, A&A, 564, A11Ainsworth R. E., Scaife A. M. M., Ray T. P., Taylor A. M., Green D. A.,Buckle J. V., 2014, ApJ, 792, L18Akima H., 1970, J. ACM, 17, 589Allen M. G., Groves B. A., Dopita M. A., Sutherland R. S., Kewley L. J.,2008, ApJS, 178, 20Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ARAA, 47, 481Bacciotti F., Eisl¨offel J., 1999, A&A, 342, 717Bacciotti F., Mundt R., Ray T. P., Eisl¨offel J., Solf J., Camezind M., 2000,ApJ, 537, L49Bacciotti F., Ray T. P., Mundt R., Eisl¨offel J., Solf J., 2002, ApJ, 576, 222Bachiller R., 1996, ARAA, 34, 111Beck T. L., McGregor P. J., Takami M., Pyo T.-S., 2008, ApJ, 676, 472Bicknell G. V., 1984, ApJ, 286, 68Binette L., Cabrit S., Raga A. C., Cant´o J., 1999, A&A, 266, 260Birch S. F., Eggers J. M., 1972, in Free Turbulent Shear Flows, Vol. I. pp11–40Blandford R. D., Payne D. G., 1982, MNRAS, 199, 883Blanksby S. J., Ellison G. B., 2003, Acc. Chem. Res., 36, 255Bonito R., Orlando S., Miceli M., Peres G., Micela G., Favata F., 2011, ApJ,737, 54Burrows C. J., et al., 1996, ApJ, 473, 437Cant´o J., Raga A. C., 1991, ApJ, 372, 646Caratti o Garatti A., Garcia Lopez R., Weigelt G., Tambovtseva L. V., GrininV. P., Wheelwright H., Ilee J. D., 2013, A&A, 554, A66Chernin L. M., Masson C. R., de Gouveia Dal Pino E. M., Benz W., 1994,ApJ, 426, 204Coffey D., Bacciotti F., Woitas J., Ray T. P., Eisl¨offel J., 2004, ApJ, 604,758Coffey D., Bacciotti F., Ray T. P., Eisl¨offel J., Woitas J., 2007, ApJ, 663,350Coffey D., Bacciotti F., Podio L., 2008, ApJ, 689, 1112Dash S. M., Wolf D. E., Seiner J. M., 1985, AIAAJ, 23, 505Dopita M. A., Sutherland R. S., 2003, Astrophysics of the Diffuse Universe,3rd edn. Astronomy & Astrophysics Library, Springer-Verlag, Berlin,HeidelbergDopita M. A., Sutherland R. S., Nicholls D. C., Kewley L. J., Vogt F. P. A.,2013, ApJS, 208, 10Downes T. P., Ray T. P., 1998, A&A, 331, 1130Eisl¨offel J., Mundt R., 1998, AJ, 115, 1554Elias J. H., 1978, ApJ, 224, 857Favre A., 1969, in Sedov L. I., ed., , Problems of Hydrodynamics and Con-tinuum Mechanics. Society for Industrial and Applied Mathematics,Philadelphia, p. 231Ferreira J., Dougados C., Cabrit S., 2006, A&A, 453, 785G¨udel M., Skinner S. L., Briggs K. R., Audard M., Arzner K., Telleschi A.,2005, ApJ, 626, L53G¨udel M., Skinner S. L., Audard M., Briggs K. R., Cabrit S., 2008, A&A,478, 797G¨unther H. M., Matt S. P., Li Z.-Y., 2009, A&A, 493, 579Hartigan P. M., Raymond J. C., Pierson R., 2004, ApJ, 614, L69Hirth G. A., Mundt R., Solf J., 1997, A&AS, 126, 437Jenkins E. B., 2009, ApJ, 700, 1299Jenkins E. B., 2013, in The Life Cycle of Dust in the Universe: Obser-vations, Theory, and Laboratory Experiments. Proceedings of Science,online, p. 5Kuncic Z., Bicknell G. V., 2004, ApJ, 616, 669Larson R. B., 2003, Rep. Prog. Phys., 66, 1651Lavalley-Fouquet C., Cabrit S., Dougados C., 2000, A&A, 356, L41Lavalley C., Cabrit S., Ferruit P., Dougados C., Bacon R., 1997, A&A, 327,671L´opez-C´amara D., Raga A. C., 2010, ApJ, 723, 449Lynch C., Mutel R. L., G¨udel M., Ray T. P., Skinner S. L., Schneider P. C.,MNRAS000
Agra-Amboage V., Dougados C., Cabrit S., Reunanen J., 2011, A&A, 532,A59Agra-Amboage V., Cabrit S., Dougados C., Kristensen L. E., Ibgui L., Re-unanen J., 2014, A&A, 564, A11Ainsworth R. E., Scaife A. M. M., Ray T. P., Taylor A. M., Green D. A.,Buckle J. V., 2014, ApJ, 792, L18Akima H., 1970, J. ACM, 17, 589Allen M. G., Groves B. A., Dopita M. A., Sutherland R. S., Kewley L. J.,2008, ApJS, 178, 20Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ARAA, 47, 481Bacciotti F., Eisl¨offel J., 1999, A&A, 342, 717Bacciotti F., Mundt R., Ray T. P., Eisl¨offel J., Solf J., Camezind M., 2000,ApJ, 537, L49Bacciotti F., Ray T. P., Mundt R., Eisl¨offel J., Solf J., 2002, ApJ, 576, 222Bachiller R., 1996, ARAA, 34, 111Beck T. L., McGregor P. J., Takami M., Pyo T.-S., 2008, ApJ, 676, 472Bicknell G. V., 1984, ApJ, 286, 68Binette L., Cabrit S., Raga A. C., Cant´o J., 1999, A&A, 266, 260Birch S. F., Eggers J. M., 1972, in Free Turbulent Shear Flows, Vol. I. pp11–40Blandford R. D., Payne D. G., 1982, MNRAS, 199, 883Blanksby S. J., Ellison G. B., 2003, Acc. Chem. Res., 36, 255Bonito R., Orlando S., Miceli M., Peres G., Micela G., Favata F., 2011, ApJ,737, 54Burrows C. J., et al., 1996, ApJ, 473, 437Cant´o J., Raga A. C., 1991, ApJ, 372, 646Caratti o Garatti A., Garcia Lopez R., Weigelt G., Tambovtseva L. V., GrininV. P., Wheelwright H., Ilee J. D., 2013, A&A, 554, A66Chernin L. M., Masson C. R., de Gouveia Dal Pino E. M., Benz W., 1994,ApJ, 426, 204Coffey D., Bacciotti F., Woitas J., Ray T. P., Eisl¨offel J., 2004, ApJ, 604,758Coffey D., Bacciotti F., Ray T. P., Eisl¨offel J., Woitas J., 2007, ApJ, 663,350Coffey D., Bacciotti F., Podio L., 2008, ApJ, 689, 1112Dash S. M., Wolf D. E., Seiner J. M., 1985, AIAAJ, 23, 505Dopita M. A., Sutherland R. S., 2003, Astrophysics of the Diffuse Universe,3rd edn. Astronomy & Astrophysics Library, Springer-Verlag, Berlin,HeidelbergDopita M. A., Sutherland R. S., Nicholls D. C., Kewley L. J., Vogt F. P. A.,2013, ApJS, 208, 10Downes T. P., Ray T. P., 1998, A&A, 331, 1130Eisl¨offel J., Mundt R., 1998, AJ, 115, 1554Elias J. H., 1978, ApJ, 224, 857Favre A., 1969, in Sedov L. I., ed., , Problems of Hydrodynamics and Con-tinuum Mechanics. Society for Industrial and Applied Mathematics,Philadelphia, p. 231Ferreira J., Dougados C., Cabrit S., 2006, A&A, 453, 785G¨udel M., Skinner S. L., Briggs K. R., Audard M., Arzner K., Telleschi A.,2005, ApJ, 626, L53G¨udel M., Skinner S. L., Audard M., Briggs K. R., Cabrit S., 2008, A&A,478, 797G¨unther H. M., Matt S. P., Li Z.-Y., 2009, A&A, 493, 579Hartigan P. M., Raymond J. C., Pierson R., 2004, ApJ, 614, L69Hirth G. A., Mundt R., Solf J., 1997, A&AS, 126, 437Jenkins E. B., 2009, ApJ, 700, 1299Jenkins E. B., 2013, in The Life Cycle of Dust in the Universe: Obser-vations, Theory, and Laboratory Experiments. Proceedings of Science,online, p. 5Kuncic Z., Bicknell G. V., 2004, ApJ, 616, 669Larson R. B., 2003, Rep. Prog. Phys., 66, 1651Lavalley-Fouquet C., Cabrit S., Dougados C., 2000, A&A, 356, L41Lavalley C., Cabrit S., Ferruit P., Dougados C., Bacon R., 1997, A&A, 327,671L´opez-C´amara D., Raga A. C., 2010, ApJ, 723, 449Lynch C., Mutel R. L., G¨udel M., Ray T. P., Skinner S. L., Schneider P. C.,MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows Gayley K. G., 2013, ApJ, 766, 53Maurri L., Bacciotti F., Podio L., Ray T. P., Mundt R., Locatelli U., CoffeyD., 2014, A&A, 565, A110Micono M., Massaglia S., Bodo G., Rossi P., Ferrari A., 1998, A&A, 333,989Micono M., Bodo G., Massaglia S., Rossi P., Ferrari A., Rosner R., 2000,A&A, 360, 795Miura A., Pritchett P. L., 1982, J. Geophys. Res., 87, 7431Nicholls D. C., Dopita M. A., Sutherland R. S., 2012, ApJ, 752, 148Nicholls D. C., Dopita M. A., Sutherland R. S., Kewley L. J., Palay E.,2013, ApJS, 207, 21Nisini B., Caratti o Garatti A., Giannini T., Lorenzetti D., 2002, A&A, 393,1035Papamoschou D., Roshko A., 1988, J. Fluid. Mech., pp 453–477Pesenti N., Dougados C., Cabrit S., O’Brien D., Garcia P. J. V., Ferreira J.,2003, A&A, 410, 155Podio L., Medves S., Bacciotti F., Eisl¨offel J., Ray T. P., 2009, A&A, 506,779Podio L., Eisl¨offel J., Melnikov S. Y., Hodapp K. W., Bacciotti F., 2011,A&A, 527, A13Pudritz R. E., Norman C. A., 1983, ApJ, 274, 677Pyo T.-S., et al., 2003, ApJ, 590, 340Raga A. C., Cant´o J., Binette L., Calvet N., 1990, ApJ, 364, 601Raga A. C., Cabrit S., Cant´o J., 1995, MNRAS, 273, 422Ray T. P., Ershkovich A. I., 1983, MNRAS, 204, 821Reipurth B., Bachiller R., 1997, in Latter W. B., ed., IAU Symp. Vol.170, CO: Twenty-Five Years of Millimeter-Wave Spectroscopy. KluwerAcademic Publishers, Dordrecht, pp 165–174Reipurth B., Bally J., 2001, ARAA, 39, 403Reipurth B., Heathcote S., Morse J. A., Hartigan P. M., Bally J., 2002, AJ,123, 362Rodr´ıguez-Gonz´alez A., Esquivel A., Raga A. C., Cant´o J., Riera A., CurielS., Beck T. L., 2012, AJ, 143, 60Romanova M. M., Ustyugova G. V., Koldoba A. V., Lovelace R. V. E., 2009,MNRAS, 399, 1802Schneider P. C., Eisl¨offel J., G¨udel M., G¨unther H. M., Herczeg G. J., Ro-brade J., Schmitt J. H. M. M., 2013, A&A, 550, L1Shakura N. I., Sunyaev R. A., 1973, A&A, 24, 337Shu F. H., Najita J. R., Ostriker E. C., Wilkin F., Ruden S. P., Lizano S.,1994, ApJ, 429, 781Slessor M. D., Zhuang M., Dimotakis P. E., 2000, J. Fluid. Mech., 414, 35Sutherland R. S., Dopita M. A., 1993, ApJS, 88, 253Takami M., Chrysostomou A., Ray T. P., Davis C. J., Dent W. R. F., BaileyJ., Tamura M., Terada H., 2004, A&A, 416, 213Trussoni E., 2008, in Massaglia S., Bodo G., Mignone A., Rossi P., eds,Lect. Notes. Phys., Vol. 754, Jets from Young Stars III. Springer-Verlag,Berlin, Heidelberg, pp 105–130, doi:10.1007/978-3-540-76967-5Vreman A. W., Sandham N. D., Luo K. H., 2006, Journal of Fluid Mechan-ics, 320, 235White M. C., McGregor P. J., Bicknell G. V., Salmeron R., Beck T. L.,2014a, MNRAS, 441, 1681White M. C., Bicknell G. V., McGregor P. J., Salmeron R., 2014b, MNRAS,442, 28Woitas J., Ray T. P., Bacciotti F., Davis C. J., Eisl¨offel J., 2002, ApJ, 580,336Zanni C., Ferrari A., Rosner R., Bodo G., Massaglia S., 2007, A&A, 469,811MNRAS , 000–000 (0000) M. C. White et al.
APPENDIX A: SUPPLEMENTARY CALCULATIONSA1 Characteristic Equations of MHD
For magnetohydrodynamic fluids with density ρ , velocity v , pressure p , viscous stress tensor t vij , magnetic field B immersed in a gravitationalpotential φ G , the equation of mass continuity can be written in Cartesian coordinate notation thus: ∂ρ∂t + ∂ ( ρv i ) ∂x i = 0 . (A1)Similarly, the equation of momentum conservation is written as ∂ ( ρv i ) ∂t + ∂ ( ρv i v j ) ∂x j = − ρ ∂φ G ∂x i − ∂p∂x i + ∂t Bij ∂x j + ∂t vij ∂x j (A2)(Kuncic & Bicknell 2004).Time-averaging equation (A1) is trivial, yielding equation (6). Averaging the momentum conservation equation, equation (A2), is morecomplex. The viscous stress tensor, t vij , is disregarded, as it is unimportant to the transfer of momentum on large scales. The Reynolds stresstensor, t Rij , appears in equation (7) as a result of averaging the second term on the left-hand side of equation (A2). This may then be combinedwith the magnetic stress tensor, t Bij , to form a single term encapsulating the total stress in the system, t ij = t Rij + t Bij . A2 Order of Magnitude Calculations
In this section, we conduct an order of magnitude analysis of the characteristic equations in § t xx term from equation (11);(ii) Demonstrate the constancy of total (thermal plus turbulent) pressure across the mixing layer.Consider the orders of magnitude of the terms in the equation of mass continuity, equation (6), as shown by the expressions beneath thebraces: ∂ (¯ ρ ˜ v x ) ∂x (cid:124) (cid:123)(cid:122) (cid:125) ¯ ρ ˜ v x /L + ∂ (¯ ρ ˜ v y ) ∂y (cid:124) (cid:123)(cid:122) (cid:125) ¯ ρ ˜ v y /h = 0 ⇒ ˜ v y ∼ hL ˜ v x . (A3)It can be shown from our definition of ξ , equation (15), and the transformed equation of mass continuity, equation (23), that h/L ∼ h (cid:48) .It then follows that the orders of magnitude in the equation of conservation of streamwise momentum, equation (11), are: ∂ (¯ ρ ˜ v x ) ∂x (cid:124) (cid:123)(cid:122) (cid:125) ρv /L + ∂ (¯ ρ ˜ v x ˜ v y ) ∂y (cid:124) (cid:123)(cid:122) (cid:125) ρv /L = − ∂ (cid:104) ρv (cid:48) x (cid:105) ∂x (cid:124) (cid:123)(cid:122) (cid:125) (R1): ρv t /L − ∂ (cid:104) ρv (cid:48) x v (cid:48) y (cid:105) ∂y (cid:124) (cid:123)(cid:122) (cid:125) (R2): ρv t /h + ∂∂x (cid:18) (cid:104) B (cid:48) x (cid:105) − (cid:104) B (cid:48) y (cid:105) − (cid:104) B (cid:48) z (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R3): B t / πL + ∂∂y (cid:18) (cid:104) B (cid:48) x B (cid:48) y (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R4): B t / πh . (A4)Note that we have expanded the stress terms, and for the purposes of expressing orders of magnitude, have replaced ¯ ρ with ρ , ˜ v x with v , andhave denoted the generic turbulent velocity, typically v (cid:48) , as v t for clarity. Because h (cid:48) (cid:28) ⇒ L (cid:29) h , terms (R2) and (R4) of equation (A4)are the most important terms on the right-hand side. Therefore, the t xx component of the turbulent stress (encapsulated in the third term onthe right-hand side of equation A4) is unimportant ( § II ] 1.644 µ m emission-line component from DG Tau is (cid:38) km s − (Paper I), indicating that the turbulent velocity in the mixing layer, v t ,is of that order. Recalling that we assume that the magnetic field is dominated by its turbulent component within the mixing layer ( § v A,t = B/ √ πρ , using the inferred magnetic field strength in the DG Tau jet, which is in therange 30-100 µ G (Lavalley-Fouquet et al. 2000; Ainsworth et al. 2014). Assuming an equipartition magnetic field of ∼ µ G, this yields aturbulent Alfv´en velocity of ∼ − (cid:28) v t ; hence, the turbulent velocity term, (R2), is dominant. This implies that ρv L ∼ ρv t h ⇒ v t ∼ hL v . (A5)Consider the equation for the conservation of transverse momentum, i.e., the i = y expansion of equation (7): ∂ (¯ ρ ˜ v x ˜ v y ) ∂x (cid:124) (cid:123)(cid:122) (cid:125) L1: ( ρv /L )( h/L ) + ∂ (¯ ρ ˜ v y ) ∂y (cid:124) (cid:123)(cid:122) (cid:125) L2: ( ρv /L )( h/L ) = − ∂ ¯ p∂y (cid:124) (cid:123)(cid:122) (cid:125) (R1): p/h − ∂ (cid:104) ρv (cid:48) x v (cid:48) y (cid:105) ∂x (cid:124) (cid:123)(cid:122) (cid:125) (R2): ( ρv /L )( h/L ) − ∂ (cid:104) ρv (cid:48) y (cid:105) ∂y (cid:124) (cid:123)(cid:122) (cid:125) (R3): ρv /L + ∂∂x (cid:18) (cid:104) B (cid:48) x B (cid:48) y (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R4): B t / πL + ∂∂y (cid:18) (cid:104) B (cid:48) y (cid:105) − (cid:104) B (cid:48) x (cid:105) − (cid:104) B (cid:48) z (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R5): B t / πh . (A6)On the right-hand side of this equation, terms (R3) and (R5) dominate terms (R2) and (R4); then, as in equation (A4), the turbulent velocityterm (R3) dominates the turbulent magnetic field term (R5). It can also be seen that term (R3) also dominates both of the terms on theleft-hand side of the equation because h/L (cid:28) . Therefore, the two dominant terms of equation (A6) are terms (R1) and (R3).Discarding the unimportant terms of equation (A6), and integrating over y , yields the ‘total’ pressure, p tot ≡ ¯ p + (cid:104) ρv (cid:48) y (cid:105) = const. (A7) MNRAS000
In this section, we conduct an order of magnitude analysis of the characteristic equations in § t xx term from equation (11);(ii) Demonstrate the constancy of total (thermal plus turbulent) pressure across the mixing layer.Consider the orders of magnitude of the terms in the equation of mass continuity, equation (6), as shown by the expressions beneath thebraces: ∂ (¯ ρ ˜ v x ) ∂x (cid:124) (cid:123)(cid:122) (cid:125) ¯ ρ ˜ v x /L + ∂ (¯ ρ ˜ v y ) ∂y (cid:124) (cid:123)(cid:122) (cid:125) ¯ ρ ˜ v y /h = 0 ⇒ ˜ v y ∼ hL ˜ v x . (A3)It can be shown from our definition of ξ , equation (15), and the transformed equation of mass continuity, equation (23), that h/L ∼ h (cid:48) .It then follows that the orders of magnitude in the equation of conservation of streamwise momentum, equation (11), are: ∂ (¯ ρ ˜ v x ) ∂x (cid:124) (cid:123)(cid:122) (cid:125) ρv /L + ∂ (¯ ρ ˜ v x ˜ v y ) ∂y (cid:124) (cid:123)(cid:122) (cid:125) ρv /L = − ∂ (cid:104) ρv (cid:48) x (cid:105) ∂x (cid:124) (cid:123)(cid:122) (cid:125) (R1): ρv t /L − ∂ (cid:104) ρv (cid:48) x v (cid:48) y (cid:105) ∂y (cid:124) (cid:123)(cid:122) (cid:125) (R2): ρv t /h + ∂∂x (cid:18) (cid:104) B (cid:48) x (cid:105) − (cid:104) B (cid:48) y (cid:105) − (cid:104) B (cid:48) z (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R3): B t / πL + ∂∂y (cid:18) (cid:104) B (cid:48) x B (cid:48) y (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R4): B t / πh . (A4)Note that we have expanded the stress terms, and for the purposes of expressing orders of magnitude, have replaced ¯ ρ with ρ , ˜ v x with v , andhave denoted the generic turbulent velocity, typically v (cid:48) , as v t for clarity. Because h (cid:48) (cid:28) ⇒ L (cid:29) h , terms (R2) and (R4) of equation (A4)are the most important terms on the right-hand side. Therefore, the t xx component of the turbulent stress (encapsulated in the third term onthe right-hand side of equation A4) is unimportant ( § II ] 1.644 µ m emission-line component from DG Tau is (cid:38) km s − (Paper I), indicating that the turbulent velocity in the mixing layer, v t ,is of that order. Recalling that we assume that the magnetic field is dominated by its turbulent component within the mixing layer ( § v A,t = B/ √ πρ , using the inferred magnetic field strength in the DG Tau jet, which is in therange 30-100 µ G (Lavalley-Fouquet et al. 2000; Ainsworth et al. 2014). Assuming an equipartition magnetic field of ∼ µ G, this yields aturbulent Alfv´en velocity of ∼ − (cid:28) v t ; hence, the turbulent velocity term, (R2), is dominant. This implies that ρv L ∼ ρv t h ⇒ v t ∼ hL v . (A5)Consider the equation for the conservation of transverse momentum, i.e., the i = y expansion of equation (7): ∂ (¯ ρ ˜ v x ˜ v y ) ∂x (cid:124) (cid:123)(cid:122) (cid:125) L1: ( ρv /L )( h/L ) + ∂ (¯ ρ ˜ v y ) ∂y (cid:124) (cid:123)(cid:122) (cid:125) L2: ( ρv /L )( h/L ) = − ∂ ¯ p∂y (cid:124) (cid:123)(cid:122) (cid:125) (R1): p/h − ∂ (cid:104) ρv (cid:48) x v (cid:48) y (cid:105) ∂x (cid:124) (cid:123)(cid:122) (cid:125) (R2): ( ρv /L )( h/L ) − ∂ (cid:104) ρv (cid:48) y (cid:105) ∂y (cid:124) (cid:123)(cid:122) (cid:125) (R3): ρv /L + ∂∂x (cid:18) (cid:104) B (cid:48) x B (cid:48) y (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R4): B t / πL + ∂∂y (cid:18) (cid:104) B (cid:48) y (cid:105) − (cid:104) B (cid:48) x (cid:105) − (cid:104) B (cid:48) z (cid:105) π (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) (R5): B t / πh . (A6)On the right-hand side of this equation, terms (R3) and (R5) dominate terms (R2) and (R4); then, as in equation (A4), the turbulent velocityterm (R3) dominates the turbulent magnetic field term (R5). It can also be seen that term (R3) also dominates both of the terms on theleft-hand side of the equation because h/L (cid:28) . Therefore, the two dominant terms of equation (A6) are terms (R1) and (R3).Discarding the unimportant terms of equation (A6), and integrating over y , yields the ‘total’ pressure, p tot ≡ ¯ p + (cid:104) ρv (cid:48) y (cid:105) = const. (A7) MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows This relationship allows us to define a new quantity, W : W ≡ ˜ T thermal + µm (cid:104) ρv (cid:48) y (cid:105) k ¯ ρ , (A8)which is related to density and total pressure thus: p tot = ¯ ρkµm W . (A9)We approximate W as being linear across the mixing layer, with bounds given by equation (19). This is a reasonable approximation, given:(i) The constancy of total pressure across the mixing layer (equation A7);(ii) The fact that turbulent velocities (that is, velocity dispersion), are present in at least both the mixing layer and the jet (Paper I, fig. 6therein). A3 Dimensionless Functions
The dimensionless function D ( η, ξ ) in the mixing layer transverse velocity profile, equation (24), is given by D ( η, ξ ) = (cid:90) ξξ ξ (cid:48) dd ξ (cid:48) (cid:18) U ( ξ (cid:48) ) η + (1 − η ) S ( ξ (cid:48) ) (cid:19) d ξ (cid:48) (A10) = η ( η − (cid:26) ( η − (cid:20) ξ ( η − ξ − ξ ) + 1 − ξ (cid:21) + log [( η − ξ − ξ ) + 1] (cid:27) . (A11)The dimensionless function F ( η, ξ ) in the mixing layer transverse turbulent stress profile, equation (26), is given by F ( η, ξ ) = (cid:90) ξξ − ξ (cid:48) ξ (cid:48) − ξ + 1 η + (1 − η )( ξ (cid:48) − ξ + 1) + D ( η, ξ (cid:48) ) d ξ (cid:48) . (A12) = 12( η − (cid:110) η [(1 − η )( ξ − ξ ) −
1] log [( η − ξ − ξ ) + 1] − [( η − η − ξ − η ( ξ −
2) + ξ ] ( η −
1) ( ξ − ξ ) (cid:111) . (A13) A4 Mixing Layer Transverse Velocity and Turbulent Stress Profiles
The transverse velocity profile across the mixing layer, ˜ v y ( η, ξ ) , may be found by substituting the expression for ξ ( η ) , equation (27), intoequation (24): ˜ v y ( η, ξ ) = v j h (cid:48) ( x ) η η − (cid:40) (cid:2) η (2 log( η ) −
1) + 1 (cid:3) (cid:2) η − η log( η ) + 2( η − ξ − η + 1 (cid:3) − η − (cid:2) η − η log( η ) + 2( η − ξ − (cid:3) log (cid:18) − η (1 − η )) + 2( η − ξ − η − (cid:19) (cid:41) . (A14)The turbulent stress profile across the mixing layer, t xy ( η, ξ ) , may be found in the same way, using equation (26): t xy ( η, ξ ) = ρ j v j h (cid:48) ( x )8( η − (cid:40) η (cid:104) ( η − (cid:2) η (1 − η )) + 2( η − ξ − (cid:3) × log (cid:18) − η (1 − η )) + 2( η − ξ − η − (cid:19) + η log( η ) (cid:16) η (cid:0) − η − η − ξ + 7 η − (cid:1) + (2 η − η log( η ) + 1 (cid:17)(cid:105) + ( η − (cid:104) η − ξ + 4 η ( η + 1)( η − ξ + η ( η (6 η −
5) + 4) − (cid:105)(cid:41) . (A15) MNRAS , 000–000 (0000) M. C. White et al. − − − Jet-to-wind density ratio, η . . E n e r gy r e qu i r e d f o r H d i ss o c i a ti on / t o t a lt u r bu l e n t e n e r gyp r odu ce d Figure B1.
The ratio of the energy required to dissociate all the H molecules entrained into the turbulent mixing layer to the total turbulent energy producedwithin the mixing layer. The ratio is plotted as a function of the jet-to-ambient-wind density ratio, η . The ratio is computed assuming a jet velocity of km s − , a mixing layer growth rate of 0.05, and a jet density of . × − g cm − , as for DG Tau ( § A5 Calculation of the Mass Entrainment Rate
The mass entrainment rate is given by taking the x -derivative of equation (32), ˙ M (cid:48) ent = dd x (cid:90) ξ ξ ρ ( ξ ) v x ( x, ξ ) h ( x ) d ξ − ρ j v j h (cid:48) ( x ) ξ ( η ) (A16) = ¯ ρ ( ξ )˜ v x ( ξ ) h ( x ) d ξ d x − ¯ ρ ( ξ )˜ v x ( ξ ) h ( x ) d ξ d x + (cid:90) ξ ξ dd x [ ρ ( ξ ) v x ( x, ξ ) h ( x )] d ξ − ρ j v j h (cid:48) ( x ) ξ ( η ) . (A17)The first two terms of the above are zero in ( x, ξ ) -space. The transformed equation of continuity, equation (23), may be written as ∂∂x (¯ ρ ˜ v x h ( x )) = h (cid:48) ( x ) ∂∂ξ ( ξ ¯ ρ ˜ v x ) − ∂∂ξ (¯ ρ ˜ v y ) , (A18)which reduces equation (A17) to ˙ M (cid:48) ent = h (cid:48) ( x ) [ ρ j v j ξ − ρ w v w ξ ] − ρ J ˜ v y ( x, ξ ) + ρ w ˜ v y ( x, ξ ) − ρ j v j ξ h (cid:48) ( x ) . (A19)Most of these terms are zero, or cancel, leaving ˙ M (cid:48) ent = ρ w ˜ v y ( x, ξ ) (A20) = ρ w v ent by definition. (A21) APPENDIX B: DISSOCIATION OF ENTRAINED MOLECULAR HYDROGEN WITHIN THE MIXING LAYER
The calculation presented in this paper does not account for the energy required to dissociate hydrogen molecules (H ) that are entrainedinto the mixing layer. In this Appendix, we show that the energy required to dissociate the entrained H is of the order of 1 per cent of allturbulent energy produced in the mixing layer, and hence has a negligible effect on the comparison to the observed mixing layer luminosityin § ˙ M ent , may be computed in a similar fashion to the total rate-of-change of mixing layer mass flux ( § ˙ M (cid:48) ent , equation 32) by the mixinglayer length L , and by πR mix , where R mix is an approximation of the mixing layer radius. Therefore, the total mass of molecular windmaterial entering the mixing layer at any one time is M ent = ρ jet η (cid:124)(cid:123)(cid:122)(cid:125) ρ w v jet h (cid:48) ( x ) η ( η − η log( η ) − η − (cid:124) (cid:123)(cid:122) (cid:125) v ent πR mix L . (B1)The bond dissociation energy of molecular hydrogen is 4.52 eV per molecule = 2 . × erg g − (Blanksby & Ellison 2003). Thetotal energy required to dissociate all of the molecular hydrogen present in the mixing layer may then be calculated by multiplying this valueby equation (B1). The ratio of this required energy is compared to the total turbulent energy produced within the mixing layer, equation MNRAS000
The calculation presented in this paper does not account for the energy required to dissociate hydrogen molecules (H ) that are entrainedinto the mixing layer. In this Appendix, we show that the energy required to dissociate the entrained H is of the order of 1 per cent of allturbulent energy produced in the mixing layer, and hence has a negligible effect on the comparison to the observed mixing layer luminosityin § ˙ M ent , may be computed in a similar fashion to the total rate-of-change of mixing layer mass flux ( § ˙ M (cid:48) ent , equation 32) by the mixinglayer length L , and by πR mix , where R mix is an approximation of the mixing layer radius. Therefore, the total mass of molecular windmaterial entering the mixing layer at any one time is M ent = ρ jet η (cid:124)(cid:123)(cid:122)(cid:125) ρ w v jet h (cid:48) ( x ) η ( η − η log( η ) − η − (cid:124) (cid:123)(cid:122) (cid:125) v ent πR mix L . (B1)The bond dissociation energy of molecular hydrogen is 4.52 eV per molecule = 2 . × erg g − (Blanksby & Ellison 2003). Thetotal energy required to dissociate all of the molecular hydrogen present in the mixing layer may then be calculated by multiplying this valueby equation (B1). The ratio of this required energy is compared to the total turbulent energy produced within the mixing layer, equation MNRAS000 , 000–000 (0000) urbulent mixing layers in protostellar outflows Table C1.
Pre-shock ionization fractions based on collisional ionization equilibrium models. χ = 0 . a χ = 0 . χ = 0 . T = . × K b T = . × K T = . × KSpecies I II III I II III I II IIIH 0.991 0.009 — 0.951 0.049 — 0.899 0.101 —He 1.000 0.000 0.000 1.000 0.000 0.000 0.997 0.003 0.000Fe 0.368 0.631 0.001 0.291 0.706 0.003 0.247 0.746 0.007 χ = 0 . χ = 0 . χ = 0 . T = . × K T = . × K T = . × KI II III I II III I II IIIH 0.699 0.301 — 0.502 0.498 — 0.301 0.699 —He 1.000 0.000 0.000 1.000 0.000 0.000 1.000 0.000 0.000Fe 0.170 0.797 0.033 0.125 0.789 0.085 0.083 0.712 0.205 a Shown are the fractions of each species in the neutral (I), singly ionized (II) and doubly ionized (III)states as a function of hydrogen ionization fraction, χ . b Temperatures determined from the collisional ionization equilibrium model. (41), in Fig. B1. The energy required to dissociate the entrained molecular material is of the order of 1 per cent of the total turbulent energyproduced. Given that our shock models indicate that at most 1 per cent of the total turbulent energy produced is radiated away as [Fe II ] 1.644 µ m emission, the effect on our predicted [Fe II ] mixing layer luminosity is of the order of 0.01 per cent, i.e. negligible. APPENDIX C: SHOCK MODELLING SUPPLEMENTARY INFORMATIONC1 Pre-Shock Gas Parameters
A collisional ionization equilibrium (CIE) model was used to determine the relative abundances of various ions in the pre-shock gas. Theresults of these models are shown in Table C1. Conversion of these fractions to number densities was then achieved by multiplying bythe relevant elemental abundances (Asplund et al. 2009 for solar abundances, and Jenkins 2009, 2013 for the depleted case). The weakdependence of [Fe II ] 1.644 µ m luminosity fraction on pre-shock ionization fraction (Figs. 6, 7) can be seen in this Table; the relative fractionof iron in the singly ionized state does not vary significantly with hydrogen ionization. C2 Effect of Varying Pre-Shock Density and Iron Depletion
Shown in Figs. C1 and C2 are the dependencies of the fraction of total emission generated by [Fe II ] 1.644 µ m line emission as a functionof pre-shock ionization and pre-shock number density (cf. Figs. 6 and 7, where we show the [Fe II ] 1.644 µ m line emission fraction asa function of pre-shock ionization and plasma β ), for solar abundances (Fig. C1) and for iron depleted by a factor of 100 below solarabundance (Fig. C2). As for variations with pre-shock plasma β , for the purposes of our study, there is relatively little effect on the fractionof total luminosity emitted in [Fe II ] 1.644 µ m; the fraction is ∼ − for solar abundances, and ∼ − for an iron depletion factor of 100.The total emission from these partially-ionized slow shocks arises from a non-linear, time-dependent ionization state. This situationcomes about via competition between collisional ionization caused by the post-shock temperature jump, and recombination that follows asthe gas cools quickly. The ionization balance is generally far from equilibrium and is somewhat dependent on the initial ionization fractions,as the shocks never reach the hot equilibrium post-shock phase seen in faster shocks (e.g. Allen et al. 2008). This makes the initial ionizationstate important in determining the final shock emission spectra. The cooling arises primarily from the collisional excitation of lines (e.g. H α ,Ly α ) and two-photon continuum emission of hydrogen, as well as a large number of forbidden and fine-structure transitions in metal speciessuch as neutral oxygen or Fe II . As the pre-ionization state changes, the availability of electrons and the changing mean molecular weight(which alters the post-shock temperature), combined with the changing initial values for ion abundances and magnetic pressure, means thatthe integrated contribution of a particular line in a particular species (such as [Fe II ] 1.644 µ m) can vary both up and down as the mix ofcompeting processes and species change.For example, the density reached when an important ion is most abundant may vary between models, depending on the temperatureprofile and the contribution of magnetic pressure support. Density ‘quenching’ via collisional de-excitation can change the emission fromsome lines but not others, depending on their critical densities (cf. Dopita & Sutherland 2003). Therefore, in any given model, the emissioncontribution from species competing with [Fe II ] 1.644 µ m emission will change, and so the Fe II cooling fraction will also change. Fortu-nately, the overall integrated variations in the [Fe II ] fraction range over only factors of a few for the range of pre-shock conditions consideredhere. The relatively small variations in the initial CIE Fe II fractions in Table C1 may contribute to the rough stability of the [Fe II ] 1.644 µ memission efficiency.Additionally, the collisional excitation of neutral hydrogen in particular, enhancing H α , Ly α and the two-photon continuum are inmany models a dominant coolant, and as the neutral hydrogen species is common throughout all the models, the overall cooling is not asvariable as it would be if the cooling were dominated by fleeting species. Likewise, the Fe II ion is often the most common Fe ion in themodels. However, in any given model, the detailed efficiency outcome is difficult to predict in this highly non-linear system, so the variationsseen in Figs. 6, 7, C1, and C2 are best interpreted in general terms, as specifics are difficult to isolate and prove. MNRAS , 000–000 (0000) M. C. White et al. − − − Ionization fraction, χ . . . . [ F e II ] . m l u m i no s it y fr ac ti on , L . / L t o t × − Pre-shock number density, n H Post-shock temperature ( × K)2.0 3.0 4.0 5.0 6.0Post-shock temperature ( × K)2.0 3.0 4.0 5.0 6.0
Figure C1.
Contribution of [Fe II ] 1.644 µ m line emission to total shock emission from a grid of MAPPINGS
IV models, where the iron abundance is notdepleted from solar abundance values and the pre-shock hydrogen ionization χ = 0 . . The percentage of shock emission which is radiated as [Fe II ] 1.644 µ mline emission is shown as a function of the pre-shock ionization fraction, χ . Models have been computed for a range of post-shock temperatures and pre-shocknumber densities as indicated. Points are joined using a one-dimensional Akima spline interpolation (Akima 1970). The vertical dotted line denotes the knownjet ionization, χ = 0 . . − − − Ionization fraction, χ [ F e II ] . m l u m i no s it y fr ac ti on , L . / L t o t × − Pre-shock number density, n H Post-shock temperature ( × K)2.0 3.0 4.0 5.0 6.0Post-shock temperature ( × K)2.0 3.0 4.0 5.0 6.0
Figure C2.
As for Fig. C1, but with an iron depletion factor of 100. MNRAS000