© ESO 2019

1. Introduction

The formation of stars in the mass interval 1–20 M_⊙ involves accretion disks and fast collimated outflows or jets. Atacama Large Millimeter Array (ALMA) observations at very high angular resolution (≤0′′.1) have recently provided clear examples of disk-jet systems around young stellar objects (YSO) with luminosities from low-mass to late O-type zero-age-main-sequence (ZAMS) stars (Lee et al. 2017; Ginsburg et al. 2018; Sanna et al. 2019). From a theoretical point of view, the role played by turbulence, magnetic field, gravitational instability and radiation feedback in defining the properties of the system (disk size and velocity profile; jet orientation and collimation) is still highly debated (Matsumoto et al. 2017; Tan et al. 2014; Kuiper & Hosokawa 2018). Detailed studies of disk-jet systems are therefore essential to constrain current star-formation models.

Since all star-formation models state that accretion and ejection are intimately related, studying the properties of the outflows from YSOs, which are much more extended and easier to observe than the small accretion disks, can be an effective way of probing the disk-jet systems. Observations reveal different characteristics for solar-type and B-type YSOs: (1) the molecular outflows tend to be less collimated with increasing YSO mass (Beuther & Shepherd 2005), and (2) only a few radio jets are known in B-type YSOs, whereas they are commonly observed towards low-mass protostars (Moscadelli et al. 2016; Sanna et al. 2018). If the larger distances and the more rapid evolution of B-type YSOs could in part explain these differences, we know that additional processes are at work in the formation of the most massive stars. In comparison with low-mass protostars, the formation of B-type YSOs involves more energetic stellar radiation, which can ionize the surrounding gas and exert radiation pressure on it. The combined action of the magnetic field, the thermal pressure of the ionized gas, and the radiation pressure has been modeled by several authors (see, e.g., Peters et al. 2011; Vaidya et al. 2011), and results in a lower collimation of the outflow. Another notable difference is that more massive B-type stars form in richer clusters (Carpenter et al. 1997; Hillenbrand & Hartmann 1998), and dynamical and radiative interactions among cluster members can significantly affect the properties of the disk-jet systems, causing disk fragmentation, jet precession, or even disruption of the system (see, e.g., Farias & Tan 2018).

In this paper we report on recent ALMA observations of the star-forming region (SFR) G16.59–0.05, also known as IRAS 18182−1433. This region has been the target of several millimeter (mm) and centimeter (cm) interferometric (Submillimeter Array, Owens Valley Radio Observatory, Plateau de Bure Interferometer, Very Large Array – VLA) and Very Long Baseline Interferometry (VLBI) observations (Beuther et al. 2006; Zapata et al. 2006; Furuya et al. 2008; Sanna et al. 2010; Moscadelli et al. 2013). The bolometric luminosity of the region is ~10⁴ L_⊙ (Moscadelli et al. 2013) at a distance of 3.6 ± 0.3 kpc, as determined via maser trigonometric parallax observations (Sato et al. 2014). Inside a molecular clump of size ≈0.5 pc and mass ≈1900 M_⊙ (Beuther et al. 2002), the VLA C-Array observations at 3.6, 1.3, and 0.7 cm by Zapata et al. (2006) identified two compact cm sources, named “a” and “b”, separated by ≈2′′ along the SE–NW direction. From previous mm interferometric observations we know that, while source “a” coincides with the peak of the dust emission, indicating that this source is the most embedded, source “b” is found in good positional correspondence with intense, high-density molecular tracers, marking a hot molecular core (HMC).

European VLBI Network (EVN) observations of the 6.7-GHz methanol masers demonstrate that this maser emission is associated with the HMC and traces an elongated structure of ≈2000 au, with the three-dimensional (3D) maser velocity pattern suggesting rotation (Sanna et al. 2010, see their Fig. 7) about a central mass of ≈12 M_⊙ (Moscadelli et al. 2013). The sensitive Jansky Very Large Array (JVLA) A-Array observations at 5, 2.3, and 1.3 cm by Moscadelli et al. (2013, see their Fig. 1) resolve the compact cm source “b” into a radio jet, elongated ≈3′′ along the E–W direction. The spectral index of the jet is negative, indicating nonthermal emission over most of the jet, except the peak close to the maser disk, where thermal free-free emission is observed (Moscadelli et al. 2013, see their Fig. 2). Water masers, monitored with the Very Long Baseline Array (VLBA) by Sanna et al. (2010), are distributed close to the HMC, and appear to trace a wide, fast bow-shock at the head of the western lobe of the radio jet.

The combination of a rotating disk, as suggested by the 6.7-GHz maser proper motions, and a nonthermal radio jet, one of the few observed in high-mass SFRs (see, e.g., Rodríguez-Kamenetzky et al. 2017), makes G16.59−0.05 a promising target to constrain the formation mechanism of B-type stars. For this reason, we conducted new ALMA observations towards G16.59−0.05, as described in Sect. 2. In Sect. 3, we report our observational results presenting a general view of the physical conditions and kinematics of the gas at the center of the mm clump. Section 4 derives the properties of the disk and the YSO, and Sect. 5 provides a more general discussion of the geometrical and physical relations among the YSO, the disk, and the jet. Finally, our conclusions are presented in Sect. 6.

2. ALMA observations

ALMA observed G16.59−0.05 during Cycle 3 in September 2016. Fourty-seven 12-m antennae were used in an extended configuration with baselines ranging from 15 to 3300 m. The observing time on the target source was about 1.8 h. The bandpass calibrator was the strong quasar J1924−2914, and the phases were calibrated from interleaved observations (every 8–9 min) of the quasar J1832−2039, separated from the target on the sky by ≈7°. To check the astrometric accuracy, ALMA also observed the quasar J1830–1606, separated from the target by ≈3°. We have verified that the peak emission of the phase-calibrated image of J1830–1606 differs from the nominal position by less than 10 mas, which provides an estimate of the astrometric uncertainty of our ALMA images.

The correlator frequency setup consisted of six spectral windows (SPWs), one broad 2-GHz spectral unit to obtain a sensitive continuum measurement at ≈242 GHz, and five narrower SPWs to cover a large number of lines, in particular the CH₃ CN (J = 14–13), CH₃OH, and SiO rotational transitions. For each SPW, Table 1 reports the frequency coverage, spectral resolution, and sensitivity, and also indicates the most prominent molecular species.

Data calibration was performed using the pipeline for ALMA data analysis in the Common Astronomy Software Applications (CASA, McMullin et al. 2007) package, version 4.7. For each SPW, looking at the plots of (baseline-averaged) uv-amplitude versus channels, we selected the most intense emission line (always ≳10 times stronger than the continuum level of the uv-amplitude spectra), self-calibrated its channel-averaged emission, and, before imaging, applied the self-calibration phase solutions to all the channels of the SPW. Self-calibration improves the signal to noise ratio (S/N) of the final images of the SPWs by a factor of 1.1–1.5, with the larger gain factors obtained in correspondence with the stronger self-calibrated lines. The dynamical range of the line images varies in the range 20–50, depending on the considered SPW. The images for the continuum and line emissions were produced manually using the CLEAN task, with the robust parameter of Briggs (1995) set to 0.5, as a compromise between resolution and sensitivity to extended emission. The clean beam FWHM of the resulting images varies in the range 0′′. 13–0′′.17. The mm continuum image of G16.59−0.05 has a 1σ rms noise level of 0.08 mJy beam⁻¹, limited by the dynamic range. The 1σ rms noise in a single spectral channel varies in the interval 1–3 mJy beam⁻¹, depending on the considered SPW/channel.

We have employed a specific procedure to determine the continuum level of the spectra and subtract it from the line emission. Since all the observed SPWs present a “line forest”, identifying the channels where no line emission is present is a very difficult task. For each SPW, we use STATCONT¹ (Sanchez-Monge et al. 2018), a statistical method to estimate the continuum level at each position of the map from the spectral distribution of the intensity at that position.

3. Results

This section presents the main results of our ALMA observations, first reporting on the physical properties of the parental core and then focusing on the kinematics of the embedded high-mass YSO, the main target of our study.

3.1. G16.59−0.05 B: the prominent molecular source in a cluster

The ALMA 1.2-mm continuum map presented in Fig. 1 reveals four main sources, distributed from SE to NW over a region of ~10⁴ au. The source B corresponds to the HMC identified in previous interferometric observations at 3 and 1 mm (Furuya et al. 2008; Beuther et al. 2006), named source “b” in Moscadelli et al. (2013). The compact “mm core” from previous, lower-angular-resolution, interferometric observations, named source “a” in Moscadelli et al. (2013), is now resolved into three distinct sources, labeled A1, A2, and A3, separated by 2000–3000 au. Considering emission down to 10% of the peak, the sources B, A2, and A3 are relatively compact with size ≲0′′.5, while A1 has a significantly flatter spatial distribution extending up to ≈1′′.5.

Figure 2 shows the spectra of the CH₃CN and CH₃OH lines prominent in our ALMA frequency setup, extracted at the position of the 1.2-mm peak for each of the four continuum sources. In agreement with previous findings, the molecular emission is significantlystronger towards B (by a factor of 2.5–5). While in sources A1 and A2 the intensity of the CH₃CN and CH₃OH lines is comparable, in A3 only weak CH₃OH emission is detected. Table 2 lists all the molecular lines analyzed in this paper. The intense lines of CH₃CN and CH₃OH are mostly employed to trace the gas kinematics inside B. While the optical depth of these lines close to the peak approaches 1 in some cases, the line wings are always optically thin and allow us to map the gas kinematics around the high-mass YSO in B with good S/N. Using the SLIM (Spectral Line Identification and Modeling) tool of MADCUBA², we surveyedthe spectra of the four mm sources to search for molecular species with a relatively large number of unblended and optically thin lines suitable for deriving the gas physical conditions. We selected the transitions of the CH₃OCH₃ and CH₃OH³ molecules reported in boldface characters in Table 2, which have optical depths of less (or much less) than 0.5, and cover a sufficiently wide range in excitation energy to guarantee a good estimate of the excitation temperature. To derive the physical parameters we used the tool AUTOFIT of MADCUBA, which compares the observed spectra with the LTE synthetic spectra, taking into account all transitions and the line opacities. Leaving four parameters free to vary, that is, column density (N_col), excitation temperature (T_ex), Local Standard of Rest (LSR) velocity (V_LSR), and line width (FWHM), AUTOFIT provides the best nonlinear least-squared fit using the Levenberg-Marquardt algorithm. Figure 3 shows that the emission profiles of the unblended CH₃OCH₃ and CH₃OH transitions are reasonably well fitted with MADCUBA, determining the values of column density, temperature, velocity and line width listed in Table 3. The masses in molecular gas of the mm sources are derived from the mm fluxes (integrating inside the contours at 10% of the peak marked in Fig. 1) and the fitted gas temperatures, assuming a dust opacity of 1 cm² g⁻¹ at 1.4 mm (Ossenkopf & Henning 1994) and a gas-to-dust mass ratio of 100. We used the mm continuum peak flux and the fitted gas temperature to estimate the optical depth of the dust emission inside the mm sources, ranging from 0.25 for source B to 0.45 for source A3. For each mm source, the mass value reported in Table 3 was corrected using the corresponding dust opacity. The flux error is evaluated by considering the contributions of the rms noise of the dust continuum map close to the intense sources, ≈0.4 mJy beam⁻¹, and the ALMA flux calibration uncertainty of 5%, as recently estimated from the analysis of calibrators in bands 3 and 6 by Bonato et al. (2018, see also references therein). The uncertainty in the mass of the mm sources is derived by taking into account both the flux and gas temperature uncertainties.

Fig. 1 ALMA 1.2-mm continuum (grayscale map and white contours). The color scale at the topgives the intensity of the map. The plotted levels are from 10 to 90% of 0.032 Jy beam−1 in steps of 10%. The main continuum sources are indicated with their corresponding labels. The level at 10% of the peak around each source is plotted in magenta. The velocity-integrated intensity of the CH3 OH 50,5 –40,4 line is shown with dashed cyan contours, plotting levels from 20 to 90% of 1.1 Jy beam−1 km s−1 in steps of 10%. The ALMA beams for the continuum and the CH3OH 50,5 –40,4 line are reported on the bottom left and right corners, respectively. The big white star marks the position of source “a” (Moscadelli et al. 2013).

Open with DEXTER

Fig. 2 Spectra of the prominent lines of CH3CN (J = 14–13; left panels) and CH3OH (right panels) extracted at the position of the 1.2-mm peak for each of the four continuum sources B, A1, A2 and A3 (from top to bottom panels, respectively). The spectra are shown as brightness temperature vs. rest frequency. The transitions of the two molecular species (see Table 2) are labeled in the upper panels.

Open with DEXTER

3.2. The massive YSO Bm

We now focus on the gas kinematics around the massive YSO inside the mm source B. In the following we refer to thisYSO as Bm. This is probably the most massive object in source B. A well-defined, SW–NE V_LSR gradient is detected at the position of Bm in all the observed high-density gas tracers. Figures 4–6 present results from low-E_u lines of CH₃OH, CH₃CN, and C³⁴S, respectively. We note that the JVLA 22-GHz continuum, pinpointing the YSO, falls exactly between the SW blue-shifted and the NE red-shifted line emission, and that the direction of the V_LSR gradient forms a large (≈70°) angle with the radio jet traced by the extended JVLA 13-GHz continuum. These findings suggest that the V_LSR gradient could be due to rotation of the gas in an envelope and/or disk surrounding the YSO.

The P–V plots produced along the axis (at PA = 18°) of the V_LSR gradient confirm that we are observing envelope–disk rotation. The plots have a butterfly-like shape (mostly clear in Fig. 6 for the C³⁴S line), with well-defined spurs at high absolute velocities and small offsets in the second and fourth quadrants. These spurs correspond to gas whose line-of-sight velocity increases with decreasing radius and could be consistent with Keplerian rotation. In Sect. 4, we fit the velocity profile and estimate the disk radius and mass, and the YSO mass. Here, we constrain the YSO V_LSR and position by inspecting the P–V plots. The positional offset (along the axis at PA = 18°) of the YSO is delimited by the offsets of the high-velocity spurs in the second and fourth quadrants of the P–V plots, falling inside the interval 0′′–0′′.15 (see Fig. 6). The spurs at low (line of sight) velocities and large offsets are less clear, because they are contaminated by emission of the ambient gas at the systemic velocity. We make use of the P–V plot of the C³⁴ S line (see Fig. 6) to confidently set a range of values for the YSO V_LSR within 60–61 km s⁻¹. These constraints for the YSO position and the V_LSR are used in the following analysis.

Looking at Figs. 4–6 (left panels), in all the three molecular tracers the blue-shifted emission extends over an area significantly smaller than that of the red-shifted emission. Moreover, while the red-shifted emission is distributed approximately symmetrically about the major axis of the envelope–disk, the barycenter of the blue-shifted gas is offset to W. These asymmetries in the velocity distribution of the gas around Bm are further discussed in Sects. 5.1 and 5.2.

3.3. The surroundings of Bm

We now consider a larger region encompassing the jet emerging from Bm. Figures 7 and 8 show the channel maps of two intense CH₃CN lines with quite different excitations, the J_K = 14₂–13₂ with E_u∕k = 121 K and the J_K = 14₇–13₇ with E_u∕k = 442 K. These two lines are representative of the velocity distribution of the large majority of the low- and high-excitation molecular transitions detected towards Bm. At high velocities (channels with V_LSR ≤ 56 km s⁻¹ or ≥63 km s⁻¹) both low- and high-excitation molecular lines trace the compact disk close to the high-mass YSO. At V_LSR ≈ 55 km s⁻¹, a compact source is detected only in low-excitation lines at ≈0′′.6 towards N–NW of the mm peak of source B, in proximity to the weak continuum emission and the cluster of water masers located at the largest distance from the radio jet. In the following, we refer to this compact source as B-NW. At central velocities (56 km s⁻¹ ≤ V_LSR≤ 63 km s⁻¹), low-excitation lines show extended emission emerging from (only) the western side of the jet at blue-shifted velocity, and from an arc-like/linear bridge at blue-/red-shifted velocities connecting the massive YSO Bm with the compact source B-NW. High-excitation lines have a quasi-compact structure also at central velocities, with some notable features: (1) at V_LSR ≈ 59 km s⁻¹, a slightly resolved funnel-like structure is traced on top of the western lobe of the jet near Bm, and (2) at V_LSR ≈ 57–58 km s⁻¹, an emission spot is observed adjacent to the radio jet towards N. This latter feature is referred to as HS in the following.

The radio jet from Bm is one of the few known cases of nonthermal jets from high-mass YSOs. In particular, the eastern lobe of the radio jet is well detected with the JVLA at 6 GHz, and is very weak and/or undetected at frequencies ≥10 GHz (Moscadelli et al. 2013, see their Figs. 1 and 2). Our ALMA observations confirm the presence of a strong asymmetry between the eastern and western lobes of the jet, with the western lobe being much more intense in both continuum and molecular lines. Considering also the distribution of the dust emission (see Fig. 1), predominant towards W-NW of Bm, the simplest explanation for all these facts is an E-W density gradient at the position of Bm. Integrating the ALMA 1.2-mm continuum inside the area of each jet lobe, we estimate that the average gas density in the western lobe is about one order of magnitude larger than in the eastern lobe.

Fig. 3 Observed spectrum and MADCUBA LTE best fit shown in black and red, respectively, in each panel. The spectra are shown as brightness temperature vs. rest frequency. In the upper-left corner and at the bottom of each panel, the molecular species and the corresponding line employed in the fit are indicated (see Table 2). Each raw refers to one of the four mm sources: B, A1, A2 and A3, from top to bottom panels, respectively.

Open with DEXTER

4. The Keplerian disk around the high-mass YSO Bm

In this section, we derive the properties of the YSO Bm and its disk. The angular resolution of the P–V plots shown in Figs. 4–6 varies in the range 0′′.13–0′′.17, reflecting the FWHM of the observing beam at the frequency of the molecular line employed to produce the P–V plot. The emission at high blue- (≤56 km s⁻¹) and red-shifted (≥63 km s⁻¹) velocities is approximatively compact (see Figs. 7 and 8), which allows us to fit a 2D Gaussian profile and determine the peak position at each velocity channel. The positional accuracy is equal to (see, e.g., Reid et al. 1988), where I and σ are the peak intensity and the rms noise, respectively, of a given channel. With typical values of the ratio , the error on the position of the fitted channel peaks is as small as a few mas. For channels where the emission is sufficiently compact, this method increases the (relative) positional accuracy by a large factor, while retaining the essential information of the P–V plots.

In Fig. 9 we show that the spatial distribution of the high-velocity emission of CH₃OH and CH₃CN is elongated along the same PA = 18 ± 3°, which we take as the direction of the major axis of the molecular disk around the high-mass YSO Bm. As already noted in Sect. 3.2 for the SW–NE V_LSR gradient, the red-shifted side of the disk is significantly longer than the blue-shifted side. Along the red-shifted side, ≈0′′.15 or ≈500 au in extent, the gas V_LSR increases monotonically approaching the YSO, ranging from ≈62 km s⁻¹ up to ≈69 km s⁻¹. On the contrary, the blue-shifted side of the disk is traced only at the lowest velocities, 50–55 km s⁻¹; the mildly blue-shifted channels, in the range 55–60 km s⁻¹, show more extended structure, which tends to be elongated to the W of Bm along the radio jet (see Fig. 7).

Comparing the two molecular tracers in Fig. 9, the red-shifted side of the disk appears to be better traced by the CH₃OH lines, whose spatial distribution is flatter and more extended than that of the CH₃CN lines. Table 2 shows that the upper energy level (E_u) of the CH₃OH transitions is on average lower than that of the CH₃CN lines. Since the gas temperature in the central layers of the disk is expected to be relatively low, it is reasonable that the CH₃OH transitions are found to be better tracers of the disk midplane. Accordingly, examining the different CH₃OH transitions (see Table 2), we find that the ones with E_u ≤ 150 K present a significantly flatter and more elongated spatial distribution. Considering only the CH₃OH transitions with 20 ≤ E_u ≤ 150 K, Fig. 10 plots V_LSR versus positions projected along the disk major-axis. While the radial profile of the red-shifted velocities is well reproduced with a Keplerian curve (as described below), the blue-shifted emission concentrates within a small offset range (≤0′′.04) and cannot be adequately fitted. However, the average position and V_LSR of the blue-shifted emission is useful to constrain the position and the V_LSR of the YSO.

The free parameters of the Keplerian fit are the YSO’s positional offset, S_⋆, LSR velocity, V_⋆, and mass, M_⋆ sin²(i). Indicating with i the inclination of the disk rotation axis with the line of sight, the formulation of the latter parameter takes into account that the disk plane is actually seen at an angle 90 − i from the line of sight, and the fitted V_LSR correspond to the actual rotation velocities multiplied by the factor sin(i). We have minimized the following χ² expression: (1)

where V_j and S_j are the channel V_LSR (in km s⁻¹) and corresponding peak positions (inarcsecond), the index j runs over all the fitted channels, and the + and − symbols hold for red- and blue-shifted velocities, respectively. The YSO mass M_⋆ is given in M_⊙. In order to take into account both the uncertainty on the velocity and that on the position, the global velocity error ΔV_j was obtained by summing in quadrature two errors: that on the velocity (taken equal to half of the channel width) and that obtained by converting the error on the offset into a velocity error through the function fitted to the data. In Sect. 3.2, based on the symmetrical patterns of the P–V plots, we constrained the ranges for the YSO positional offset and V_LSR to: 0′′ ≤ S_⋆ ≤ 0′′. 15, 60.0 km s⁻¹ ≤ V_⋆ ≤ 61.0 km s⁻¹. The mass of the YSO is searched over the range 5–15 M_⊙, consistent with the estimated bolometric luminosity of 10⁴ L_⊙. Figure 11 reports the plots of the distribution of the χ² as a function of the free parameters. The white contour in these two plots draws the 1-σ confidence level for the three free parameters, following Lampton et al. (1976). These plots show that we do find an absolute minimum of the χ², and the determined 1-σ intervals for the parameter values are: S_⋆ = −0′′.084 ± 0′′.020, V_⋆ = 60.4 ± 0.5 km s⁻¹, and M_⋆ sin²(i) = 10 ± 2 M_⊙.

Following the analysis by Moscadelli et al. (2013, see in particular their Fig. 5), source B, the HMC, dominates the IR emission of the region, suggesting that our estimate of 10⁴ L_⊙ for the bolometric luminosity of B is reasonable and not too high. This luminosity corresponds to a single ZAMS star of ≈13 M_⊙ (see Davies et al. 2011), which is an upper limit for the mass of Bm. Comparing this upper limit with the YSO mass from the Keplerian fit, the value of the disk inclination angle is constrained within the interval 60° ≤ i ≤ 120°, that is, the plane of the disk is within 30° of the line of sight. However, if a sizeable fraction of the bolometric luminosity is due to gravitational energy released in the accretion process rather than nuclear burning, the mass of Bm could be less than 13 M_⊙, and the disk would be almost edge-on. Inside the disk radius of ≈0′′.15, the 1.2-mmcontinuum flux is 42 mJy, corresponding to ≈1 M_⊙ (making thesame assumptions as in Sect. 3.1). Since this value is much less than the YSO mass ≈10 M_⊙, our choice of fitting a Keplerian velocity profile appears to be well justified a posteriori.

Figure 12 (lower panel) compares the kinematics of the disk, traced with thermal CH₃OH, with that of the 6.7-GHz CH₃OH masers observed with VLBI by Sanna et al. (2010). The 6.7-GHz masers trace a slightly elongated structure oriented along a PA (approximately −40°) quite different from that (18°) of the YSO disk. The maser bipolar V_LSR distribution (red-(blue-)shifted to NW(SE)) and the proper motion pattern were interpreted by Sanna et al. (2010) in terms of rotation, and the 6.7-GHz masers were the first indication of the existence of a rotating disk around Bm. However, the new ALMA data indicate now that the 6.7-GHz masers are not tracing the disk midplane, as originally assumed. While the positions and V_LSR of the NW red-shifted 6.7-GHz cluster are in reasonable agreement with the red-shifted side of the YSO disk, the SE blue-shifted cluster is placed 0′′.1–0′′.2 SE of the disk midplane. Looking at Fig. 12, the two groups of maser features with measured proper motions are projected on the sky close to the disk rotation axis (i.e., along the perpendicular to the disk major axis crossing Bm, at the peak of the compact 22-GHz continuum) and the proper motions are mainly directed perpendicular to this axis. These findings still suggest that the 6.7-GHz masers are mainly tracing rotation about Bm, although they may be partly offset from the disk midplane.

Fig. 4 Gas kinematics towards the high-mass YSO Bm. Left panel: black contours and the color map show, respectively, the velocity-integrated intensity and the intensity-averaged velocity of the CH3 OH 50,5 –40,4 line (Eu = 35 K). Plotted contours are from 10 to 90% of 1.1 Jy beam−1 km s−1 in steps of 10% and the scale on the right gives the color-velocity conversion. The grayscale-filled and white contours show the JVLA 22-GHz and 13-GHz continuum emission, respectively. Contours are from 50 to 90% of 0.28 mJy beam−1 in steps of 10% at 22 GHz, and from 20 to 90% of 0.1 mJy beam−1 in steps of 10% at 13 GHz. The JVLA beams at 22 and 13 GHz are reported in the bottom-left and bottom-right corners, respectively. Right panel: P–V plot of the CH3 OH 50,5 –40,4 line. The cut(at PA = 18°) along which positions are evaluated is indicated with the dashed black line in the left panel. To produce the P–V plot, we averaged the emission inside a strip parallel to the cut and 0.′′ 4 in width to include only the disk–envelope around the YSO. The intensity scale is shown on the right. The spatial and velocity resolutions are indicated with horizontal white and vertical cyan error bars, respectively, in the bottom-left corner. The horizontal and vertical white continuous lines mark the YSO VLSR and positional offset, respectively, estimated with a Keplerian fit (see Sect. 4). The yellow dashed vertical line indicates the positional offset of the peak of the 1.2-mm continuum emission of source B. The four quadrants are labeled.

Open with DEXTER

	Fig. 5 As in Fig. 4 but for the CH3CN JK = 142–132 line (Eu = 121 K). Left panel: black contours show the velocity-integrated intensity with levels ranging from 10 to 90% of 0.9 Jy beam−1 km s−1 in steps of 10%.
Open with DEXTER

Fig. 6 As in Fig. 4 but for the C34S J = 5–4 line (Eu = 28 K). Left panel: black contours show the velocity-integrated intensity with levels ranging from 10 to 90% of 0.8 Jy beam−1 km s−1 in steps of 10%. Right panel: horizontal and vertical white continuous lines, and the blue curves mark, respectively, the VLSR and positional offset of the YSO, and the Keplerian profile around a YSO of 10 M⊙. The two horizontal and vertical, white dashed lines indicate the maximum interval of variation for the YSO VLSR and positional offset, respectively, estimated by eye whilst looking at the symmetrical patterns of the P–V plot. We note that the plotted offset range is larger than for the P–V plots of Figs. 4 and 5.

Open with DEXTER

Fig. 7 Color map of the intensity of the CH3CN JK = 142–132 line (Eu = 121 K) at a differentVLSR (km s−1), indicated in the upper-left corner of each panel. The color-intensity scale is shown on the right. The grayscale-filled and white contours show the JVLA 22-GHz and 13-GHz continuum, respectively, plotting the same levels as in Fig. 4. The green contours reproduce the ALMA 1.2-mm continuum, showing levels at 10, 20, 40, and 80% of the peak emission of 0.032 Jy beam−1. The magenta-edged white triangles mark the positions of the 22-GHz water masers derived through VLBI observations by Sanna et al. (2010). Upper first and second panels: high-mass YSO Bm, and the compact molecular source B-NW are labeled.

Open with DEXTER

	Fig. 8 As in Fig. 7 but for the CH3CN JK = 147–137 line (Eu = 442 K). Upper last and middle first panels: emission spot HS, adjacent to the radio jet, is labeled.
Open with DEXTER

Fig. 9 The disk around the high-mass YSO Bm. Colored dots indicate the peak positions of the most blue- and red-shifted velocity channels for the emission of the nine CH3 OH (left panel) and seven CH3CN (right panel) lines listed in Table 2. Colors represent VLSR as coded on the right of each panel. The dashed black lines, at PA = 18° and 17° for CH3 OH and CH3CN, respectively, show the linear fits to the spatial distribution of the channel peaks. The grayscale-filled and black contours represent the JVLA 22-GHz and 13-GHz continuum, respectively, with the same levels as in Fig. 4.

Open with DEXTER

	Fig. 10 Red and black error bars: major-axis projected positions and VLSR (together with the corresponding errors), respectively, for the highest-velocity emission peaks of the six CH3 OH lines with 20 K ≤ Eu ≤ 150 K listed in Table 2. The blue curve is the best Keplerian fit to the data. The horizontal and vertical dashed lines indicate the fitted YSO VLSR and position, respectively.
Open with DEXTER

Fig. 11 χ2 distribution from the Keplerian fit (see Fig. 10) as a function of the free parameters YSO mass and position (left panel), and YSO VLSR and mass (right panel). To produce these plots, the third free parameter of the Keplerian fit, that is, the YSO VLSR and position in the left and right plot, respectively, is taken equal to the best-fit value. The color scale at the top of each panel gives the value of the χ2 in a logarithmic scale. In each of the two plots, the position of the minimum χ2, = 6.6, is labeled, and the white contour indicates the level +3.5 = 10.1, which corresponds to the 1-σ confidence level for three free parameters (Lampton et al. 1976).

Open with DEXTER

5. The YSO, the disk, and the jet

In this section, we discuss in more detail the interactions among the YSO, the disk, and the jet.

5.1. Jet structure and orientation

Looking at Fig. 9, one can note that, near Bm, the slightly resolved 13-GHz continuum is elongated to NW along a direction roughly perpendicular to the disk. As shown in Figs. 7 and 8, the axis of the funnel-like structure on top of the radio jet traced by the CH₃CN emission at V_LSR ≈ 59 km s⁻¹ also has a similar NW orientation. This structure can be interpreted in terms of the jet cavity, and together with the shape of the 13-GHz emission near Bm, could mark the direction of ejection of the jet close to the YSO. At larger separation from the YSO, the 13-GHz continuum seems to bend and approach the E-W direction, which could suggest that precession or recollimation is taking place.

The nonthermal continuum emission in the two lobes of the radio jet, shown in the upper panel of Fig. 12, emerges from shocks responsible for the acceleration of the electrons to relativistic velocities. These shocks originate near Bm where the jet impinges on high-density material and propagate away along the jet direction at a speed, assuming momentum conservation, of V_sh = V_jet, where n_jet and V_jet are the jet density and ejection velocity, respectively, and n_amb is the ambient density (see, e.g., Masson & Chernin 1993). Looking at the upper panel of Fig. 12, one can note that the western lobe is significantly closer to Bm, which could indicate that the nonthermal shock has moved more slowly towards W because of the higher ambient density, as already discussed in Sect. 3.3. From the JVLA 6-GHz map we derive that the eastern lobe is about three times more distant from Bm than the western lobe, which, assuming the jet is ejected symmetrically, agrees well with our previous estimate of a density contrast in the ambient gas by a factor of approximately ten between the western and eastern side. We denote with ΔS^W and ΔS^E the separations (projected onto the line connecting the lobes – see Fig. 12) between the position of Bm and the emission peaks of the western and eastern lobes, respectively, and with and the shock propagation velocities towards W and E, respectively. The time elapsed since the episode of ejection responsible for the present emission in the jet lobes is then ΔT = ΔS^W/ = ΔS^E/. From the 6-GHz map we obtain ΔS^W ≈ 0′′.33 or ≈1200 au, and ΔT ≈ 6000/ yr, where is expressed in km s⁻¹.

Fig. 12 Lower panel, colored dots: positions and VLSR of the thermal CH3OH emission as in Fig. 9, with the black dashed line marking the major axis of the disk around the high-mass YSO Bm. The grayscale map and the white contours reproduce the ALMA 1.2-mm continuum. The intensity scale for the map is shown on the top, and the plotted contours are the same as in Fig. 1. The JVLA 22-GHz continuum is shown with cyan contours (same levels as in Fig. 4). Colored squares denote positions and VLSR of the 6.7-GHz CH3OH masers (Sanna et al. 2010). The maser proper motions are represented by black arrows, with the amplitude scale given in the bottom-right corner. The black dotted line marks the major axis of the distribution of the CH3 OH masers. The dashed magenta line has the same meaning as in the upper panel. Upper panel, black contours: JVLA 6-GHz continuum with levels from 20 to 90%, in steps of 10%, and 95% of 0.17 mJy beam−1. The dashed magenta line connects the peak emission of the eastern and western lobes of the radio jet. The grayscale map, and the white and cyan contours have the same meaning as in the lower panel. VLBI positions and proper motions of the water masers by Sanna et al. (2010) are indicated with magenta triangles and black arrows, respectively. The amplitude scale for the proper motions is given in the bottom-left corner. The beam of the JVLA 6 GHz continuum image is reported in the bottom-right corner.

Open with DEXTER

5.2. Proper motion and precession

Comparing the upper to the lower panel of Fig. 12, it is clear that the line connecting the 6-GHz lobes of the radio jet crosses source B near the 1.2-mm peak but is displaced from Bm (pinpointed by the 22-GHz emission) by ≳100 mas towards N. This offset is significantly larger than the expected positional error of a few tens of milliarcseconds between the 6- and 22-GHz JVLA images⁴ (Moscadelli et al. 2013, 2016). In Sect. 3.2 we noted the clear asymmetry between the NE red-shifted and SW blue-shifted sides of the YSO envelope–disk, with the NE side being much more extended than the SE one (see Figs. 4 and 12). We argue now that both the displacement of Bm from the jet axis and the disk asymmetry could result from the relative motion of Bm with respect to the disk. We assume that Bm has moved from its original position, at the disk center, towards S-SW, along the disk plane.

What is causing the motion of Bm? One possibility is that Bm is orbiting in a gravitationally bound, multiple stellar system. This would naturally explain its motion relative to the disk. In this case, the observed disk would rotate around the multiplesystem and mediate the accretion onto it. Recent observations of bound systems of low-mass YSOs surrounded by a disk show that the regular distribution in mass and velocity of the disk is truncated at an inner radius comparable to the separation among the most massive stellar members (Tobin et al. 2016; Artur de la Villarmois et al. 2018). At radii comparable to these stellar orbits, the disk fragments and gaps, rings, and spirals appear, like those observed in protoplanetary disks. Looking at Fig. 9, one can note that the disk around Bm presents a regular (Keplerian) velocity pattern until reaching the position of the YSO, marked by the compact 22-GHz continuum. The mass responsible for the observed velocity pattern appears to be confined inside a radius ≤200 au from Bm, which implies that other (putative) YSOs embedded in the parental mm source B are not dynamically relevant for the disk or the motion of Bm.

An alternative explanation for the motion of Bm is that the YSO has been ejected from its original position at the disk center by dynamical interaction with one or more companions. The YSO moved through the disk plane until it reached the current position, traced by the 22-GHz continuum. The distribution of the gas in the disk remained essentially unaffected by the movement of the YSO. The fact that we see Keplerian rotation around Bm implies that the material of the disk has been rotating fast enough to re-adjust its velocity field to the new position of the YSO. This sets a lower limit to the crossing time, T_cr, needed for the star to go from the original position to the current position. For a stellar mass of 10 M_⊙ and a radius of the Keplerian disk of 0′′.15 (540 au), we derive a rotational period T_rot = 4000 yr. Thus, the condition T_cr ≥ T_rot = 4000 yr must be satisfied. Using this estimate of T_cr and the distance of ~0′′.1 traveled by Bm, the average speed of Bm must be ≤1 km s⁻¹. This is in agreement with the result that the LSR velocity of Bm of V_⋆ = 60.4 ± 0.5 km s⁻¹ determined from the Keplerian fit (see Sect. 4) corresponds to the systemic velocity of the molecular cloud, V_sys = 60 ± 0.04 km s⁻¹, obtained from single-dish observations (see, e.g., Hill et al. 2010).

For a plausible mass accretion rate of ~10⁻³ M_⊙ yr⁻¹, T_cr is sufficiently long that most of the mass of Bm could be accreted after its ejection from the original position. In this scenario, the structure of the disk has not been significantly affected by the motion of the YSO through it, because the star has become sufficiently massive to dominate the gravitational field of the region, only when it was already close to the current position.

We still have to explain why the jet is offset from the current position of the YSO. One possibility is that the jet was ejected before the YSO left its original position (marked by the intersection between the jet axis and the disk). Following the analysis of Sect. 5.1 and assuming ΔT ~ T_cr ≥ 4000 yr, we derive that the velocities of the shocks must be ≤ 1 km s⁻¹. Such a low velocity, much less than the typical jet speeds ≥ 100 km s⁻¹ (see, e.g., Masqué et al. 2015), implies a very high density ratio (~ 10⁴) between the ambient gas and the jet. Another possibility is that the jet is ejected from the present location of Bm, but is then recollimated on a larger scale along the axis of the magnetic field, which, on that scale, has not been affected by the movement of the YSO. This latter interpretation is supported by the morphology of the 13-GHz continuum and the funnel-like molecular emission observed near Bm, which, as noted in Sect. 5.1, suggests that the jet emerges from Bm perpendicular to the Keplerian disk.

Looking at Fig. 12, it is also interesting to note that the 1.2-mm peak is slightly offset from the jet axis in the direction of the 22-GHz continuum. This offset, 0′′.03–0′′.04, is comparable with the error in the relative astrometry between the JVLA and ALMA images, but it could have a physical explanation, too. In fact, it could result from a combination of two components of dust emission with similar intensities: one characterized by high temperature and low column density at the position of the 22-GHz continuum, the other with low temperature and large column density at the disk center. In order to resolve these two putative components, we have redone the 1.2-mm continuum image uniformly weighting the visibilities, thus increasing the angular resolution from 0′′.15 to 0′′.12, but the image still presents a single peak. We conclude that, if two components are present, their separation must be ≤0′′.1.

6. Conclusions

The main results of our Cycle 3 high-angular resolution (beam FWHM ≈ 0′′.15) ALMA observations towards the high-mass SFR G16.59−0.05 can be resumed as follows. At the center of the main clump, the dust emission is resolved into four small (size ~1′′) mm sources, among which the one (source B) harboring the high-mass YSO (Bm) is the most prominent in molecular emission. Fitting unblended, optically thin lines of CH₃OCH₃ and CH₃OH, we determined temperatures for the gas inside the mm sources in the range 42–131 K, and estimated masses of between 1 and 5 M_⊙.

A well-defined V_LSR gradient (prominent in CH₃OH, CH₃CN and C³⁴ S) is traced in many high-density molecular tracers at the position of the high-mass YSO Bm, and is oriented along a direction forming a large (≈70°) angle with the radio jet previously revealed through sensitive JVLA observations. The P–V plots of this gradient present butterfly-like shapes, and the emission peaks of the molecular lines at high velocity draw linear patterns, indicating that we are observing rotation of the disk-envelope surrounding Bm. The disk radius is estimated to be ≈500 au, and the V_LSR radial distribution is well reproduced by Keplerian rotation around a central mass of 10 ± 2 M_⊙.

The position of Bm, pinpointed by the compact 22-GHz emission, is found to be offset by ≳0′′.1 S of the jet axis and the mm peak of source B. We explain these findings assuming that, following a multiple stellar event, Bm was ejected from the center of the parental mm source and moved through the disk plane until it reached the current position. While the material of the disk had enough time to re-adjust its velocity field to the new position of the YSO, the distribution of matter and magnetic field, responsible for collimating the jet, remained basically unaffected by the movement of the YSO.

Acknowledgements

V.M.R. is funded by the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 664931. We thank the referee Todd Hunter for useful comments, which improved the paper. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2015.1.00408.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

References

All Tables

All Figures

	Fig. 1 ALMA 1.2-mm continuum (grayscale map and white contours). The color scale at the topgives the intensity of the map. The plotted levels are from 10 to 90% of 0.032 Jy beam−1 in steps of 10%. The main continuum sources are indicated with their corresponding labels. The level at 10% of the peak around each source is plotted in magenta. The velocity-integrated intensity of the CH3 OH 50,5 –40,4 line is shown with dashed cyan contours, plotting levels from 20 to 90% of 1.1 Jy beam−1 km s−1 in steps of 10%. The ALMA beams for the continuum and the CH3OH 50,5 –40,4 line are reported on the bottom left and right corners, respectively. The big white star marks the position of source “a” (Moscadelli et al. 2013).
Open with DEXTER
In the text

	Fig. 2 Spectra of the prominent lines of CH3CN (J = 14–13; left panels) and CH3OH (right panels) extracted at the position of the 1.2-mm peak for each of the four continuum sources B, A1, A2 and A3 (from top to bottom panels, respectively). The spectra are shown as brightness temperature vs. rest frequency. The transitions of the two molecular species (see Table 2) are labeled in the upper panels.
Open with DEXTER
In the text

	Fig. 3 Observed spectrum and MADCUBA LTE best fit shown in black and red, respectively, in each panel. The spectra are shown as brightness temperature vs. rest frequency. In the upper-left corner and at the bottom of each panel, the molecular species and the corresponding line employed in the fit are indicated (see Table 2). Each raw refers to one of the four mm sources: B, A1, A2 and A3, from top to bottom panels, respectively.
Open with DEXTER
In the text

	Fig. 4 Gas kinematics towards the high-mass YSO Bm. Left panel: black contours and the color map show, respectively, the velocity-integrated intensity and the intensity-averaged velocity of the CH3 OH 50,5 –40,4 line (Eu = 35 K). Plotted contours are from 10 to 90% of 1.1 Jy beam−1 km s−1 in steps of 10% and the scale on the right gives the color-velocity conversion. The grayscale-filled and white contours show the JVLA 22-GHz and 13-GHz continuum emission, respectively. Contours are from 50 to 90% of 0.28 mJy beam−1 in steps of 10% at 22 GHz, and from 20 to 90% of 0.1 mJy beam−1 in steps of 10% at 13 GHz. The JVLA beams at 22 and 13 GHz are reported in the bottom-left and bottom-right corners, respectively. Right panel: P–V plot of the CH3 OH 50,5 –40,4 line. The cut(at PA = 18°) along which positions are evaluated is indicated with the dashed black line in the left panel. To produce the P–V plot, we averaged the emission inside a strip parallel to the cut and 0.′′ 4 in width to include only the disk–envelope around the YSO. The intensity scale is shown on the right. The spatial and velocity resolutions are indicated with horizontal white and vertical cyan error bars, respectively, in the bottom-left corner. The horizontal and vertical white continuous lines mark the YSO VLSR and positional offset, respectively, estimated with a Keplerian fit (see Sect. 4). The yellow dashed vertical line indicates the positional offset of the peak of the 1.2-mm continuum emission of source B. The four quadrants are labeled.
Open with DEXTER
In the text

	Fig. 5 As in Fig. 4 but for the CH3CN JK = 142–132 line (Eu = 121 K). Left panel: black contours show the velocity-integrated intensity with levels ranging from 10 to 90% of 0.9 Jy beam−1 km s−1 in steps of 10%.
Open with DEXTER
In the text

	Fig. 6 As in Fig. 4 but for the C34S J = 5–4 line (Eu = 28 K). Left panel: black contours show the velocity-integrated intensity with levels ranging from 10 to 90% of 0.8 Jy beam−1 km s−1 in steps of 10%. Right panel: horizontal and vertical white continuous lines, and the blue curves mark, respectively, the VLSR and positional offset of the YSO, and the Keplerian profile around a YSO of 10 M⊙. The two horizontal and vertical, white dashed lines indicate the maximum interval of variation for the YSO VLSR and positional offset, respectively, estimated by eye whilst looking at the symmetrical patterns of the P–V plot. We note that the plotted offset range is larger than for the P–V plots of Figs. 4 and 5.
Open with DEXTER
In the text

	Fig. 7 Color map of the intensity of the CH3CN JK = 142–132 line (Eu = 121 K) at a differentVLSR (km s−1), indicated in the upper-left corner of each panel. The color-intensity scale is shown on the right. The grayscale-filled and white contours show the JVLA 22-GHz and 13-GHz continuum, respectively, plotting the same levels as in Fig. 4. The green contours reproduce the ALMA 1.2-mm continuum, showing levels at 10, 20, 40, and 80% of the peak emission of 0.032 Jy beam−1. The magenta-edged white triangles mark the positions of the 22-GHz water masers derived through VLBI observations by Sanna et al. (2010). Upper first and second panels: high-mass YSO Bm, and the compact molecular source B-NW are labeled.
Open with DEXTER
In the text

	Fig. 8 As in Fig. 7 but for the CH3CN JK = 147–137 line (Eu = 442 K). Upper last and middle first panels: emission spot HS, adjacent to the radio jet, is labeled.
Open with DEXTER
In the text

	Fig. 9 The disk around the high-mass YSO Bm. Colored dots indicate the peak positions of the most blue- and red-shifted velocity channels for the emission of the nine CH3 OH (left panel) and seven CH3CN (right panel) lines listed in Table 2. Colors represent VLSR as coded on the right of each panel. The dashed black lines, at PA = 18° and 17° for CH3 OH and CH3CN, respectively, show the linear fits to the spatial distribution of the channel peaks. The grayscale-filled and black contours represent the JVLA 22-GHz and 13-GHz continuum, respectively, with the same levels as in Fig. 4.
Open with DEXTER
In the text

	Fig. 10 Red and black error bars: major-axis projected positions and VLSR (together with the corresponding errors), respectively, for the highest-velocity emission peaks of the six CH3 OH lines with 20 K ≤ Eu ≤ 150 K listed in Table 2. The blue curve is the best Keplerian fit to the data. The horizontal and vertical dashed lines indicate the fitted YSO VLSR and position, respectively.
Open with DEXTER
In the text

	Fig. 11 χ2 distribution from the Keplerian fit (see Fig. 10) as a function of the free parameters YSO mass and position (left panel), and YSO VLSR and mass (right panel). To produce these plots, the third free parameter of the Keplerian fit, that is, the YSO VLSR and position in the left and right plot, respectively, is taken equal to the best-fit value. The color scale at the top of each panel gives the value of the χ2 in a logarithmic scale. In each of the two plots, the position of the minimum χ2, = 6.6, is labeled, and the white contour indicates the level +3.5 = 10.1, which corresponds to the 1-σ confidence level for three free parameters (Lampton et al. 1976).
Open with DEXTER
In the text

	Fig. 12 Lower panel, colored dots: positions and VLSR of the thermal CH3OH emission as in Fig. 9, with the black dashed line marking the major axis of the disk around the high-mass YSO Bm. The grayscale map and the white contours reproduce the ALMA 1.2-mm continuum. The intensity scale for the map is shown on the top, and the plotted contours are the same as in Fig. 1. The JVLA 22-GHz continuum is shown with cyan contours (same levels as in Fig. 4). Colored squares denote positions and VLSR of the 6.7-GHz CH3OH masers (Sanna et al. 2010). The maser proper motions are represented by black arrows, with the amplitude scale given in the bottom-right corner. The black dotted line marks the major axis of the distribution of the CH3 OH masers. The dashed magenta line has the same meaning as in the upper panel. Upper panel, black contours: JVLA 6-GHz continuum with levels from 20 to 90%, in steps of 10%, and 95% of 0.17 mJy beam−1. The dashed magenta line connects the peak emission of the eastern and western lobes of the radio jet. The grayscale map, and the white and cyan contours have the same meaning as in the lower panel. VLBI positions and proper motions of the water masers by Sanna et al. (2010) are indicated with magenta triangles and black arrows, respectively. The amplitude scale for the proper motions is given in the bottom-left corner. The beam of the JVLA 6 GHz continuum image is reported in the bottom-right corner.
Open with DEXTER
In the text