The ALMA-PILS survey: Gas dynamics in IRAS 16293$-$2422 and the connection between its two protostars

TL;DR

This study uses ALMA data to analyze the gas dynamics and interactions between the two protostars in IRAS 16293-2422, revealing outflows, a quiescent bridge, and possible evolutionary differences.

Contribution

First detailed kinematic analysis of the IRAS 16293-2422 system at 60-1000 au scales using molecular line tracers from ALMA-PILS survey.

Findings

Detection of multiple outflows from protostar A.

Identification of a quiescent gas and dust bridge between the protostars.

Evidence suggesting protostar B is at an earlier evolutionary stage.

Abstract

[Abridged] The majority of stars form in binary or higher order systems. The Class 0 protostellar system IRAS16293-2422 contains two protostars, 'A' and 'B', separated by ~600 au and embedded in a single, 10^4 au scale envelope. Their relative evolutionary stages have been debated. We aim to study the relation and interplay between the two protostars A and B at spatial scales of 60 to ~1000 au. We selected molecular gas line transitions of CO, H2CO, HCN, CS, SiO, and CCH from the ALMA-PILS spectral imaging survey (329-363 GHz) and used them as tracers of kinematics, density, and temperature in the IRAS16293-2422 system. The angular resolution of the PILS data set allows us to study these quantities at a resolution of 0.5 arcsec (60 au [..]). Line-of-sight velocity maps of both optically thick and optically thin molecular lines reveal: (i) new manifestations of previously known outflows emanating from protostar A; (ii) a kinematically quiescent bridge of dust and gas spanning between the two protostars, with an inferred density between 4 10^4 and 3 10^7 cm^-3; and (iii) a separate, straight filament seemingly connected to protostar B seen only in CCH, with a flat kinematic signature. Signs of various outflows, all emanating from source A, are evidence of high-density and warmer gas; none of them coincide spatially and kinematically with the bridge. We hypothesize that the bridge arc is a remnant of filamentary substructure in the protostellar envelope material from which protostellar sources A and B have formed. One particular morphological structure appears to be due to outflowing gas impacting the quiescent bridge material. The continuing lack of clear outflow signatures unambiguously associated to protostar B and the vertically extended shape derived for its disk-like structure lead us to conclude that source B may be in an earlier evolutionary stage than source A.

Tables

Table 1. Table 1: Coordinates of submillimeter continuum peaks in IRAS 16293 − - 2422.

Component	Right Ascension	Declination
	(J2000)	(J2000)
source A	16h32m22$\stackrel{s}{.}$873	$-$24°28′36$\stackrel{\prime \prime }{.}$54
source B	16h32m22$\stackrel{s}{.}$6147	$-$24°28′32$\stackrel{\prime \prime }{.}$566

Table 2. Table 2: Selected line transitions.

Molecule	Transition	Rest frequencya𝑎aa𝑎a$a$Rest frequencies taken from the Cologne Database for Molecular Spectroscopy, CDMS (Müller et al. 2005; Endres et al. 2016). Original sources of spectroscopic data: C17O, Klapper et al. (2003); CO, Winnewisser et al. (1997); H13CN, Cazzoli & Puzzarini (2005); H2CO and H213CO, Cornet & Winnewisser (1980); C34S Gottlieb et al. (2003); SiO and 29SiO, Müller et al. (2013); C2H, Padovani et al. (2009).	$E_{\mathrm{up}}$/$k$	$n_{\mathrm{crit}}$b𝑏bb𝑏b$b$Critical density is calculated as the ratio of the Einstein $A$ coefficient and the collision rate at a temperature of 100 K, which are accessed through the LAMDA database (Schöier et al. 2005), with the following original sources: CO (adopting o-H2 as dominant collision partner), Yang et al. (2010); H13CN, Dumouchel et al. (2010) for main isotopologue; H2CO and adopting identical rates for H213CO, Wiesenfeld & Faure (2013); C34S adopting rates for the main isotopologue C32S, Lique et al. (2006); SiO and adopting identical rates for 29SiO, Dayou & Balança (2006); C2H, Spielfiedel et al. (2012).	Traced componentf𝑓ff𝑓f$f$source A and source B are abbreviated as ‘A’ and ‘B’; ‘PA’ denotes position angle east of north, following Fig. 1.
		(GHz)	(K)	(${\mathrm{cm}}^{-3}$)
Dust	(continuum)	344 (0.87 mm)	-	-	Quiescent bridge between A and B.
C17O	$J$ = 3–2	337.061c𝑐cc𝑐c$c$The hyperfine splitting of C17O 3–2 and H13CN 4–3 is discussed in Appendix A.	32.5	4.1$\times$104	Quiescent bridge between A and B.
12C16O d𝑑dd𝑑d$d$For 12C16O, we study morphological features using only velocity channels at $|$ $V_{\mathrm{lsr}}$$-$ $v_{\mathrm{systemic}}$ $|>5$ $\mathrm{km}{s}^{-1}$, as shown in Fig. 4.	$J$ = 3–2	345.796	33.2	4.1$\times$104	Redshifted, outflow cone accelerating westward from A;
					onset of eastern counterpart of westward outflow;
					blueshifted knot south of B.
H13CN	$J$ = 4–3	345.33977	41.4	1.2$\times$108	Redshifted northern wall of PA 270° outflow;
					Redshifted outflow southeast (PA 135°) from A.
o-H2CO	$J_{K_a,K_c}$= 51,5 – 41,4	351.76864	62.5	7.5$\times$107	Redshifted northern wall of PA 270° outflow;
					Blueshifted interface wall with bridge.
					Redshifted outflow southeast (PA 135°) from A.
o-H213CO	$J_{K_a,K_c}$= 51,5 – 41,4	343.32571	61.3	7.5$\times$107	Redshifted northern wall of PA 270° outflow;
C34S	$J$ = 7–6	337.39646	50.2	1.5$\times$108	Blueshifted gas in same direction as ‘axis’ of bridge.
SiO	$J$ = 8–7	347.33058	75.0	6.1$\times$107	Axes of of redshifted (PA 135°) and
					blueshifted (315°) outflows from A;
					Mixed-velocity pockets north and south of B;
29SiO	$J$ = 8–7	342.98084	74.1	5.9$\times$107	Redshifted outflow southeast (PA 135°) from A.
					Mixed-velocity pockets north and south of B;
C2H	$N_J$ = $4_{7/2}$ – $3_{5/2}$	349.39997e𝑒ee𝑒e$e$The rest frequency for C2H $4_{7/2}$ – $3_{5/2}$ is the average of the tabulated rest frequencies for the individual hyperfine components, $F$=4–3 at 349.39928 GHz and $F$=3–2 at 349.40067 GHz.	41.9	1.4$\times$107	Quiescent, narrow filament across B.

Table 3. Table 3: LIME model parameters.

Molecule	$X$/H2	$T_{\mathrm{subl}}$a𝑎aa𝑎a$a$Ice mantle sublimation temperature for HCN follows those of other CN-bearing species (e.g., Noble et al. 2013). That of H2CO follows Aikawa et al. (1997), Ceccarelli et al. (2001) and Rodgers & Charnley (2003).	Depletion	$b$c𝑐cc𝑐c$c$Doppler parameter as defined in LIME, that is, the 1/e half-width of a thermally broadened line.
	range	[K]	factorb𝑏bb𝑏b$b$When depletion is applied in the model, the gas phase abundance of the molecule in question is divided by the depletion factor.	[$\mathrm{km}{s}^{-1}$]
C17O	5.58$\times$10-8	30	100	1
o-H2CO	[$1\times {10}^{-8}$, $3\times {10}^{-6}$]	50	100	1
H13CN	[$2\times {10}^{-12}$, $2\times {10}^{-10}$]	100	100	1

Equations

$$M/L=\frac{2k{T}_{\mathrm{kin}}}{\mu m_{0}G}.$$

Full text

*Aims. * We aim to study the relation and interplay between the two protostars A and B at spatial scales of 60 au up to \sim$$10^{3} au.

*Methods. * We selected molecular gas line transitions of the species CO, H2CO, HCN, CS, SiO, and C2H from the ALMA-PILS spectral imaging survey (329–363 GHz) and used them as tracers of kinematics, density, and temperature in the IRAS 16293$-$2422 system. The angular resolution of the PILS data set allows us to study these quantities at a resolution of 0.5″ (60 au at the distance of the source).

Results. * Line-of-sight velocity maps of both optically thick and optically thin molecular lines reveal: (i) new manifestations of previously known outflows emanating from protostar A; (ii) a kinematically quiescent bridge of dust and gas spanning between the two protostars, with an inferred density between 4$\times$104* $\mathrm{cm}^{-3}$ and $\sim$3$\times$107 $\mathrm{cm}^{-3}$; and (iii) a separate, straight filament seemingly connected to protostar B seen only in C2H, with a flat kinematic signature. Signs of various outflows, all emanating from source A, are evidence of high-density and warmer gas; none of them coincide spatially and kinematically with the bridge.

*Conclusions. * We hypothesize that the bridge arc is a remnant of filamentary substructure in the protostellar envelope material from which protostellar sources A and B have formed. One particular morphological structure appears to be due to outflowing gas impacting the quiescent bridge material. The continuing lack of clear outflow signatures unambiguously associated to protostar B and the vertically extended shape derived for its disk-like structure lead us to conclude that source B may be in an earlier evolutionary stage than source A.

1 Introduction

The majority of currently forming stars are part of multiple systems (Chen et al. 2013; Duchêne & Kraus 2013; Tobin et al. 2016b). Two main scenarios have been proposed to form stellar systems of binary and higher order: disk fragmentation (Adams et al. 1989; Bonnell & Bate 1994), leading to close binaries (up to a few hundred au separation), and turbulent fragmentation of a natal protostellar envelope (Offner et al. 2010; Pineda et al. 2015), leading to wide separation binaries ($\gtrsim$1000 au). The difference in initial separation is often obfuscated at later stages by the effects of tidal interactions after the initial formation of protostars. Another observable difference that can help to distinguish between formation scenarios is that binaries formed through disk instability may be more prone to exhibit aligned rotation axes, whereas rotation axes of those formed from turbulent fragmentation should be randomly distributed. The rotation axis of a protostar can be inferred through studying the morphology and kinematics of its disk-outflow system (Lee et al. 2016). High angular resolution observations of embedded protostellar systems suggest that both mechanisms of multiple star formation occur in nature (e.g., Tobin et al. 2016a; Brinch et al. 2016; Lee et al. 2017). In both cases, remnant material of the natal body may manifest as a bridge between recently formed binary companions. The object of study in this paper is one of the most well-studied protostellar binary systems, IRAS 16293$-$2422 (hereafter IRAS 16293), which has such a bridge and exhibits kinematic signatures of active outflow and infall (see below).

IRAS 16293 is a nearby, young, Class 0 protostellar system in the Ophiuchus cloud complex (see Sect. 2 of Jørgensen et al. 2016 for a review of the source). The distance to three other young stellar objects in the same cloud complex has been determined at 147$\pm$3 pc (Ortiz-León et al. 2017), and trigonometric parallax measurements of water maser spots in IRAS 16293 itself put it at a distance of 141${}^{+30}_{-21}$ pc (Dzib et al. 2018). However, to be consistent with earlier work and with detailed numerical modeling of the IRAS 16293 system conducted by our team, we adopt a distance estimate of 120 pc based on extinction measurements and VLBI parallax measurements of two stars in the core of the complex (Loinard et al. 2008). This value is within the 1$\sigma$ uncertainty margin of Dzib et al. (2018). In general, an 18% increase in the source distance estimate (with an uncertainty margin of 15–20% of its own) inflates projected linear scales by a factor 1.18, while quantities such as mass and luminosity, determined from broadband flux and column density, would scale by the distance squared, that is, a factor 1.39. With $d$=120 pc, the angular distance of 5.3″ between the two submillimeter sources A and B corresponds to a projected distance of 636 au. The exact three-dimensional geometry of the two sources, their potential disks and the filament structure is unknown. Based on a total of 515 spectral lines from 54 molecular species detected in a Submillimeter Array (SMA) spectral survey, Jørgensen et al. (2011) derived LSR velocities of sources A and B of $+3.2$ and $+2.7$ $\mathrm{km}\,\mathrm{s}^{-1}$, respectively (throughout this work, IRAS 16293A will be abbreviated as “source A” and IRAS 16293B as “source B”). Different centroid velocities have been reported based on individual molecular lines: for source A, $+3.8$ $\mathrm{km}\,\mathrm{s}^{-1}$ from C17O, C34S (Favre et al. 2014), $+3.6$ $\mathrm{km}\,\mathrm{s}^{-1}$ from HCN (Takakuwa et al. 2007); for source B, $+3.4$ $\mathrm{km}\,\mathrm{s}^{-1}$ from CH3OCHO, H2CCO (Pineda et al. 2012). Differences in the derived centroid velocity may stem from different angular resolution of the data sets, the different molecular tracers used, and/or varying optical depth of individual transitions of the same species, so that each traces a somewhat different ensemble of gas. All of the centroid velocities listed above fall well within the $\sim$3–4 $\mathrm{km}\,\mathrm{s}^{-1}$ wide distributions of $V_{\mathrm{lsr}}$ values fitted by Jørgensen et al. (2011).

Caux et al. (2011) used single-dish observations to estimate the masses of the two protostellar sources. Assuming that the line broadening of the profiles is due to infalling motions, these authors arrive at $\sim$0.8 $M_{\sun}$ for source A and 0.1 $M_{\sun}$ for source B. The same mass of source A was found by Bottinelli et al. (2004), also assuming infalling motion to explain the observed line broadening, but using interferometric observations in which source A is spatially separated from source B. Alternatively, the line broadening could be caused by rotating motion instead of infall, as assumed for source A by Pineda et al. (2012), who used observations of methyl formate (CH3OCHO) and ketene (H2CCO) with a 2.2″$\times$1.0″ beam to infer Keplerian rotation around a central object of 0.53 $M_{\sun}$. Oya et al. (2016) found a mass varying between 0.5 and 1.0 $M_{\sun}$, depending on the assumed inclination and centrifugal barrier radius. These studies show that the source masses are in the regime of low-mass protostars, but the various methods yield masses differing by up to a factor of two.

The binary protostellar system IRAS 16293 is embedded in an envelope with a radius of (6–8)\times$$10^{3} au (e.g., Schöier et al. 2002; Crimier et al. 2010). Rotating, infalling motions of this envelope have been inferred from spectral line profiles (e.g., Menten et al. 1987; Zhou 1995; Ceccarelli et al. 2000; Schöier et al. 2002; Takakuwa et al. 2007). In a detailed velocity model posited by Oya et al. (2016), the large-scale envelope could transition into a Keplerian disk within the centrifugal barrier around source A, as indicated by observations presented in Favre et al. (2014). The latter study found rotation on 50–400 au scales, which could not be explained by simple Keplerian rotation around a point-mass, but needed to take into account the extra material of the enclosed mass at these scales.

At least two major outflows have been observed from IRAS 16293A: an east-west bipolar outflow (Yeh et al. 2008) and a northwest-southeast outflow pair (Kristensen et al. 2013; Girart et al. 2014). On larger scales of $>$5000 au, a northeast outflow has also been observed (Mizuno et al. 1990; Stark et al. 2004), the origin of which is speculated to be IRAS 16293A, based on high-resolution CO images revealing a collimated structure near source A (see Fig. 1). Combining single-dish observations of a broad range of HCO+ spectral lines ($J$=1–0 up to $J$=13–12), Quénard et al. (2018) conclude that the spatially unresolved spectral line profiles are dominated by outflow contributions, explained by these authors using a model of the northwest-southeast outflow emanating from source A. In contrast, any signs of outflows associated with IRAS 16293B have long escaped detection. Loinard et al. (2013) argued that the blueshifted emission found southeast of source B is from a young monopolar outflow from source B, but Kristensen et al. (2013) demonstrated, using the same data, that it could be a bow-shock from the northwest-southeast outflow from IRAS 16293A. Oya et al. (2018) presented SiO velocity maps indicative of a pole-on pair of outflows from source B, but these authors also recognize that interaction with the outflow from source A is not ruled out as a possible scenario.

The quadruple outflow structure from source A has led to speculations that source A itself is a multiple system. Indeed, continuum observations at centimeter waves resolved IRAS 16293A into two (Wootten 1989; Chandler et al. 2005) or even three components (Loinard et al. 2007; Pech et al. 2010), with a separation of up to $\sim$0.5″. The 0.5″ resolution maps from the ALMA PILS program appear to indicate a singularly peaked source at the position of source A, both in continuum (Jørgensen et al. 2016) and in optically thin C17O (Jacobsen et al. 2018). However, the integrated C17O 3–2 intensity is affected by contamination from more complex molecules in the dense, warm regions close to the two protostellar sources, which modifies the apparent morphology of emission integrated over channels bracketing the C17O line frequency (see Fig. 8 in Appendix A). Moreover, the 0.87 mm dust continuum is optically thick in the disk domains, which could hide any intrinsically multipeaked structure (Calcutt et al. 2018). In contrast, the PILS maps of various isotopologues of methyl cyanide in Calcutt et al. (2018), for which the emission is integrated over a moving velocity interval to avoid incorporating contributions from other species, do show two clearly separated emission peaks. The methyl cyanide peak positions are consistent with those of radio continuum (3.6 cm) positions A1 and A2 from Loinard et al. (2007) after correction for proper motion (Pech et al. 2010). The slight offset of $\sim$0.2″ is likely due to pointing errors and uncertainties in the proper motion correction coefficients (Calcutt et al. 2018).

Material spanning the interbinary region between IRAS 16293A and B was first identified in millimeter observations from the BIMA array111BIMA stands for Berkeley-Illinois-Maryland Association, which operated and funded the nine-telescope array at Hat Creek in California until 2004, when its hardware was relocated and merged into the Combined Array for Research in Millimeter-wave Astronomy (CARMA). by Looney et al. (2000), who interpreted it as a circumbinary envelope. The same bridge material is also clearly detected in sensitive Atacama Large Millimeter/Submillimeter Array (ALMA) observations (Pineda et al. 2012; Jørgensen et al. 2016). Observations of polarized dust emission indicate that the magnetic field in the filamentary bridge is oriented along the long axis of the bridge (Rao et al. 2009, 2014; Sadavoy et al. 2018). Jacobsen et al. (2018) (hereafter J2018) show that a dust filament model can match the observed spatial emission from high-resolution submillimeter continuum emission and C17O 3–2 gas line emission. Similar structures observed in other multiple protostellar systems are interpreted as fragmentation of large-scale circum-multiple envelopes (Lee et al. 2017) or (Keplerian) rotating disks that have become locally gravitationally unstable (Tobin et al. 2016a; Fernández-López et al. 2017; Dutrey et al. 2014, 2016).

The exact evolutionary stage of the two protostars is unknown, partly because of the different inclination angles of sources A (edge on; Pineda et al. 2012; Oya et al. 2016) and B (face on; Oya et al. 2018). Given the outflow and infall signatures near source A, its status as a protostar is firmly established. In contrast, the apparent quiescence of source B has led to speculations that it is not a young protostar, but rather harbors a late-stage T Tauri disk (Stark et al. 2004). This interpretation was later called into question when spectral line signatures of infalling material were observed toward source B (Chandler et al. 2005; Pineda et al. 2012; Jørgensen et al. 2012). In conclusion, it is as yet unknown if the two protostars are at the same evolutionary stage, or that one may be more evolved than the other.

The outflows from source A, along with the striking arc of dust and gas connecting sources A and B, make IRAS 16293 a very complex system, where the observed gas line emission can only be explained by a combination of multiple physical components. We aim to disentangle these structures in our new observations and to map the physical origins of the observed molecular emission lines.

This paper is based on observations from the ALMA band 7 segment of the Protostellar Interferometric Line Survey (PILS222http://youngstars.nbi.dk/PILS) targeting IRAS 16293. The many thousands of line detections in the PILS data set facilitate the discovery of complex organic molecules ($>$6 atoms, including carbon; Herbst & Van Dishoeck 2009) and the study of their interstellar chemistry. Relatively weak lines of complex molecules are most easily separable in spectra extracted toward spatial positions near the kinematically simple and spatially compact source B, with narrow ($\sim$1 $\mathrm{km}\,\mathrm{s}^{-1}$), single-moded spectral line shapes. Most analysis based on the PILS data cubes has therefore been restricted to spatial positions close to the two protostars (e.g., Jørgensen et al. 2016; Coutens et al. 2016; Lykke et al. 2017; Ligterink et al. 2017; Fayolle et al. 2017; Coutens et al. 2018). In contrast, in this work, we fully exploit the spatial dimensions of the PILS data set, but focus on a few restricted frequency ranges containing lines of well-known, simple molecular species. The aim is to study the kinematic structure of the binary system and the intervening and surrounding gaseous material.

We summarize the characteristics of the observational data in Sect. 2, and describe the selection of molecular tracers and the morphology and dynamics observed in each of them in Sect. 3. Section 4 presents analysis of the physical characteristics of the bridge filament between the two protostars and the kinematics of outflow motions. Discussion of the results is given in Sect. 5, and Sect. 6 summarizes the main conclusions.

2 Observations

The PILS program was conducted with ALMA between June 2014 and May 2015, using both the main array of 12-m dishes and the shorter baselines available in the array of 7-m dishes in the Atacama Compact Array (ACA). Its end product comprises a three-dimensional data set covering the uninterrupted spectral range 329–363 GHz (all in ALMA band 7) in spectral channels of 0.244 MHz in width. The sky area of $\sim$15″ in diameter is covered in a single pointing with the primary beam of the 12-m dishes, sampling the two binary components of IRAS 16293 and their surroundings with a synthesized beam FWHM of 0.5″. Self-calibration of the visibility phases was applied, based on the continuum signal. While the main array provides the unprecedented sensitivity, which is roughly uniform at $\sim$8 mJy beam*-1* channel*-1* across the full spectral range, the addition of ACA baselines ensures that large scale emission up to $\sim$13″ is recovered. We refer to Jørgensen et al. (2016) for additional details of observing conditions, data characteristics, and the data processing strategy.

3 Results: observed dynamics and morphology

The morphology of the dust continuum emission in IRAS 16293 as observed with ALMA at 0.5″ resolution at 0.87 mm (Jørgensen et al. 2016; J2018) can be broken up into three distinct components: a nearly circular structure that hints at a face-on disk surrounding source B; an elliptically shaped structure related to the inclined disk-like structure around source A; and a ridge of material stretching from B to A and even beyond, to the southeast of A. The centroid coordinates of the continuum peaks of components A and B, measured using two-dimensional Gaussian fits on the PILS 0.826–0.912 mm continuum image, are listed in Table 1. These positions of the 0.87 mm continuum measured from our PILS data products are consistent with the 1.15–1.30 mm continuum peak locations reported by Oya et al. (2018) to within 0.05″, that is, one tenth of the synthesized beam size of either observation.

We select spectral signatures of C17O, H13CN, H2CO, H213CO, C34S, SiO, 29SiO, and C2H. The motivation for selecting these tracer transitions is that their emission is spatially extended across the inter- and circumbinary region, the molecules are chemically simple, and are chosen because they probe a range of densities, temperatures (see $E_{\mathrm{up}}$ and $n_{\mathrm{crit}}$ in Table 2) and physical processes such as shocks. C18O 3–2 largely follows the spatial distribution of C17O (J2018), and the former is not further discussed in this paper. Although we recognize the additional value of studying other tracer molecules also covered in the 329–363 GHz PILS data, the selection in this paper is deliberately restricted to its current scope. The rest frequencies for the relevant transitions of these species are listed in Table 2. For the systemic velocities of each of the protostars, we adopt the statistically averaged values based on the multiline analysis by Jørgensen et al. (2011). Velocity channel maps of each of these species are shown in Figs. 16–23, with a velocity range in $V_{\mathrm{lsr}}$ between 0.0 and $+6.5$ $\mathrm{km}\,\mathrm{s}^{-1}$, bracketing the systemic velocities of $+3.2$ and $+2.7$ $\mathrm{km}\,\mathrm{s}^{-1}$ for sources A and B in the IRAS 16293 binary (Jørgensen et al. 2011). As our aim is to study material nearby, but not in the dense protostellar sources, two masks are placed on the disks of each source in all maps presented in this work. This procedure ensures that high-intensity, partially (self-)absorbed spectral line signals do not skew the intensity scaling on the integrated intensity maps and the calculation of the weighted velocity (moment 1) maps. Contamination of the selected molecular line transitions by nearby transition of other species are discussed in detail in Appendix A. Spectral line shapes deviate from Gaussian profiles in some cases, particularly in directions where line optical depth is high, and/or outflowing or infalling gas motions contribute to the gas emission.

We construct velocity maps of the selected species by considering only channels with $V_{\mathrm{lsr}}$ in the range [$-4$,$+11$] $\mathrm{km}\,\mathrm{s}^{-1}$ (except for C2H, where we choose [$-1$,$+7$] $\mathrm{km}\,\mathrm{s}^{-1}$ to limit contamination by other species), and calculating an intensity-weighted mean velocity for each pixel, excluding all flux density values below 40 mJy beam*-1* (i.e., roughly 4–6 $\sigma$ depending on the noise level in the particular section of frequency coverage). The result is shown in Fig. 2. The velocity structure of each tracer is also displayed in the form of contour maps for three integrated ranges of velocity (red, systemic, blue) in Fig. 3, and split out for each 0.2 $\mathrm{km}\,\mathrm{s}^{-1}$ spectral channel in Appendix B.

As seen in Fig. 2, the C17O gas has an intensity-weighted mean velocity restricted to $\pm$ 1 $\mathrm{km}\,\mathrm{s}^{-1}$ of the systemic velocity in the bridge region between source A and B. The bridge is therefore kinematically quiescent, that is, we find no evidence that the bridge is part of an outflow motion nor that it is a rotating structure. The spatial distribution of C17O coincides with the dust bridge traced by submillimeter continuum, while none of the other (higher density) tracers in this study have intensity peaks that are cospatial with the bridge, at 0.5″ resolution (60 au in projection, at the distance of the source). The bridge morphology is also recovered in the PILS C18O 3–2 map (J2018), while the 13CO 3–2 starts looking markedly different due to high line optical depth (J2018). The C34S 7–6 distribution, discussed in more detailed in Sect. 4.3, may seem to overlap partly with the dust and C17O, but we do not regard it as tracing the bridge.

A tentatively axisymmetric structure is seen in the channel maps of H2CO and H13CN (Figs. 2, 3, 17, 18): two arcs emanating from source A, one on the NE side at velocities between $+2$ and $+3.5$ $\mathrm{km}\,\mathrm{s}^{-1}$, and its counterpart on the SW side at velocities between $+4$ and $+5.5$ $\mathrm{km}\,\mathrm{s}^{-1}$ (e.g., Fig. 17). These arcs bracket the axis defined by the NW outflow seen in CO 6–5 by Kristensen et al. (2013), and are kinematically symmetric around $V_{\mathrm{lsr}}$ = $+3$ $\mathrm{km}\,\mathrm{s}^{-1}$. Examining the molecular tracers considered in this work: this velocity gradient (roughly perpendicular to the main axis of the bridge) is spatially more compact in H13CN than in H2CO; there may be evidence of it in C34S; but it is absent in C17O, SiO, and C2H (Fig. 2, Fig. 3). We therefore conclude that the striking symmetry observed in H2CO is not the result of a bulk rotation of the gas about the axis of the bridge or the axis of the northwest outflow from source A.

While none of the tracer molecules listed above have signals bright enough to study faint line wings at velocity offsets more than 4 $\mathrm{km}\,\mathrm{s}^{-1}$, the main isotopologue of CO provides sufficiently bright signal to probe higher line-of-sight velocities. In fact, in Fig. 4 we use only high-velocity ($>$5 $\mathrm{km}\,\mathrm{s}^{-1}$) channels of CO 3–2 emission, thereby avoiding complications in interpreting extremely optically thick emission at lower velocities. The western lobe of the east-west outflow pair is prominently visible in CO 3–2 in a cone-shaped distribution, with a much smaller region of blueshifted counterpart appearing on the eastern side of source A. At somewhat coarser angular resolution (1.5″ beam, using the SMA), Yeh et al. (2008) have previously studied the east-west outflow from source A in CO 2–1 and 3–2. Compared with their work, our 0.5″ resolution map in CO 3–2 reveals a more cone-like shape of the westward, redshifted outflow lobe, and there is an apparent acceleration taking place with increased distance from source A. See Sect. 4.3 for a discussion of the structures described here. There is also bright CO 3–2 emission with an emission peak 1″ south of source B, with velocities ranging up to 10 $\mathrm{km}\,\mathrm{s}^{-1}$ blueshifted with respect to source B. Its location is consistent with the position ‘b2’ described in Yeh et al. (2008). In our map (Fig. 4), this blueshifted knot is clearly morphologically connected to source B.

We observe that the SiO velocity map in Fig. 2 traces both sides of the NW-SE outflow pair extending up to $\sim$8″ from source A, blueshifted on the northwest side and redshifted on the southeast side. The lower optical depth line of 29SiO only shows emission on the redshifted side southeast of source A. In contrast, a redshifted southeast outflow appears absent in the high-temperature tracer CO 6–5 (Kristensen et al. 2013), although it is seen in colder gas (CO 2–1; Girart et al. 2014). Likewise, the PILS CO 3–2 data (Fig. 4) do show redshifted emission about 6″ southeast of source A, and it is very prominent in the SiO velocity map (Fig. 2) at V_{\mathrm{lsr}}$$\geq $+5$ $\mathrm{km}\,\mathrm{s}^{-1}$, as well as in H2CO and H13CN. In addition, a significantly redshifted patch of emission is observed in SiO about 3$\aas@@fstack{\prime\prime}$5 due west of A (and 4″ south of B), overlapping with the ‘kink’ location in the CO outflow, as highlighted in Fig. 5 (see also Sect. 4.3). SiO 8–7 is also seen prominently 1–2″ north and south of source B, with mixed red and blue velocity components along the line of sight (Figs. 3, 21), as well as in the fainter 29SiO line (Fig. 22). This structure is consistent with that observed in SiO 6–5 by Oya et al. (2018).

Finally, a filamentary structure connected to source B with a position angle of $\sim$170° is visible in C2H, but not in any other tracer. Its velocity structure is flat, restricted to $V_{\mathrm{lsr}}$ values between $+2$ and $+4$ $\mathrm{km}\,\mathrm{s}^{-1}$ (see Sect. 4.3 for further discussion). The redshifted (¿$+5$ $\mathrm{km}\,\mathrm{s}^{-1}$) emission in the velocity map of C2H (Fig. 2) is due to a bright CH3CN transition at 349.393 GHz. In positions near source A, the CH3CN line also contaminates channels closer to the systemic velocity of C2H (see Appendix A). In other positions in the map with integrated intensity of C2H above 0.07 Jy beam*-1* $\mathrm{km}\,\mathrm{s}^{-1}$, the spectral profile shows a characteristic structure with two roughly equal intensity peaks spaced by 1.4 MHz, the separation between the $F$=4–3 and $F$=3–2 hyperfine components (see Table 2, Fig. 15 in this work, and Fig. 4 of Murillo et al. 2018). The interpretation of the observed morphology and velocity structure is addressed in Sect. 4.3 and Sect. 5.

4 Analysis

4.1 Radiative transfer modeling of the interbinary bridge

The bridge between protostars A and B, as traced by the submillimeter dust continuum and cold C17O gas (Sect. 3), is kinematically quiescent (Sect. 3) and its central velocity at $V_{\mathrm{lsr}}$=+3 $\mathrm{km}\,\mathrm{s}^{-1}$ is consistent with the velocities of the individual protostars ($+2.7$ and $+3.1$ $\mathrm{km}\,\mathrm{s}^{-1}$). It is therefore unlikely to be associated with any of the outflows, and it is treated as a separate entity in this work. Besides C17O, none of the other selected molecules match, simultaneously, the morphology of the dust bridge filament and the narrow line-of-sight velocity distribution of C17O.

In principle, an assessment of density and temperature conditions throughout the inter- and circumbinary region of IRAS 16293 could be set up using spectral line intensity ratios of various molecules and transitions highlighted in this work (see Table 2). However, given the multitude of partially resolved dynamical components, it is not at all certain that line intensities measured toward a particular position emanate in the same gas for one transition of a species and another transition of another species. We therefore refrain from embarking on such an analysis.

To study the quiescent interbinary bridge, we use a curved cylinder filament model, spanning 636 au between and enveloping the two sources in an arc-like structure. The density profile of the bridge arc is a power-law which decreases with radius to the power $-0.25$ and is independent of the distance from either protostellar source. The enclosing envelope extends to beyond our field of view, and has a purely spherical geometry, with a density profile decreasing radially with a power of $-1.7$ and a density plateau imposed at radii within 600 au. This is the same filament model used in J2018, who show that, when using luminosities of 18 $L_{\odot}$ and 3 $L_{\odot}$ for sources A and B, respectively, the model yields a qualitative match with the morphology and peak emission levels seen in the 0.87 mm dust continuum emission and C18O/C17O gas line emission. The fiducial model adopted in this paper is equivalent to ‘Rotating Toroid model 1’ in J2018. The temperature structure of each cell in the model is derived self-consistently, based on the adopted luminosities of both sources A and B (Fig. F.4 of J2018). We do not scale distances in the model to the newly published distance of IRAS 16293 (see Sect. 1), to allow inclusion of and direct comparison with results from the modeling work of J2018. An 18% increase in source distance would have little impact on the model outcomes, compared with other assumptions with considerably larger uncertainties, such as dust opacity and sublimation thresholds (see below).

The line radiative transfer code LIME (Brinch & Hogerheijde 2010) is used to obtain synthetic line emission cubes of the H2CO 51,5–41,4 and H13CN 4–3 gas line emission, analogous to the strategy used in J2018 for CO and its isotopologues. We use an abundance jump model to emulate freeze-out onto the dust grains. When the dust temperature is below the sublimation temperature, the gas phase abundance is decreased by a depletion factor. See Table 5 for more information on the individual molecules. Uncertainties in the sublimation thresholds could be up to $\sim$30% (see references in Table 5). This study only considers the inter- and circumbinary region and excludes the warm ‘disk regimes’ at $\gtrsim$50 K (cf. the model temperature structure in Fig. F.4 of J2018). Therefore, adjusting the adopted sublimation temperature by 30% would only impact the modeled emission morphology of CO species, and only those that are optically thin (C17O, C18O), while those of H2CO and H13CN would be left unaffected. The number density of H2 is above 10^{8}\,$$\mathrm{cm}^{-3} in the modeled bridge filament, sufficiently high to validate the occurrence of depletion of gas-phase molecules onto dust grain surfaces. Since we are investigating the appropriateness of a static filament structure to the observation of the quiescent bridge only, we keep the velocity structure of the model static, with only random velocity dispersions included, such as from turbulence. We fix the built-in velocity dispersion parameter (Doppler $b$) in LIME to 1 $\mathrm{km}\,\mathrm{s}^{-1}$. This value is the 1/e half-width, equivalent to a FWHM of 1.4 $\mathrm{km}\,\mathrm{s}^{-1}$ for a Gaussian line profile, which is representative of the measured line widths in the region under study. The only free parameter is the molecular abundance. Twenty different abundance values are run for each molecule in the range given in Table 5. After LIME produces a synthetic line emission cube, it is convolved with a $0.5^{\prime\prime}\times 0.5^{\prime\prime}$ beam in the image analysis software package MIRIAD. We refer to J2018 for further details of the model definition and the radiative transfer approach.

In Fig. 6, integrated intensity maps of the radiative transfer calculations described above are juxtaposed with their observed counterparts. The latter are produced by integrating over the narrow velocity range $V_{\mathrm{lsr}}$=[+2.0,+4.0] $\mathrm{km}\,\mathrm{s}^{-1}$ (i.e., all channels with contours in black in Figs. 16–18 in the Appendix), with the aim of filtering out dynamically active gas associated to the various outflow components which are not included in the density distribution of the model (J2018). As already concluded by J2018, the model in C17O provides a qualitative match to the morphology of the bridge connecting the two protostars, although the observed C17O filament is laterally $\sim$3 times wider than the modeled curved filament. Absolute flux values produced by the model are a factor $\sim$2 higher than observed. H2CO is also present and detected in its 51,5–41,4 transition across the bridge area outlined by the continuum contours, but its observed morphological shape is different from that of C17O and the continuum. The difference is partly due to the inclusion of outflow-impacted gas, despite the attempt to exclude this by using a narrow velocity integration range. The modeled H2CO emission extends radially over many hundreds of au surrounding the binary, an effect of H2CO molecules being present in the gas phase in the spherical envelope component of the model. The superposed density structure of the disks-bridge-envelope model is dominated by that of the spherical envelope model at such large radii of $\sim$700 au (from the midpoint between A and B). Outside this radial distance of 700 au, the number density of the modeled envelope component drops below $10^{7}$ $\mathrm{cm}^{-3}$ (the outermost gray, circular contour in Fig. F.4 of J2018). Along with radially decreasing temperatures, this leaves conditions insufficient to excite o-H2CO 51,5, with a critical density of 7.5$\times$107 $\mathrm{cm}^{-3}$ and an upper energy level of $E_{\mathrm{up}}$/$k$=62.5 K (Table 2). In addition to tracing the bridge filament, the H2CO model map shows enhanced H2CO emission up to $\sim$2″ above the disk of source A (northwest and southeast of A), where the temperature is sufficiently high to sublimate H2CO from dust grains (cf. 55 K contour in Fig. F.4 of J2018). Finally, the observed morphology of H13CN in the inter-binary region does not follow the shape of the bridge filament. Most of the H13CN 4–3 emission shown in Fig. 6, bottom left, is due to contamination from outflow-impacted components (despite our effort to exclude these using a narrow velocity integration interval), none of which are included in the physical model. As expected from its high critical density, the radiative transfer model does not populate the $J$=4 level of HCN species beyond the disk domains. The exception is a small region immediately above the disk plane, where T_{\mathrm{kin}}$$>100 K (again, see Fig. F.4 of J2018), where the observed enhancement in H13CN emission is reproduced by the model.

In conclusion, the bridge morphology is observed in dust, C17O, and partly also in o-H2CO, and the physical model constructed by J2018 roughly reproduces the morphology of the dust continuum and C17O emission. The partial overlap between the observed H2CO emission and that produced by the model (Fig. 6, middle panels) should not be over-interpreted. As shown in Figs. 2 and 18, the H2CO observed in the direction of the bridge filament domain is largely at line-of-sight velocities that are offset from the quiescent C17O bridge filament at the systemic velocity. In fact, the lower optical depth tracer H213CO in the same transition (Fig. 19) shows no detectable emission in the bridge domain in the $V_{\mathrm{lsr}}$ range from +2.0 to +4.0 $\mathrm{km}\,\mathrm{s}^{-1}$. This means that the observed H2CO emission is not tracing the bulk gas in the quiescent bridge filament, but rather surface layers of dynamically stirred components. In our interpretation, this indicates that the number densities in the modeled bridge component are higher than in reality. A more realistic model would not produce any H2CO emission from the bridge component. H13CN 4–3 emission is not seen outside of the disk domains in the modeled map, and does not trace the bridge morphology in the observed map (Fig. 6, bottom panels). The explanation for the lack of H2CO and H13CN emission in the bridge filament is provided by the excitation balance governed by density and temperature. Looking at the critical densities and $E_{\mathrm{up}}$ values in Table 2, the upper levels of the relevant H13CN and H2CO transitions do not get populated sufficiently. The freeze-out abundance drop in the model (Table 5) is not the cause, which we confirm with a separate radiative transfer model run in which freeze-out is neglected completely, keeping all molecules in the gas phase. In the resulting emission maps (shown in Fig. 24 of Appendix C), there is still no H13CN 4–3 emission outside of the disk regions. Without freeze-out, the model o-H2CO 51,5–41,4 emission map is dominated by the contribution from the cold, extended envelope (Fig. 24). In contrast with H13CN, H2CO still does show a small intensity enhancement in part of the bridge arc region, probably reflecting its slightly lower critical density when compared with HCN 4–3. The total line-of-sight optical depth of the superposed components, particularly for H13CN, becomes much higher in the ‘no freeze-out’ case, as expected when all H13CN molecules remain in the gas phase even in the colder, extended envelope that surrounds the binary.

The fiducial model (RTM1 from J2018) is not necessarily the only and best fit to the data. The only observational constraints that were taken into account were the multiwavelength SED and the distribution of CO 3–2 isotopologues and submillimeter dust emission; the J2018 work did not include an assessment of molecular lines tracing a range of densities. The number density at the axis of the modeled bridge arc is 7.5$\times$108 $\mathrm{cm}^{-3}$. The density structure that was adopted by J2018 was mainly driven by the wish to reproduce the dust emission morphology. However, dust emission strength depends on several factors which are all assumed to be fixed in the model setup: dust opacity, gas-to-dust ratio, and a dust temperature being coupled rigidly to that of the gas. For example, if a different grain size distribution is assumed, higher submillimeter grain opacities could arise, lowering the peak (column) density required to match the observed dust emission. Similar arguments hold for a lower gas-to-dust ratio and a higher dust temperature. A combination of these effects can easily bring the model peak number density down by a factor of 10–100 to well below the critical densities of the CS, H2CO and HCN transitions selected in this work (Table 2). In this work, a suite of molecular tracers with different critical densities are studied, among which C17O 3–2 (4$\times$104 $\mathrm{cm}^{-3}$), o-H2CO 51,5–41,4 (7.5$\times$107 $\mathrm{cm}^{-3}$) and H13CN 4–3 (1.2$\times$108 $\mathrm{cm}^{-3}$). Of these, C17O 3–2 is the only species that traces the bridge morphology observed in dust continuum (see above), and, additionally, the distribution of two lines of 34SO2 observed in the PILS data cube (critical densities of 3.8$\times$107 $\mathrm{cm}^{-3}$ and 8.3$\times$107 $\mathrm{cm}^{-3}$) is also confined to the disk domains of the two protostellar sources (Drozdovskaya et al. 2018). Put together, this indicates that, indeed, the real density of the bridge is likely between 4$\times$104 $\mathrm{cm}^{-3}$ and 3$\times$107 $\mathrm{cm}^{-3}$.

Besides the low-density tracers C17O 3–2 and 0.87 mm dust, all emission observed from other molecular species in the vicinity of the bridge domain can be ascribed to components with a measurable velocity along the line of sight, or their impact on more quiescent regions (see below in Sect. 4.3).

4.2 Fragmentation of the bridge

It is conceivable that both protostellar sources A and B formed from the bridge filament through gravitational instability. In this context, we follow the analytic approach by Ostriker (1964) and treat the bridge as an isothermal cylinder. In this idealized scenario (far more simplified than the self-consistently derived temperature structure modeled by J2018 and used in our Sect. 4.1), the thermal pressure in the cylinder keeps it from collapsing radially if the mass per unit length, $M/L$, is below a certain threshold,

[TABLE]

This quantity depends linearly on kinetic temperature $T_{\mathrm{kin}}$ and is inversely proportional to the mass of the typical particle $\mu m_{0}$, with $m_{0}$ the mass of a proton. It is the equivalent of Jeans mass in cylindrical geometry. Using $\mu=2.0$ for H2 and three different temperature values of [25, 30, 40] K, we arrive at a stability threshold of [0.14, 0.17, 0.23] $M_{\sun}$ per 600 au (the approximate projected length of the bridge).

We compare the mass per length stability criterion derived above with an observationally derived mass of bridge material, calculated as follows. We take the average C17O 3–2 intensity of 600 mJy beam*-1* (integrated over the [+2,+4] $\mathrm{km}\,\mathrm{s}^{-1}$ velocity range), and use RADEX (Van der Tak et al. 2007) non-LTE radiative transfer calculations to convert to a C17O column density. Such a calculation requires a kinetic temperature and a number density of collision partners as input. To represent a natal, isothermal, cylindrical filament, we choose a temperature between the minimum temperature of the modeled bridge arc (24 K) and the average temperature in the bridge domain (41 K) (see Sect. 4.1 and J2018). For number density, $n_{\mathrm{H_{2}}}$, we adopt $10^{5}$ $\mathrm{cm}^{-3}$, which is between the critical density of C17O 3–2 and those of HCN 3–2, H2CO 51,5–41,4 and SO2 174,14–173,15 (Table 2, Sect. 4.1). With $T_{\mathrm{kin}}$=[25, 30, 40] K, we obtain a C17O column density of [6.2, 5.4, 4.8]\times$$10^{15} $\mathrm{cm}^{-2}$, only mildly sensitive to temperature. Using a 16O/17O ratio666The adopted 16O/17O ratio is consistent with that measured in H2CO toward IRAS 16293B by Persson et al. (2018): (2.0–3.0)\times$$10^{3}. of 2.005$\times$103 (Wilson 1999) and CO/H2=$10^{-4}$, this corresponds to $N$(H2)=[1.24, 1.08, 0.96]\times$$10^{23} $\mathrm{cm}^{-2}$. Taking a rough half-power width of the dust bridge of 200 au, the total mass encompassed by the bridge ‘cylinder’ area of 200$\times$600 au between sources A and B is [0.0056, 0.0049, 0.0044] $M_{\sun}$.

On the other hand, integrating the mass of all dust and gas in the bridge component of the model described in Sect. 4.1 yields a tenfold higher total mass of 0.055 $M_{\sun}$. Either way, the total bridge mass values are below the stability threshold at the assumed temperatures of 25, 30, and 40 K. The discrepancy between these two bridge mass calculations is not surprising, given the assumptions in the dust model from J2018. The dust opacity model, from Ossenkopf & Henning (1994), uses a micron-sized dust distribution. In the hypothetical case that grain growth has already taken place and that the true size distribution of dust grains in the bridge peaks at sizes of 100-1000 $\mu$m, then opacities at 345 GHz (expressed per gram of material) can be up to a factor of 10 larger, which would mean that the mass determination from our radiative transfer calculations is overestimated by up to an order of magnitude. However, there is no evidence to support that the size of the dust grains in the interbinary bridge of IRAS 16293 would deviate from the roughly micron-sized populations representative of the unprocessed interstellar medium. Moreover, the canonical gas-to-dust mass ratio of 100 was used, and any deviation from this value would also change the total mass estimate of the modeled filament.

In the case presented above in which the filament mass is below the Ostriker stability criterion, the filament would be supported against radial collapse. This scenario, in fact, is the only way in which a filament may persist long enough for it to fragment along its vertical (length) direction (Inutsuka & Miyama 1992) and form additional dense cores. However, the total mass budget available in the current bridge filament (see above) seems insufficient to form even a brown dwarf. Again in an idealized nearly isothermal scenario, Inutsuka & Miyama (1997) have shown that separation between fragmented cores within a radially stable filament are on the order of four times its width. This ratio is consistent with the bridge observed in IRAS 16293, with an observed diameter of $\sim$200 au, and separation between the cores of at least 600 au, knowing that the vector connecting sources A and B may be inclined with respect to the plane of the sky. While the bridge filament present little to no gas motion along the line of sight ($\lesssim$1 $\mathrm{km}\,\mathrm{s}^{-1}$), observations presented in this work and in other papers in the literature do not rule out gas flow vectors in the plane of the sky along the length of the filament, possibly regulated by the highly ordered magnetic field in the bridge (Sadavoy et al. 2018).

Finally, we emphasize that, if indeed the bridge filament has been host to sources A and B, the filament must have been considerably more massive in the past. After adding even half of the current mass of sources A and B (0.5 $M_{\sun}$) to the current bridge filament mass calculated above, it is lifted to a factor few above the gravitational instability threshold. We emphasize that the interpretation of stability of the filament, sketched in this section, only holds for the filament in its current state. Radial collapse in an earlier evolutionary stage, with a different mass distribution and under different temperature conditions, is not ruled out. To conclude, the observed filament is obviously more complex in nature than a static, cylindrically symmetric shape characterized by a single temperature. A true assessment of its gravitational stability should include, for example, turbulent and magnetic support. Any additional mechanisms to provide support against gravitational collapse would only increase the stability threshold, and the conclusions drawn in this section about the dynamical stability of the filament would thus remain unchanged.

4.3 Kinematics

Besides the bridge filament (Sects. 4.1, 4.2), almost all other morphological components discussed in this paper are dynamic and can be attributed to protostellar outflow activity. Distinct components are discussed in this subsection and are summarized in Fig. 1 and Table 2.

The CO 3–2 velocity map in Fig. 4 clearly traces the western lobe of the known east-west outflow pair. However, its eastern counterpart and the separate southeast-northwest pair are all much less pronounced, possibly because the eastern outflow is intrinsically dimmer (see also Yeh et al. 2008) and the southeast-northwest pair may have lower line-of-sight velocity. In fact, the southeast outflow lobe from source A is detected in larger spatial extent at more modestly redshifted velocities ($V_{\mathrm{lsr}}$=[+6,+8] $\mathrm{km}\,\mathrm{s}^{-1}$) than those shown in Fig. 4. The line-of-sight velocity in the western lobe increases with distance from source A, implying an acceleration with time. This is consistent with a centrally driven outflow that would naturally accelerate as the density of the ambient envelope, and therefore the pressure, decreases with radius (e.g., Moriarty-Schieven & Snell 1988; Arce et al. 2007) . In addition, the CO outflow in Fig. 4 exhibits a relatively collimated morphology up to $\sim$3″ ($\sim$350 au in projection) from source A, at which point it abruptly fans out into a cone-like structure. This position coincides with a marked enhancement in SiO 8–7 emission (Fig. 5), possibly related to a simultaneous acceleration of the outflow (leading to SiO enhancement through dust grain erosion) and decrease in ambient pressure (leading to a wider angle CO structure). The position of this kink may be related to the radial extent of the rotating ‘inner envelope’ described by Oya et al. (2016). However, these authors infer a radius of the inner envelope of 200 au, whereas a density drop would be needed beyond 300 au to explain the location of the kink in the CO outflow morphology described above. Alternatively, the opening angle of the western outflow may be restricted by higher pressure inside the quiescent interbinary bridge described in Sect. 4.1 of this work. Its radius (as observed in dust continuum and C17O, Fig. 2, top left; Fig. 6, top left) appears to be more congruent with the 350 au radial separation from source A.

Closer in to the source, at about 120 au projected distance west of source A, a H2O maser spot detected by Dzib et al. (2018) seems to coincide with the interface between the dense, rotating, disk-like structure surrounding source A (blue on the west, red on the east) and the emerging redshifted outflow in the western direction (see Figs. 1, 4). The direction of proper motion of the maser spot (position angle $\sim$290°, Dzib et al. 2018) roughly follows the flow direction of the western outflow as seen in CO 3–2 PILS data (Fig. 4). This overlap supports the suggestion made by Dzib et al. (2018) that the redshifted maser spot may be due to a launching shock at the base of the outflow. A similar argument may hold for a blueshifted maser spot detected to the east of source A, although the spatial association is more ambiguous than for the redshifted counterpart (Dzib et al. 2018).

In contrast to 12C16O, the other molecular tracers selected from the PILS data set mainly highlight the collimated northwest-southeast pair. Kristensen et al. (2013) detected only the northwest side of this outflow in CO 6–5, while Girart et al. (2014) detected both sides in CO 3–2 and SiO 8–7. In addition to previously detected tracers, in this work, we also detect both symmetric counterparts in o-H2CO, H13CN: redshifted with respect to source A to the southeast, blueshifted to the northwest (Fig. 2). Blueshifted emission in o-H2CO, H213CO, and H13CN appears to trace the interface of the northwest outflow with the dust bridge. This observed blueshifted emission at the southwestern edge (‘inside bend’) of the bridge could be an indication of the outflow impacting the quiescent bridge and sweeping up material in its wake. Pockets of enhanced H2CO emission may be the result of molecules being sputtered from the mantles of mildly shocked grains associated to one of the outflow lobes. Model calculations show that such effects can be attained by shock velocities of $\sim$10–15 $\mathrm{km}\,\mathrm{s}^{-1}$ (Caselli et al. 1997; Suutarinen et al. 2014). More energetic shocks would produce increased amounts of SiO in the gas phase through erosion of the grain cores, for which simulations indicate that shock velocities should exceed $\sim$25 $\mathrm{km}\,\mathrm{s}^{-1}$ (e.g., Gusdorf et al. 2008a, b). The lack of correlation between the observed spatial extent of gas-phase SiO and H2CO suggests that it is not the same shocks responsible for both molecules; other aspects, such as local density and temperature, must also play a role in determining the morphology of the H2CO emission map. The morphology of SiO emission, particularly on the southeast side of source A, outlines a more collimated, redshifted outflow (jet) structure than the maps of H2CO and H13CN do (Fig. 2).

The blueshifted structure of C34S is seen deeper into the bridge than H2CO and H13CN, closer to the central axis of the bridge. Despite the overlap between the bridge and the C34S, the blueshifted velocities of the latter suggest that it is not part of the quiescent bridge. The redshifted emission component seen in H2CO, H13CN and C34S appears to be at the northernmost section of the redshifted outflow at position angle $\sim$270° (Yeh et al. 2008). In summary, the two outflow-like emission shapes at position angles 335° and 285° (Fig. 1) are interpreted as two unrelated features. In our data, there is insufficient symmetry and alignment to assign them to the rotating structures around source A at position angle (rotation axis) 326–360° (Favre et al. 2014; Oya et al. 2016). We refer to Fig. 1 for an illustration of the misalignment.

The picture drawn in Fig. 1 also addresses the coupling of small scales (60–1500 au), as seen in interferometric observations, with larger scales ($\sim$2000–30 000 au) measured in single-dish observations. Firstly, we confirm the finding by Yeh et al. (2008) that the east-west outflow pair at a few-hundred au scales is consistent with the much larger outflow lobes extending out to thousands of au in the same direction, likely both driven by component A in the IRAS 16293 system. Secondly, the arcminute-scale outflow pair which is redshifted to the northeast and blueshifted to the southwest (Mizuno et al. 1990; Stark et al. 2004) does not appear to have any counterpart at scales below $\sim$1000 au. This points to an outflow driven by IRAS 16293A at some time in the past, but which has been quenched in recent times. Adopting an inclination angle of 65° for the northeast-southwest outflow lobes and a line-of-sight velocity offset of 10 $\mathrm{km}\,\mathrm{s}^{-1}$ (Stark et al. 2004), the absence of a northeast-southwest flow at scales $\lesssim$1000 au translates into a timeframe for the quenching of a few hundred years.

Finally, the strongly blueshifted components apparent in CO 3–2 (Fig. 4) near source B are interpreted using the position-velocity diagram in Fig. 7. We choose a south-north slice in order to cover the C2H emission on either side of source B as well as the 12CO blob south of source B. In the diagram, displacement values below $-3$″ correspond to the regime of the northern edge of the outflow cone emanating westward from source A, where redshifted velocities are observed. All weighted velocity values for the main CO isotopologue within 3″ on either side of source B, however, are significantly blueshifted with respect to the systemic velocity of source B (vertical dotted black line in Fig. 7). This must mean that surface layers on the front side (facing the observer), probed by optically thick CO emission, are moving away from the center of mass of protostellar source B. Already noted in the 2008 publication by Yeh et al. (their position ‘b2’), the origin of this kinematic structure is still under debate. It may be a compact (and therefore young) outflow feature driven by source B (Loinard et al. 2013; Oya et al. 2018), or alternatively, a bow shock feature related to the northwest outflow driven by A (Kristensen et al. 2013). If driven by B, it remains to be explained why the symmetrical counterpart is hidden from sight, even at displacements 2–3″ north of source B, where dust optical depth would not be sufficient to absorb line emission from the background. If driven by source A, the relatively modest velocity of the NW outflow axis, $\lesssim$5 $\mathrm{km}\,\mathrm{s}^{-1}$ with respect to either source (H13CN and SiO panels in Fig. 2), has somehow been translated into much higher velocities for the dense material in the outskirts of source B (up to 10 $\mathrm{km}\,\mathrm{s}^{-1}$ with respect to the protostar, see dark blue contours in Fig. 4). One conceivable scenario is that a dense stream of gas emanating from source A impacts dense material around source B, which deflects part of the gas such that its velocity vector becomes more aligned with the line of sight direction, which leads to higher line-of-sight velocity offsets. It may also be that some other, yet unidentified dynamical process is contributing. In the same region where 12CO is blueshifted, the bulk mass traced by optically thin C17O has velocities consistent with the systemic velocity of source B. In the disk domain, masked in Fig. 7 with a box spanning a width identical to the diameter of the masks in Fig. 2, even C17O is affected by the high optical depth of line and continuum photons. The C2H emission, which also stretches north-south across source B (Figs. 2, 23), shows a velocity trend in Fig. 7 overlapping that of C17O. These two species show a modest velocity gradient across source B,: $V_{\mathrm{lsr}}$ is +1.8 to +2.7 $\mathrm{km}\,\mathrm{s}^{-1}$ north of source B, and +2.6 to +4.1 $\mathrm{km}\,\mathrm{s}^{-1}$ south of B, whichever physical component it traces.

Two pockets of SiO 8–7 emission are seen $\sim$1″ north and south of source B, at position angle 15° east of north (Fig. 21), not overlapping with the C2H filament (Sect. 5.4) at position angle $-15$°. Line of sight velocities of SiO in these pockets are shifted toward both blue and red sides of the systemic velocity, which makes their association to either outflow or inflow motions ambiguous. The much weaker, but optically thin emission from 29SiO (Fig. 22) shows the same morphology, again with mixed blue and red velocities. The morphology and kinematics of SiO 8–7 in our map is consistent with that observed in SiO 7–6 by Oya et al. (2018), who interpret this structure as a signature of a pole-on pair of outflows.

5 Discussion

To highlight the different physical and dynamical components studied in this work, the molecular gas observed between and around the binary protostellar system IRAS 16293 at 60–1500 au scales is divided into three distinct domains: (i) dense and hot ($>$100 K) gas in the disk or disk-like regions around sources A and B; (ii) more tenuous and colder gas residing in the dust bridge between the protostars; and (iii) kinematically active gas within or on the borders of outflow lobes driven by source A. Domain (i) is not examined in this work, but its kinematics and temperature structure was extensively studied, using a different set of molecular tracers, by Oya et al. (2016) for source A and by Oya et al. (2018) for source B. Domains (ii) and (iii) are discussed in Sects. 5.1 and 5.2. A structure seen in C2H, seemingly unrelated to outflow or infall dynamics, is addressed in Sect. 5.4. Finally, the relative evolutionary stage of the two sources in IRAS 16293 is addressed in Sect. 5.5.

5.1 The quiescent bridge

The arc-like bridge structure between protostars A and B is made up of material which has a moderate density (4$\times$104 $\mathrm{cm}^{-3}$– $\sim$3$\times$107 $\mathrm{cm}^{-3}$). Evidence for the density range is provided by the observation that the only molecular species clearly tracing the morphology of the dust arc is C17O, while all tracers of higher density (see Table 2) follow outflow structures and have line-of-sight velocity structures that deviate significantly from that of C17O. This interpretation is supported when our observations are compared with three-dimensional models from J2018, yielding a qualitative match in morphology of the observed bridge structure (Sect. 4.1). To finally distinguish between the effects of gas temperature and gas density, observation would be needed of the distribution of a similarly high density tracer ($n_{\mathrm{crit}}>10^{6}$ $\mathrm{cm}^{-3}$) which at the same time has a low upper level energy (E_{\mathrm{up}}$$< 30 K). For example, the $J$=3–2 or $J$=2–1 transition of H13CN both satisfy these conditions.

We hypothesize that the bridge is a remnant substructure of the circumbinary envelope (Schöier et al. 2002) or a filamentary core, from which both protostellar sources have formed in the past. With its current mass budget and temperature conditions, the IRAS 16293 bridge filament is stable against further gravitational collapse (Sect. 4.2). In addition, the bridge arc is kinematically quiescent, lying at a flat $V_{\mathrm{lsr}}$ within 0.5 $\mathrm{km}\,\mathrm{s}^{-1}$ of the systemic velocities of both protostars (Fig. 2, top left). This straightforward observational fact rules out a scenario in which the bridge arc is one segment of a large, circumbinary disk (or torus) in which both protostars would be embedded. In such a scenario, the bridge gas would have shown a line-of-sight velocity gradient following the kinematic signatures of sources A ($+$3.1 $\mathrm{km}\,\mathrm{s}^{-1}$) and B ($+$2.7 $\mathrm{km}\,\mathrm{s}^{-1}$), that is, a blue-to-red gradient in the southeast-to-northwest direction. This would only comply with the observed flat velocity distribution if the ‘disk’ would rotate entirely in the plane of the sky, which is inconsistent with the line of sight velocity difference of the two protostars. Compared to other, tighter protostellar binary systems such as GG Tau (Dutrey et al. 2014, 2016) and IRS 43 (Brinch et al. 2016), the bridge that we observe between the components of IRAS 16293 lacks a symmetric complement to close a full circumbinary disk or torus, and its velocity structure is inconsistent with disk-like rotation. If there is any velocity gradient across the bridge, it is in the transverse direction rather than along the length of the arc.

Although the 600 au separation between the two components in IRAS 16293 is somewhere in the mid-field between close and wide binaries (see Sect. 1), we conclude that the formation of protostars A and B has probably occurred through turbulent fragmentation. The reason is that the competing scenario, disk fragmentation (Adams et al. 1989), is unlikely due to: (i) the lack of evidence for a remnant of a circumbinary disk, as discussed in detail in this work (see above); and (ii) the stark misalignment between the two disk-like structures, face-on for source B (Jørgensen et al. 2016; Oya et al. 2018), and roughly edge-on for source A (Pineda et al. 2012; Girart et al. 2014).

Like IRAS 16293, the more evolved Class I binary protostellar system IRAS 04191$+$1523, with two components separated by 860 au, also shows a bridge of gas connecting the two protostellar sources (Lee et al. 2017). Based on the velocity structure of the C18O bridge, Lee et al. (2017) conclude that the two sources and the bridge are a substructure of the same natal envelope, which has given rise to protostar formation through turbulent fragmentation. The object L1521F, hosting only a very young, low-luminosity source (Bourke et al. 2006), has also revealed a filamentary arc of gas (HCO+) when studied at high resolution with ALMA (Tokuda et al. 2014). In this case, the filament itself is not detected in dust continuum emission; it connects to a compact millimeter continuum component at its western-most terminus, whereas protostars have not yet formed on its eastern end. The dense core Barnard 5, with at least one embedded protostar, also exhibits filamentary structure within which three additional, possibly prestellar condensations have formed (Pineda et al. 2015). The physical origin of the filaments in both L1521F and in Barnard 5 is unclear, but their long axis dimensions of $\gtrsim$1000 au make a disk origin seem unlikely. In contrast, the Class 0 L1448 IRS3B triple protostellar system hosts a heavily curved connecting filament seen in C18O, and its kinematics are consistent with being a remnant of a 254 au diameter disk that is hypothesized to have spawned all three cores through fragmentation of a disk or the innermost envelope (Tobin et al. 2016a). A connecting bridge was also found in the triple system SR24, with a $\sim$700 au size, by Fernández-López et al. (2017), who suggest that it, too, may be interpreted as a remnant of a disk spiral arm. In each of these case studies, the bridging filaments appear to have a relation to the nursery of the protostellar sources, either in the form of a protostellar disk undergoing fragmentation or a substructure of a larger-scale envelope. IRAS 16293 appears to belong to the latter category.

5.2 Outflow signatures

The SiO gas surrounding IRAS 16293B at projected radii up to 250 au ($J$=8–7, Sect. 4.3 of this work; $J$=7–6, Oya et al. 2018) is of mixed blue- and red-shifted velocity up to $\pm$3 $\mathrm{km}\,\mathrm{s}^{-1}$ from the systemic velocity of source B. The presence of SiO in the gas phase is often used as a tracer of outflow shocks (e.g., Flower et al. 1996; Caselli et al. 1997; Gusdorf et al. 2008a, b). While the overlapping blueshifted and redshifted SiO emission is suggestive of a geometry with approaching and receding outflows superposed in projection (Oya et al. 2018), the relatively low velocities would not be expected to lead to grain erosion. Correcting the observed velocity for projection effects has a marginal effect, since the $V_{\mathrm{lsr}}$ offset is necessarily almost completely parallel with the true velocity vectors in case of a face-on disk launching an outflow along the line of sight. However, in the scenario that the SiO pockets near source B are due to impact of the northwest outflow from source A, the currently observed velocities of 5–10 $\mathrm{km}\,\mathrm{s}^{-1}$ may not lead to sufficiently strong shocks. Before silicate material from the grain cores is released into the gas phase, relative velocities should exceed $\sim$15$\mathrm{km}\,\mathrm{s}^{-1}$ (e.g., Caselli et al. 1997; Gusdorf et al. 2008a). This leaves neither the pole-on outflow from source B nor the impact of the outflow from A onto source B as a fully viable scenario.

The morphology of two molecular lines that are sensitive to densities above $10^{7}$ $\mathrm{cm}^{-3}$, o-H2CO 51,5–41,4 and H13CN 4–3 (Table 2), are compared with a three-dimensional model in Sect. 4.1. In the observed maps (Fig. 2), these high-density tracers do not coincide with the axis of the dust/C17O bridge. The modeled structure (J2018; Sect. 4.1 of this work) does not produce emission in H2CO or H13CN from its bridge domain, consistent with the observations. Instead, the modeled emission is concentrated in the disk domains (and is extremely optically thick for o-H2CO), with an extended hour-glass at positions above the disk plane of source A, qualitatively consistent with enhanced pockets of emission seen in the observed maps southeast and northwest of the edge-on source A. This is caused by self-absorption in the high-density midplane of the disk or disk-like structure, whereas its surface regions are less opaque, but still sufficiently dense and hot to excite the H13CN into $J$=4. The absence of H13CN in the bridge in the models is due to the chosen ice sublimation temperature of 100 K, which is only exceeded in the innermost regions close to the protostars.

The high-density molecular gas observed in emission near the bridge filament is interpreted as being the effect of two outflows at position angles 315° and 270°, impacting patches of higher density in the ambient envelope. The northernmost of the two outflows is blueshifted, and is posited to be material on the southwestern edge of the bridge, swept up by the blueshifted side of the northwest-southeast outflow pair (Sect. 4.3, Fig. 1). The southern, redshifted component at position angle 285°, coincides with the northern border of the redshifted western outflow cone (Yeh et al. 2008; Fig. 4 in this work). The same molecules are not observed at the southern edge of the western outflow; we speculate that this is due to a density contrast in the ambient medium, with higher density at the northern side, toward the bridge filament. H2CO in particular is known to sublimate efficiently from dust grain surfaces not just by thermal desorption, but also by photodesorption such as occurs in classical photodissociation regions (Van der Wiel et al. 2009; Guzmán et al. 2011). Chemical models of outflow-envelope interface regions (Drozdovskaya et al. 2015) also show enhanced H2CO abundance, at early times. It is therefore natural to expect enhanced gas-phase H2CO abundances in directly irradiated walls of outflow cavities in young protostellar systems such as IRAS 16293.

5.3 Peak location for complex organic molecules

Searches for spectral signatures of relatively faint molecular line signals from new, complex (organic) molecules in the IRAS 16293 system often focus on positions close to source B, rather than source A. The primary reason is the simpler, narrower spectral line profiles near source B, owing to the difference in orientation of the disk-like structures surrounding sources A and B, and the multitude of dynamically active outflow components associated with source A (Sect. 4.3, 5.2). With intrinsically narrower, symmetric line profiles, the disambiguation of nearby, partially overlapping lines is more straightforward. Specifically, the emission strength of lines of complex molecules is seen to peak in a region west-southwest of source B at about 0.5″ from its continuum peak (e.g., Baryshev et al. 2015; Coutens et al. 2016; Jørgensen et al. 2016; Lykke et al. 2017; Ligterink et al. 2017; Calcutt et al. 2018). The offset can be explained by the effect of high optical depth of both line and continuum emission closer to the source center. The comparative lack of complex molecule emission on the north and east side of source B may be related to the bridge filament connecting to source B on that side (Sect. 5.1), where physical conditions, such as temperature and dust optical depth, are likely different. In addition, excitation conditions and spectral line shapes on the east side of source B may be affected by the northwest outflow from source A impacting there (Sect. 5.2). Finally, even to the south and west of source B, dynamical effects do play a role in shaping molecular emission, as shown for example in the SiO morphology (Sect. 3, Fig. 2).

5.4 The C2H filament

A narrow filament of material stretching across source B in the south-north direction is visible in C2H, and a similar, but slightly offset patch of high-velocity (V_{\mathrm{lsr}}$$<$$-6 $\mathrm{km}\,\mathrm{s}^{-1}$) 12CO up to 1$\aas@@fstack{\prime\prime}$5 south of source B. Compared with the C2H emission observed in single-dish JCMT observations (Caux et al. 2011), it is evident that the majority of C2H flux is filtered out in our interferometric ALMA observations, much more than for spectral lines of other species. However, regardless of missing flux from a smooth, large-scale component, our ALMA C2H map reveals compact substructures at scales \sim$$10^{2}-$10^{3}$ au. The velocity profiles of 12CO, C17O and C2H along this feature are shown in Fig. 7. C2H has a rather flat velocity structure along the length of the filament, at $V_{\mathrm{lsr}}$ within $\sim$1 $\mathrm{km}\,\mathrm{s}^{-1}$ of that of source B itself, as already apparent from its intensity-weighted velocity map in Fig. 2. The morphology of this feature, combined with its small line of sight velocity gradient, would be more consistent with a symmetric (outflow or infall) motion almost exactly in the plane of the sky, than with an outflow emanating perpendicularly from the face-on disk of source B. Moreover, the observed spatial distribution of C2H is not consistent with tracing one or multiple walls of outflow cavities carved out by any of the outflows in this system, such as has been observed in C2H morphologies in other protostellar system at similar spatial scales, for example IRAS 15398$-$3359 (Bjerkeli et al. 2016) and VLA 1623$-$2417 (Murillo et al. 2018).

The morphology of C2H in IRAS 16293 is also studied in Murillo et al. (2018), who find that it anti-correlates with that of cyclic C3H2. While non-overlapping distributions of these two chemically related species (Gerin et al. 2011) may seem surprising, there are also models and observations that indicate somewhat different physical conditions being needed for C2H and $c$-C3H2 emission lines to occur (e.g., Cuadrado et al. 2015). Lacking an obvious connection to any of the other features connected to source A, the physical, dynamical and chemical origin of the C2H emission remains an open question.

5.5 The relative evolutionary stage of sources A and B

Based on the conclusions drawn in this paper and accompanying PILS papers (J2018; Calcutt et al. 2018), it remains unclear if sources A and B in the IRAS 16293 system are coeval. Their rotation axes are not aligned and a formation scenario from a giant circumbinary disk is therefore unlikely (see Sect. 5.1). Sources A and B differ in various aspects. Here we list those of relevance to the matter of evolutionary stage. ($i$) Source A is known to drive at least two sets of bipolar outflows, and an additional one in the past, whereas source B shows signs of infall, but exhibits no unambiguous proof of outflow activity (Pineda et al. 2012; Loinard et al. 2013; Kristensen et al. 2013; Sect. 4.3 of this work). ($ii$) J2018 find that the observed dust emission of the face-on disk of source B can only be explained with a vertically extended disk or disk-like model (i.e., one that has not yet undergone significant settling), while the vertical structure of the disk of source A is not constrained by their model. ($iii$) Chemical differentiation is seen between both sources particularly when analyzing CN-bearing species. Such species have time dependent formation pathways (Garrod et al. 2017) with some, such as vinyl cyanide (C2H3CN), being formed later during the warm-up phase, as a result of the destruction of ethyl cyanide (C2H5CN). Calcutt et al. (2018) find vinyl cyanide to be at least nine times more abundant toward the B source. However, these authors stress that although this could indicate that source A is less evolved, it could also imply that the warm-up timescales for source A were very short, limiting the amount of vinyl cyanide that could be formed. Collecting the above information elements, ($i$) and ($ii$) are consistent with a view that source B is at an earlier stage of evolution than source A, whereas ($iii$) could be explained either way, appreciating that the temperature history of each source is not necessarily monotonous and smooth. We thus pose that source B is at the same evolutionary stage or less evolved than source A.

If there is a difference in evolutionary stage between sources A and B, outflowing material from source A may have provided a trigger for gravitational collapse to set in, some six hundred au (projected) further down its own natal filament. Alternatively, even if the initial collapse of source B started spontaneously, the northwest outflow from source A could now be feeding gas and dust originating in protostar A onto the disk of source B. Such a scenario would add a new element to the investigation of the different levels of various groups of chemical compounds found in sources A and B (Bisschop et al. 2008; Calcutt et al. 2018).

6 Summary and conclusions

We have presented an analysis of the distribution and kinematics of molecular gas in the protostellar binary system IRAS 16293 at scales of 60–2000 au, using ALMA-PILS observational data and three-dimensional line radiative transfer modeling. The selected tracers are 0.87 mm dust continuum, line transitions CO 3–2, C17O 3–2, H13CN 4–3, o-H2CO 51,5–41,4, H213CO 51,5–41,4, C34S 7–6, SiO 8–7, 29SiO 8–7, and C2H 47/2–35/2, together highlighting different aspect of the morphology and dynamics of the system. The main conclusions are summarized here.

The arc-shaped filament connecting the two protostellar sources is clearly seen in dust continuum and C17O 3–2 emission. The kinematic pattern of C17O indicates that it is quiescent, and consistent with the systemic velocities of both sources. Other molecular tracers considered in this work show a different spatial distribution, different kinematics, or both (Sect. 3). The gas density of the bridge filament is likely to be between 4$\times$104 $\mathrm{cm}^{-3}$ and $\sim$3$\times$107 $\mathrm{cm}^{-3}$ (Sect. 4.1). A three-dimensional radiative transfer model qualitatively matches the observed structure in dust and C17O, but not those in H13CN and H2CO emission (Sect. 4.1). 2. 2.

Using a simplified description of the bridge filament, balancing thermal pressure and gravity, the bridge is inferred to be stable against radial collapse. While a radially stable filament is a prerequisite for it to fragment along the direction of its axis, our interpretation of the observations of the bridge filament in its current state does not necessarily extend to earlier times when it was likely more massive and thus less stable (Sect. 4.2). 3. 3.

Given the physical nature of the bridge filament, we pose that it is a remnant substructure of a filamentary circumbinary envelope that has undergone turbulent fragmentation to form both protostellar sources (Sect. 5.1). The possibility that the bridge is (a remnant of) a circumbinary disk is ruled out (Sect. 5.1). 4. 4.

The western side of the east-west outflow pair is clearly traced in CO 3–2 in our observations, shown to be accelerating with distance from the launching source A. In addition, it is seen to change from a collimated to a more wide-angle structure at a projected distance of $\sim$350 au from source A, possibly related to a density and pressure gradient of one of the more or less cylindrical structures observed in this system. The launching point of the outflow may be connected to a H2O maser shock spot $\sim$120 au west of source A (Sect. 4.3). 5. 5.

Molecular lines of (isotopologues of) H2CO, HCN, SiO, and CS trace outflow motions or the impact of outflows onto more quiescent components. Dynamic structures near source B might be due to nearly pole-on outflow from source B itself, or impact of the northwest outflow from source A onto dense material in source B (Sect. 4.3). In our observations, we find no counterpart for the large scale (thousands of au) northeast-southwest outflow pair, indicating that its launching engine must have been quenched a few hundred years ago (Sect. 4.3). 6. 6.

A striking, previously unknown structure stretches straight across source B at position angle $\sim$15°. It is only seen in C2H and is kinematically flat, which rules out an origin in outflow activity of source B (Sect. 5.4). 7. 7.

Combining evidence from this work and other studies, we suggest that IRAS 16293B is not more evolved than its sister source A (Sect. 5.5). 8. 8.

When used as a target to search for complex (organic) molecules, as is often done, the physical and dynamical complexity of the IRAS 16293 protostellar binary system should be taken into account. While searches for complex molecules generally avoid the broad line regions near source A, even the often-used positions on or near the blue ‘$\times$’ mark in Fig. 2, 0.5″ offset from source B, may exhibit spectral line shapes with evidence of line-of-sight motions, which will manifest differently in different molecules (Sect. 5.3).

By abandoning spherically symmetric model descriptions of the envelopes of IRAS 16293 (e.g., Jørgensen et al. 2002; Schöier et al. 2002; Crimier et al. 2010) and other protostellar systems, J2018 have taken a significant leap by developing a three-dimensional density model that includes an envelope hosting two individual energy sources, A and B, disk-like structures around each, an interbinary bridge filament, and a self-consistently derived temperature structure of the ensemble. While the radiative transfer scheme initially included dust and CO isotopologues, it was further expanded in this work to include two additional molecular species, H2CO and H13CN. However, there are still steps to be made to make such a model more realistic. For example, chemistry, freezing out and sublimation of molecules could be included, beyond the current implementation using gas-phase abundance jumps at certain threshold temperatures. The current model also lacks various dynamic components such as the outflows, and the rotating, disk-like structure around source A, which may need to be expanded to larger radii than the 150 au used in J2018. In addition, future observational work could answer currently open questions about, for example, (i) the evolutionary state of source B, (ii) the debated small-scale multiplicity of source A, or (iii) the connection between the relatively young, small-scale outflow signatures (a few hundred to a few thousands au) observed in interferometric observations such as in this paper, and the older outflows seen at distances of tens of thousands of au from the protostellar sources. Considering the already known level of complexity of the system, and the additional aspects yet to be uncovered, adding all these elements to a (self-consistent) physical model would further increase the number of degrees of freedom, which may result in challenges when attempting to sensibly constrain the model parameters.

The protostellar stages of this and other binary systems should be studied further in order to understand how interactions with binary companions (through direct illumination or even mass transfer through an outflow) affect the evolution of the inner tens to hundreds of au of such systems. From this work and the literature so far, it appears that IRAS 16293 represents a particular stage in the evolution of binary star systems, a phase which many other protostars may go through at one point in their evolution. IRAS 16293 may not be particularly special; it is thus rather appropriate to use this system to study chemical and physical conditions in a representative analog of the young Solar system.

Acknowledgements.

We thank the anonymous referee for thoughtful comments that helped improve the quality of this paper. MHDvdW, SKJ, JKJ, LEK and HC acknowledge support from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement No. 646908) through ERC Consolidator Grant “S4F”, as well as from a Lundbeck Foundation Group Leader Fellowship. PB acknowledges support from the Swedish Research Council through contract 637-2013-472. Research at the Centre for Star and Planet Formation is funded by the Danish National Research Foundation. The postdoctoral grant of AC is funded by the ERC Starting Grant 3DICE (grant agreement 336474). MND acknowledges the financial support of the Center for Space and Habitability (CSH) Fellowship and the IAU Gruber Foundation Fellowship. CF acknowledges support from the Italian Ministry of Education, Universities and Research, through the grant project SIR (RBSI14ZRHR). This paper makes use of the following ALMA data: ADS/JAO.ALMA#2013.1.00278.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSC and ASIAA (Taiwan), and KASI (Republic of Korea), in co-operation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. This research has made use of NASA’s Astrophysics Data System Bibliographic Services. This research has made use of APLpy, an open-source plotting package for Python (Robitaille & Bressert 2012); Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration 2013), the matplotlib plotting library (Hunter 2007), and the CASA data processing package (McMullin et al. 2007).

Appendix A Line contamination

Since the PILS data set, covering 329.15–362.90 GHz, has a line density of about one line every 3 MHz (Jørgensen et al. 2016), it is important to consider the effect of contamination of other species contributing to the observed signal of a selected tracer. The potential contributions from species and transitions other than those listed in Table 2 are visualized in the spectral profiles in Figs. 8–15 below.

In general, while various contaminating contributions are identified in the spectrum extracted toward the position one beam west-southwest of source B, this signal is exclusively due to relatively complex molecular species, which are known to have spatial distributions confined to the regions in the central 1–2″ around sources A and/or B (e.g., Baryshev et al. 2015; Lykke et al. 2017; Ligterink et al. 2017), domains explicitly not under study in the current work. The only exception is an SO2 transition near H13CN 4–3 profiles. This and other cases that require a more detailed description are addressed in the following paragraphs.

C17O.

The hyperfine splitting for the C17O 3–2 rotational transition is only 0.5 MHz. We therefore indicate the rest frequency of this transition with a single vertical dashed line in Fig. 8, derived from the mean frequency of the contributing hyperfine components with Einstein $A$ coefficients >$$10^{-6} s*-1*.

H13CN 4–3.

Of the three separated hyperfine components of this rotational transition (Cazzoli & Puzzarini 2005), the central one ($\Delta F$=1) is taken as the reference for the velocity axis. There is potential contamination of the H13CN 4–3 signal by the SO2 132,12–121,11 transition (345.3385377 GHz), with a modest upper level energy of 92 K, at only 1.2 MHz from the rest frequency of H13CN. At the position 0.5″ southwest of source B (‘$\times$’ mark in Fig. 2), preferred for identifying spectral signatures of rare, complex molecules in previous PILS papers (e.g., Coutens et al. 2016; Lykke et al. 2017; Ligterink et al. 2017), lines are narrow and spectrally resolved, and we do not see evidence of SO2 contribution to the H13CN total intensity beyond a few percent. At other positions, for example in the western outflow from source A or close to source A itself, contamination cannot be ruled out.

29SiO 8–7.

The signal of the 29SiO 8–7 itself is weak and somewhat blueshifted ($\sim$1 $\mathrm{km}\,\mathrm{s}^{-1}$) toward the position 0.5″ offset from source B (Fig. 14). We are confident, however, that the signal ascribed to 29SiO in the channel maps (Fig. 22) and velocity map (Fig. 2), because the morphology of the redshifted 29SiO channel maps follows that of the brightest patches of SiO (Fig. 21) and the spectral line shape toward those positions (for example, the extension southeast of source A) is single-peaked.

C2H 47/2–35/2.

For the C2H $N_{J}$=47/2–35/2 transition, there are two hyperfine structure components ($F$=4–3 and 3–2), separated by 1.4 MHz. Such a small separation means that they lie within 1.2 $\mathrm{km}\,\mathrm{s}^{-1}$ of one another in the spectral dimension, but are resolved in our 0.2 $\mathrm{km}\,\mathrm{s}^{-1}$ spectral resolution data cubes, which complicates the interpretation of its velocity structure. In this work, the reference frequency (zero-point for $V_{\mathrm{lsr}}$) is taken as the mid-point between the two hyperfine structure lines, $F$=4–3 at 349.39928 GHz and $F$=3–2 at 349.40067 GHz (Padovani et al. 2009). In addition, most of the redshifted (¿$+5$ $\mathrm{km}\,\mathrm{s}^{-1}$) emission in the velocity map of C2H (Fig. 2, Fig. 23) is due to the bright CH3CN 193–183 transition at 349.39329 GHz. For this reason, we have adopted a narrower integration for C2H than for the other lines (cf. Fig. 15 vs. Fig. 8–14). In positions near sources A and B, the CH3CN line also contaminates channels closer to the systemic velocity of C2H.

Appendix B Channel map figures

Channel maps of the molecular lines listed in Table 2 are shown in Fig. 16–23.

Appendix C Model maps without freeze-out

In Fig. 24 we show, for completeness, the result of the radiative transfer model (Sect. 4.1) with freeze-out effects discarded for o-H2CO and H13CN. The dashed contours, representing negative intensity values, indicate that optical depth of H13CN 4–3 is high when all H13CN of all model components is kept in the gas phase, as discussed in Sect. 4.1. These ‘no freeze-out’ calculations are included in this work with the aim of disentangling effects of excitation and those of freeze-out and sublimation. They are explicitly not meant to present an alternative physical scenario; freeze-out of molecules onto dust grain surfaces should always be taken into account in conditions representative of star-forming regions.