Investigating the structure and fragmentation of a highly filamentary IRDC Free

Abstract

INTRODUCTION

Although a basic mechanism for the specific case of isolated low-mass star formation has been investigated over several decades (e.g. Shu, Adams & Lizano 1987), a more generalized model, one that incorporates a consistent description for the formation of high-mass (>8 M_⊙) stars, is still lacking. An important step in developing a holistic understanding of the star formation process is identifying and categorising the initial phases of high-mass star formation. Ultimately this requires detailed knowledge of their host molecular clouds.

Discovered as silhouettes against the bright Galactic mid-infrared background, infrared dark clouds (hereafter, IRDCs) were quickly identified as having the potential to aid our understanding of the star formation process (Pérault et al. 1996; Egan et al. 1998). Initial studies set about categorizing their physical properties, finding broad ranges in size (∼1–10 pc, with rare examples exceeding 50 pc e.g. ‘Nessie’; Jackson et al. 2010), mass (10²–10⁵ M_⊙), and column density (∼10²²–10²⁵ cm⁻²) (e.g. Carey et al. 1998; Egan et al. 1998; Rathborne, Jackson & Simon 2006; Simon et al. 2006). Subsequent investigations categorizing their temperatures (≲ 25 K; Pillai et al. 2006; Peretto et al. 2010; Ragan, Bergin & Wilner 2011; Fontani et al. 2012; Chira et al. 2013), chemistry (e.g. Sakai et al. 2008; Gibson et al. 2009; Jiménez-Serra et al. 2010; Vasyunina et al. 2011; Sanhueza et al. 2012, 2013; Pon et al. 2015; Lackington et al. 2016), and kinematics (e.g. Devine et al. 2011; Henshaw et al. 2013, 2014; Peretto et al. 2013, 2014; Tackenberg et al. 2014; Dirienzo et al. 2015; Schneider et al. 2015; Pon et al. 2016) have ensued. Although the broad range in characteristics dictates that not all IRDCs will form high-mass stars (Kauffmann & Pillai 2010), identifying massive and relatively quiescent molecular clouds (those that are yet to be globally influenced by feedback effects from massive young stellar objects) is a crucially important step in understanding the initial conditions for high-mass star formation.

The focus of this paper is G035.39–00.33, a massive (∼2 × 10⁴ M_⊙; see table 1 of Kainulainen & Tan 2013) and filamentary IRDC, with a kinematic distance of 2.9 kpc (Simon et al. 2006). Originally part of a sample selected from the Rathborne et al. (2006) study by Butler & Tan (2009), G035.39–00.33 was selected for further investigation due to its high-contrast against the Galactic mid-infrared background, and because it exhibits extended regions with no obvious tracers of star formation activity (4.5, 8, and 24 μm emission; Carey et al. 2009; Chambers et al. 2009). Since 2010, G035.39–00.33 has been the focal point of a dedicated research effort whose aim is to provide a detailed case study of the physical structure, chemistry, and dynamics of a single IRDC.

The first results of this case study were presented by Jiménez-Serra et al. (2010, Paper I), who identified faint, but widespread, SiO emission throughout G035.39–00.33. Although a currently undetected population of low-mass protostars may account for this emission, Jiménez-Serra et al. (2010) discussed the possibility that such a signature may represent a ‘fossil record’ of either the cloud formation process or a cloud merger.

The potential to use G035.39–00.33 as a laboratory for studying the early phases of the star formation process is supported by both Hernandez et al. (2011, Paper II) and Barnes et al. (2016, Paper VII), who report widespread depletion of CO and widespread emission from deuterated species (in this case, N₂D⁺), respectively. These two observations emphasize the presence of cold (<20 K) and dense gas, yet to be globally affected by stellar feedback, where CO molecules have frozen on to the surface of dust grains leading to an enhancement in the abundance of deuterated nitrogen-bearing molecules. Comparing the observed abundance of deuterated species with that predicted by chemical models (Kong et al. 2015), Barnes et al. (2016) estimate the age of the cloud to be ∼3 Myr old. This may imply, therefore, that although star formation within the cloud remains within an early evolutionary phase, the cloud itself is dynamically old. Having existed for 5–10 local free-fall times, G035.39–00.33 may have had sufficient time to settle into a state of near virial equilibrium, as concluded by Hernandez et al. (2012, Paper III).

The results of Jiménez-Serra et al. (2010) implied that G035.39–00.33 comprises multiple sub-clouds. This was investigated by both Henshaw et al. (2013, Paper IV) and Jiménez-Serra et al. (2014, Paper V), who performed systematic studies of the kinematics and structure of G035.39–00.33. On the largest scales, at least three line-of-sight kinematic features are present, with each exhibiting a unique velocity and density structure (Jiménez-Serra et al. 2014). This was confirmed by Henshaw et al. (2014, Paper VI), who performed the first high angular resolution (∼5 arcsec) study of the dense gas kinematics throughout the cloud, using observations of the J = 1 → 0 transition of N₂H⁺ taken with the Plateau de Bure Interferometer (hereafter, PdBI). It was revealed that the IRDC comprises a complex network of morphologically distinct molecular filaments. Moreover, Henshaw et al. (2014) found evidence to suggest that the kinematics of the gas are locally influenced by the presence of dense, and in some cases, starless, continuum sources.

In this paper (VIII), we revisit the PdBI 3.2 mm continuum emission data, which was first presented by Henshaw et al. (2014) for qualitative comparison with the N₂H⁺ (1−0) molecular line kinematics. Our primary aim is to investigate the structure, fragmentation process, and star formation potential of G035.39–00.33 via quantitative analysis of the dust continuum emission. Details of the observations can be found in Section 2. In Section 3, we discuss the method used to systematically identify structure within the continuum data and discuss this in the context of the complex kinematics of G035.39–00.33. Section 4 contains our quantitative analysis of the continuum data. We begin with a discussion on the spatial distribution of the identified continuum cores and how this compares to predictions from theoretical work describing the fragmentation of fluid cylinders, before turning our attention to the cores themselves. In Section 5, we discuss the implications of our findings in the context of star formation throughout G035.39–00.33. Finally, in Section 6 we summarize our findings and suggest possible avenues for future research.

OBSERVATIONS

The 3.2 mm continuum observations were carried out using the Institut de Radioastronomie Millimétrique (IRAM) PdBI, France. A 6-field mosaic was used to cover the inner area of IRDC G035.39–00.33 (the dotted circles in Fig. 1 depict the primary beam of the PdBI at 3.2 mm ∼54 arcsec). The final map size is ∼40 arcsec × 150 arcsec (corresponding to ∼0.6 pc × 2.1 pc).

Figure 1.

The combined mid- and near-infrared extinction-derived mass surface density map of Kainulainen & Tan (2013, grey-scale) overlaid with 3.2 mm continuum contours. Contour levels start at 3σ_rms and increase in steps of 2σ_rms (where σ_rms ≈ 0.07 mJy beam⁻¹). Large dotted circles indicate the 6-field mosaic obtained with the PdBI. Filled magenta circles indicate the locations of high-mass cores reported by Butler & Tan (2012). Filled cyan and yellow squares refer to the high-mass (>20M_⊙) and low-mass (<20M_⊙) dense cores identified by Nguyen Luong et al. (2011) from Herschel observations. Open red circles and red triangles refer to the 8, and 24 μm emission, respectively (Carey et al. 2009; Jiménez-Serra et al. 2010). The symbol sizes are scaled by the source flux. Filled green squares highlight the location of extended 4.5 μm emission (Chambers et al. 2009).

Open in new tabDownload slide

Observations were carried out over six days in 2011 May, June, and October, in the C and D configurations (using six and five antennas, respectively) offering baselines between 19 m and 176 m. Emission on scales larger than ∼1.2(λ/D) ∼ 42 arcsec (∼0.6 pc), where D = 19 m, is filtered out by the interferometer. The 3.2 mm continuum data were cleaned using the hogbom algorithm with natural weighting. This results in a synthesized beam of {θ_maj, θ_min} = {4.3 arcsec, 3.1 arcsec} = {0.06 pc, 0.04 pc}, with a position angle of 18|$_{.}^{\circ}$|3. Line-free channels give a total bandwidth of ∼ 3 GHz. The map noise level, σ_rms, estimated from emission free regions, is ∼0.07 mJy beam⁻¹. The reference position used to determine the relative offset positions used throughout this paper is |$\alpha \,({\rm J}2000)=18^{\rm h}57^{\rm m}08{^{\rm s}_{.}}0$|⁠, δ (J2000) = 02°10′30|${^{\prime\prime}_{.}}$|0. We refer the reader to Henshaw et al. (2014) for more details on the observations.

In addition to the 3.2 mm continuum data, complementary PdBI N₂H⁺ (1−0) observations, first presented in Henshaw et al. (2014), are also utilized throughout this work. These data have been combined with existing IRAM 30 m observations to incorporate missing short spacing information into the interferometric map. Following post-processing, the spatial and spectral resolution of the N₂H⁺ (1−0) data are ∼5 arcsec and 0.14 km s⁻¹, respectively.

OBSERVATIONAL RESULTS

Structure identification using continuum data

Fig. 1 shows the spatial extent of the PdBI 6-field mosaic (dotted circles) overlaid on the combined mid- and near-infrared extinction-derived mass surface density map of G035.39–00.33 (Kainulainen & Tan 2013). The black contours highlight the 3.2 mm continuum emission. The locations of extended 4.5 (extended green objects; Cyganowski et al. 2008), 8, and 24 μm emission, indicating locations of embedded star formation, appear as green squares, red circles, and red triangles, respectively (Carey et al. 2009; Chambers et al. 2009; Jiménez-Serra et al. 2010). Qualitatively, the continuum emission appears closely related to regions of high-mass surface density (0.15 g cm⁻² < Σ < 0.32 g cm⁻²). However, there are notable exceptions to this. For instance, there is a lack of continuum emission towards {Δα, Δδ} = { − 20.0 arcsec, 70.0 arcsec}. Such discrepancies may be due to a lack of sensitivity (our 3σ_rms column density sensitivity is N_H ≳ 10²³cm⁻²; see Section 4.2) or the result of missing short spacings in our interferometric map.

The continuum emission is highly structured. There are a number of prominent emission peaks arranged along the major axis of the IRDC. To investigate this further, we use dendrograms (Rosolowsky et al. 2008). Specifically, our analysis makes use of astrodendro, a python package used to compute dendrograms of astronomical data.¹ As well as providing a systematic approach to structure identification, dendrograms are also less sensitive to variation in the input parameters in comparison to alternative methods (Pineda, Rosolowsky & Goodman 2009). Additionally, dendrograms can be used to identify hierarchical structure, which is desirable in complex regions. The following parameters are used in computing the dendrogram: |${\rm min\_value}= 3\sigma _{\rm rms}$| (the minimum intensity considered in the analysis); |${\rm min\_delta}= 2\sigma _{\rm rms}$| (the minimum spacing between isocontours; using |${\rm min\_delta}= 1\sigma _{\rm rms}$| has no effect on the identified structure); |${\rm min\_npix}=26$| (the minimum number of pixels contained within a structure). The angular resolution of our observations is used to determine |${\rm min\_npix}=\frac{2\pi }{8{\rm ln}(2)}\frac{\theta _{\rm maj}\theta _{\rm min}}{A_{\rm pix}}$|⁠, where A_pix is the area of a single (0.76 arcsec × 0.76 arcsec) pixel.

Fig. 2 shows the result of the dendrogram analysis. Here we highlight the location of the dendrogram ‘leaves’, the highest level of the dendrogram hierarchy, representing the smallest structures identified. A total of 13 leaves are identified. Each is denoted by an ID number designated in order of increasing offset declination. Equivalent radii of the leaves range between 0.04 pc < R_eq < 0.07 pc (with a median value of R_eq = 0.05 pc; note that |$\sqrt{\theta _{\rm maj}\theta _{\rm min}}\sim 3.7$| arcsec which is ∼0.05 pc), whereby R_eq ≡ (N_pixA_pix/π)^1/2, and N_pix is the number of pixels associated with a given leaf. The equivalent radius refers to that of a circle which covers an area equivalent to N_pixA_pix. These radii are, on average, 30 per cent larger than the geometric mean radii computed from the intensity-weighted second moment output from astrodendro (see Table 1 and discussion by Rosolowsky & Leroy 2006). Selecting the equivalent radius therefore represents a conservative approach to estimating the physical properties of the leaves (e.g. number densities; Section 4.2). The median aspect ratio of the dendrogram leaves is 2.2. However, leaves #4 and #5 appear to be more filamentary, with aspect ratios >4 (although it is questionable as to whether these leaves represent single structures; see Section 3.2). There is some suggestion that certain leaves exhibit substructure, evidenced through either secondary continuum peaks or irregular boundaries (e.g. leaves #2, #6, #13). This substructure is rejected by the algorithm, following the insertion of physically motivated input parameters selected to reflect the limitations of our observations. Only when |${\rm min\_npix}=18$| (i.e. when the pixel threshold is reduced below the number of pixels contained within one beam) does the number of structures identified deviate from that presented here (leaf #13 splits into two). Even when |${\rm min\_npix}$| is reduced by a factor of >2, 11 out of the original 13 identified leaves remain unaffected. This gives us confidence that our results are robust, and relatively insensitive to the input parameters of the dendrogram algorithm. The offset locations of peak emission, semiminor and semimajor axes (R_min, R_maj) and their geometric mean, projected aspect ratios (⁠|$A\!\!{\rm R}$| ≡ R_maj/R_min), areas (N_pixA_pix), and equivalent radii (R_eq) of the leaves are listed in Table 1.

Figure 2.

Highlighting the projected location of the dendrogram leaves discussed in Section 3. Each leaf is denoted with an ID number and a coloured contour. Information relating to each leaf can be found in Table 1. The background image is the mass surface density map of Kainulainen & Tan (2013) and all symbols are defined in Fig. 1.

Open in new tabDownload slide

Comparison with molecular line observations

To complement the structure-finding algorithm employed in Section 3.1, we examine the molecular line kinematics associated with each dendrogram leaf. This investigation focuses on the isolated (F₁, F = 0, 1 → 1, 2) hyperfine component of N₂H⁺ (1−0), using the PdBI data first presented by Henshaw et al. (2014).²

To examine the kinematics, we first generate a spatially averaged spectrum from all spectra contained within the boundary defining the maximum (projected) physical extent of each leaf. Henshaw et al. (2014) find that the N₂H⁺ (1−0) emission observed throughout G035.39–00.33 can be attributed to a complex network of filamentary structures. Since the velocity separation between these sub-clouds is <1 km s⁻¹ (comparable to their mean FWHM line-widths) and velocity gradients of magnitude 1.5–2.5 km s⁻¹ pc⁻¹ are observed throughout each, analysing the N₂H⁺ (1−0) data using dendrograms is not trivial. Gaussian profiles are therefore fitted to the isolated component of all features observed within each spatially averaged spectrum.³ The line emission is then integrated over the velocity range, [v_{0, i} − (Δv_i/2)] < v_{0, i} < [v_{0, i} + (Δv_i/2)], where v_{0, i} and Δv_i are the centroid velocity and full width at half-maximum (FWHM) line-width of the ith velocity component, respectively. This enables us to examine the spatial distribution of N₂H⁺ emission associated with each identified velocity component, and relate this to the relevant continuum emission peak.

Out of the 13 spatially averaged spectra extracted from within the leaf boundaries, 8 exhibit multiple velocity components. Typically, where multiple spectral features are evident, the emission from one will dominate over the other(s). This allows us, albeit simplistically, to link the continuum sources to a particular kinematic feature. A good example of this is leaf #9, shown in the top panels of Fig. 3. The left-hand panel is the spatially averaged spectrum extracted from the dashed box shown in the centre and right-hand panels. While the continuum emission appears monolithic, two velocity components are evident in the N₂H⁺ spectrum. The centre and right-hand panels show the spatial distribution of the N₂H⁺ emission associated with each velocity component. The region covered by the integrated emission in both instances is greater than that covered by the leaf alone, indicating that the N₂H⁺(1−0) emission is extended and not exclusively associated with the continuum peak for leaf #9. Although emission from the velocity component identified at v_{0, 1} = 45.4 km s⁻¹ is spatially coincident with the leaf (centre panel), the component at v_{0, 2} = 46.6 km s⁻¹ dominates (in terms of spatial coverage; right-hand panel). We therefore speculate that the majority of flux attributed to leaf #9 is associated with the high(er)-velocity component.

Figure 3.

Top: the spectral and spatial distribution of N₂H⁺ (1−0) emission associated with leaf #9. The left-hand panel is a spatially averaged spectrum, showing the isolated (F₁, F = 0, 1 → 1, 2) hyperfine component of N₂H⁺ (1−0). The spectrum has been extracted from the black dashed box seen in the centre and right-hand panels. Two spectral components are clearly evident. The solid black and orange Gaussian profiles represent the best-fitting model solution to the data (the orange line refers to the component most closely associated with the leaf extracted from the 3.2 mm continuum data). The centre and right-hand panels display the spatial distribution of each emission feature. The plots are shown on an equivalent scale for ease of comparison, however, this has led to saturation in both central panels. The black contours are equivalent to those in Fig. 1 and the orange contour corresponds to the boundary of leaf #9. Bottom: each panel is equivalent to those shown in the top panels but for leaf #5. Note how the northern and southern portions of the dendrogram leaf appear to be attributed to two independent velocity components which overlap spatially. The black Gaussian profiles reflect the fact that the continuum flux accredited to leaf #5 cannot be attributed to a single structure (see text for details regarding the treatment of such cases).

Open in new tabDownload slide

An exception to this is presented in the bottom panels of Fig. 3. In contrast to leaf #9, the spatial distribution of each N₂H⁺ emission feature associated with leaf #5 (see the bottom-left panel) leads us to question whether this leaf represents a single structure with complex internal kinematics, or if this leaf, in fact, comprises two independent structures that overlap in projection. The southern portion of the leaf appears to be associated with the component identified at v₀ = 45.6 km s⁻¹ (centre panel). However, the northern portion, including the peak in 3.2 mm continuum emission, is associated with the higher velocity component at v₀ = 46.6 km s⁻¹(right-hand panel). As the top and bottom panels of Fig. 3 are continuous in declination, we can see that the northern tip of leaf #5 appears to be associated with a coherent structure that extends to, and beyond, leaf #9 (F3; Henshaw et al. 2014).

We stress that the inclusion of velocity information does not completely alleviate problems associated with projection effects (see e.g. Beaumont et al. 2013). This becomes particularly pertinent in an environment such as G035.39–00.33, which exhibits complex morphological structure and kinematics (e.g. Henshaw et al. 2013, 2014; Jiménez-Serra et al. 2014). However, the above analysis does highlight the importance of demonstrating caution when identifying structure in two-dimensional data. Repeating the above analysis for all identified leaves, we use the kinematic information as a rough guide to determine which leaves are to be analysed further in Section 4. The leaves are split into four categories: (1) leaves which exhibit single velocity components are retained throughout all analysis (leaves #1, #2, #3, #10, #12); (2) leaves which show multiple velocity components, but where one of these dominates (in terms of either the magnitude of, or spatial coverage of, the integrated N₂H⁺ emission), are retained, and the dynamical properties of the leaves (see Section 4.3) are estimated using only the kinematic properties of the dominant component (leaves #7, #9, #11); (3) leaves which exhibit multiple velocity components, but where the continuum emission cannot be unambiguously attributed to a single velocity component, are retained, and the dynamical properties of the leaves are estimated using the kinematic properties of all components (leaf #8); (4) leaves which exhibit multiple velocity components, but where different portions of the continuum emission within the leaf boundary may be associated with different velocity components, are rejected from all analysis that relies on intrinsic geometrical assumptions and/or assumes there is no underlying substructure (leaves #5, #6, #13 and possibly #4). We refer the reader to Appendix A for a more complete description of each dendrogram leaf and its associated kinematics.

CLIMBING THE STRUCTURE TREE: INVESTIGATING THE STRUCTURE AND FRAGMENTATION OF G035.39–00.33

Investigating the fragmentation of a filament

We begin our analysis at the foot of the structure tree, focusing on the spatial distribution of continuum emission throughout G035.39–00.33. There are several examples in the literature of filamentary molecular clouds that exhibit a quasi-regular spacing of ‘cores’ (e.g. Jackson et al. 2010; Miettinen et al. 2012; Busquet et al. 2013; Kainulainen et al. 2013; Takahashi et al. 2013; Lu et al. 2014; Wang et al. 2014; Beuther et al. 2015b; Ragan et al. 2015; Contreras et al. 2016; Teixeira et al. 2016). This regularity is often discussed in the context of predictions from theoretical work describing the fragmentation of fluid cylinders due to gravitational or magnetohydrodynamic-driven instabilities (e.g. Chandrasekhar & Fermi 1953; Nagasawa 1987; Inutsuka & Miyama 1992; Nakamura, Hanawa & Nakano 1993, 1995; Tomisaka 1995). In this theoretical framework, the characteristic spacing between fragments is defined by the wavelength of the fastest growing unstable mode of the fluid instability.

The 3.2 mm dust continuum emission associated with G035.39–00.33 is distributed along the major axis of the filamentary IRDC (see Fig. 1). To quantify the spatial separation between the dendrogram leaves identified in Section 3.1 we use the minimum spanning tree (MST) method (Prim 1957). An MST is a graph theory construct that identifies the shortest possible total path-length between a set of points where there are no closed loops. MSTs are frequently used to quantify the relative spatial distributions of both stars and gas in simulated and observed star-forming regions (e.g. Cartwright & Whitworth 2004; Allison et al. 2009; Gutermuth et al. 2009; Lomax, Whitworth & Cartwright 2011; Parker & Dale 2015). MSTs have the advantage over other methods that they are not biased by the inherent geometry of the region. The reference point for each dendrogram leaf is taken as the location of peak emission (see Table 1), and these points are then used to construct the MST. Fig. 4 includes a box plot and the corresponding cumulative distribution function of the MST lengths. The mean angular separation between dendrogram leaves according to the MST is ∼12.8 arcsec (with ∼50 per cent of all values falling within a factor of ∼2 of the mean), which corresponds to a projected physical distance of λ_obs ∼ 0.18 pc.⁴

Figure 4.

Top: a box plot of the angular leaf separation. This highlights the range in separation, the interquartile range (the box itself), the median separation (∼11.1 arcsec; thick vertical line), and the mean separation (∼12.8 arcsec; vertical dashed line). Botttom: the corresponding cumulative distribution.

Open in new tabDownload slide

Using the mass surface density map presented in Fig. 1, Hernandez et al. (2012) estimate the filament number density, n_{H, f}, of G035.39–00.33 (see their Table 1). The average value of n_{H, f}, over the region mapped with the PdBI, is n_{H, f} ∼ 0.2 × 10⁵ cm⁻³ (note that this assumes that the filament is inclined by 30°, with respect to the plane of the sky, 0°). Nguyen Luong et al. (2011) estimate the dust temperature throughout G035.39–00.33, by fitting pixel-by-pixel modified blackbody spectral energy distributions derived from Herschel observations (excluding the 70 μm emission). They find dust temperatures ranging between 13 and 16 K, with the lowest temperatures observed towards the centre of G035.39–00.33. However, the angular resolution of the Herschel temperature map is 37 arcsec. These observations are therefore sensitive to more diffuse (and warm) emission than traced by our high angular resolution PdBI observations (flux contributions from scales >42 arcsec are filtered out by the interferometer; see Section 2). Thermodynamic coupling between the gas and dust in molecular clouds is valid for densities ∼10⁵ cm⁻³ (e.g. Goldsmith 2001; Glover & Clark 2012), which is greater than the value estimated by Hernandez et al. (2012). However, in the absence of both complementary gas temperature and high angular resolution dust temperature estimates, we assume that the gas and dust temperatures are approximately equal, and adopt a fiducial value of T_f = 15 K.

One aspect not factored into the above analysis is the effect of inclination. The true core spacing throughout a filament inclined at an angle, i, with respect to the plane of the sky (i = 0°) is λ_{obs, i} = λ_obs/cos(i). Hernandez et al. (2012) assume an inclination angle of 30° when determining the number density of the filament. Following this assumption results in an ‘inclination-corrected’ leaf spacing of λ_{obs, i} = 0.21 pc. Setting λ_frag = λ_{obs, i} in equation (2), we find that a density of n_{H, f} ∼ 3.3 × 10⁵ cm⁻³, is required to reproduce the inclination-corrected spacing.

Evidently there exists a discrepancy between the observed core separation and that predicted by this particular model. When considering thermal motions only, the density required to reproduce the observed spacing is ∼ an order of magnitude greater than the value estimated by Hernandez et al. (2012). Including turbulent gas motions only exacerbates the problem, requiring densities that are, in fact, similar to the typical density of the leaves (see Section 4.2). Putting this another way, for a fixed density of n_f = 0.2 × 10⁵ cm⁻³, the theoretical fragment spacing is a factor of ∼8 greater than the observed value when considering both thermal and turbulent support (∼5, when considering only thermal support). This is similar to the conclusion of Pillai et al. (2011), who find that, in order to explain the fragment spacing in massive star-forming regions G29.96–0.02 and G35.20–1.74, gas densities similar to their estimated core densities are required. Below we discuss two alternate scenarios which may explain the discrepancy in G035.39–00.33. First of all, we discuss the possibility that the presence of dynamically important magnetic fields may influence the fragmentation length-scale (a scenario that was recently explored by Contreras et al. 2016). Secondly, we discuss how the geometric interpretation inherent to the above discussion, i.e. that G035.39–00.33 is well represented as a single cylindrical filament, may not be applicable.

The left-hand panel of Fig. 5 demonstrates how λ_{frag, f} varies as a function of density according to equation (5), for both the thermal (blue) and thermal+non-thermal (black) cases. The dotted, dashed, and dot–dashed refer to where β = 0.1, 1.0, 10.0. The solid vertical and horizontal lines indicate the observed filament density (n_{H, f} = 0.2 × 10⁵ cm⁻³; Hernandez et al. 2012) and inclination-corrected leaf spacing (∼0.21 pc), respectively. The right-hand panel of Fig. 5 shows how the density required for λ_{frag, f} = λ_{obs, i} changes as a function of the magnetic field strength according to equation (5). Solid lines indicate the locus where λ_{frag, f} = λ_{obs, i}, once again, for both the thermal (blue) and thermal+non-thermal (black) cases.

Figure 5.

Left: this figure shows how the wavelength of the fastest growing mode of the magnetohydrodynamic fluid instability discussed in Section 4.1 changes as a function of density according to equation (5) (Nakamura et al. 1993, 1995) for both the thermal (blue) and thermal+non-thermal (black) cases. The dotted, dashed, and dot–dashed refer to where β = 0.1, 1.0, 10.0. For dynamically important magnetic fields (where β = 0.1; dotted lines), the wavelength of the fastest growing mode is less than in the regime where magnetic fields are unimportant (β = 10; dot–dashed lines). The solid vertical and horizontal lines indicate the observed filament density (n_{H, f} = 0.2 × 10⁵ cm⁻³; Hernandez et al. 2012) and inclination-corrected leaf spacing, respectively. Right: this figure shows how the filament density required for λ_{frag, f} = λ_{obs, i} changes as a function of magnetic field strength. The solid line(s) indicate where λ_frag = λ_{obs, i}. The dotted, dashed, and dot–dashed refer to the loci of β = 0.1, 1.0, 10.0. For dynamically unimportant magnetic fields (β = 10; dot–dashed lines), the density required for λ_{frag, f} = λ_{obs, i} (solid lines) reverts to that derived using equation (1). This figure demonstrates that even with a dynamically important magnetic field, the reduction in the wavelength of the fastest growing mode of the instability is insufficient to explain the observed leaf spacing.

Open in new tabDownload slide

The left-hand panel demonstrates how the wavelength of the most unstable perturbation is shorter when the ratio of the magnetic pressure to the gas pressure is higher (i.e. when β is small), for a fixed density (Nakamura et al. 1993). Note however, that the effect is small. Conversely, the right-hand panel shows how the density required for λ_{frag, f} = λ_{obs, i} can be reduced if the magnetic field is dynamically important (0.1 < β < 1.0). As can be inferred from the right-hand panel, the density required for λ_{frag, f} = λ_{obs, i} in the case of thermal+non-thermal (thermal) fragmentation is 6.2 × 10⁵ cm⁻³ (2.1 × 10⁵ cm⁻³) for a magnetic field strength of 470 μG (140 μG), compared with 9.4 × 10⁵ cm⁻³ (3.3 × 10⁵ cm⁻³) for a dynamically unimportant magnetic field (see above). This figure demonstrates that the decrease in λ_{frag, f} according to equation (5), due to the presence of a dynamically important longitudinal magnetic field, on its own, cannot account for the observed discrepancy.

An alternative possibility is that the underlying assumption of the above model, that G035.39–00.33 can be described simplistically as a single cylindrical filament, may be a poor one. We stress that this is not necessarily mutually exclusive from a scenario which includes dynamically important magnetic fields (the above analysis only accounts for a very specific configuration of magnetic field). However, in this example, the discrepancy may relate to the structure of the cloud itself. Although low-angular resolution dust continuum maps may hint towards a relatively ‘simplistic’ cloud morphology, the reality is anything but simple. Henshaw et al. (2014) argue that G035.39–00.33 is organized into a serpentine network of morphologically distinct molecular sub-filaments. Each sub-filament displays not only unique kinematic properties (in terms of σ_v and a complex pattern of velocity gradients), but also its own density structure, as demonstrated by Jiménez-Serra et al. (2014, albeit these densities are derived from CO observations with ∼20 arcsec resolution).

The presence of multiple line-of-sight structures, if confirmed, influences our investigation in two key ways. First, if, as suggested by Henshaw et al. (2014), the leaves are associated with otherwise independent molecular filaments then we may underestimate the true spacing. The density required to reproduce the true spacing may therefore be lower than that stated above (assuming the velocity dispersion remains constant; equation 4). Secondly, the relevant properties in equation (4) (n_{H, f}, inclination, σ_v) should be those relating to the individual sub-filaments. The fiducial value of the density assumed above (n_{H, f} = 0.2 × 10⁵ cm⁻³) is derived from the mass surface density map (Fig. 1), assuming cylindrical geometry with a radius of R_f ≈ 30 arcsec or ∼0.4 pc at a distance of 2900 pc (Hernandez et al. 2012), which may not reflect the central density of the sub-filaments. From equation (4), a factor of 10 increase in the density would give λ_frag ∼ 0.45 pc ∼ 2λ_{obs, i}.

Estimating the physical properties of the dendrogram leaves

Initial considerations

Having focused on the distribution of continuum emission throughout G035.39–00.33, we now turn our attention to analysing the leaf properties. It is worth noting that the combination of missing continuum flux in our interferometric map and complex kinematic structure (discussed in Sections 3.2 and 4.1, and more fully in Henshaw et al. 2014), makes it difficult to unambiguously apportion flux to any given continuum source. We therefore employ two different approaches to estimating the physical properties of the dendrogram leaves.

Both methods make an underlying assumption regarding the translation of a region of emission, defined by an isosurface, into physical three-dimensional space (we refer the reader to Rosolowsky et al. 2008 for a more complete description of the philosophy behind these methods). The first approach, which is most conservative, makes the assumption that each leaf represents a discrete object superimposed on a background of flux, |$F^{\rm bg}_{\nu }$| (which is subtracted from each leaf pixel prior to the estimation of physical properties). This has been used in recent studies as a way of accounting for the fact that emission from more diffuse or larger scale structures can contaminate the flux of small-scale structures (e.g. Ragan, Henning & Beuther 2013; Pineda et al. 2015). The second approach assumes that there is no background contribution, and that all of the flux within the leaf boundary is attributed to that structure. The reality probably lies somewhere in between, and these two methods provide lower and upper bounds to the source flux, respectively. We present results from both approaches throughout the following analysis and denote the background-subtracted properties with ‘b’ (see Table 2).

It is also prudent, prior to the determination of physical properties, to estimate the contribution to the continuum flux from free–free emission originating from embedded radio sources. To estimate the free–free contribution at 93 GHz, we inspect images from The Coordinated Radio and Infrared Survey for High-mass Star Formation (CORNISH) survey (Hoare et al. 2012; Purcell et al. 2013). We identify one 5σ source (∼2 mJy at 5 GHz) at a position |$\alpha \,(J2000)=18^{\rm h}57^{\rm m}08{^{\rm s}_{.}}37$|⁠, δ (J2000) = 02°10′32|${^{\prime\prime}_{.}}$|71, corresponding to an offset location of {Δα, Δδ} = {5.8 arcsec, 2.2 arcsec} in our PdBI map. We note however, that due to artefacts in the CORNISH images, reliable source detections are limited to ≥7σ. Since the aforementioned 5σ source does not coincide with one of the 3.2 mm continuum peaks, nor is there evidence for 8, 24, 70 μm emission at this location (see Fig. 1 and Nguyen Luong et al. 2011), it is possible that this is an image artefact in the CORNISH map. In the optically thin regime, free–free emission has a frequency dependence of S_ν ∝ ν^−0.1. Based on the rms noise of the CORNISH images (0.37 mJy at 5 GHz), we estimate an upper limit to the free–free contribution of 0.28 mJy at 93 GHz. Since no other detections are made, we expect this to be a fairly generous upper limit and anticipate that the contribution to our PdBI continuum flux from free–free emission is small.

Estimating the physical properties

The dust opacity per unit mass is determined from κ_ν = κ₀(ν/ν₀)^β, assuming a dust emissivity index, β, where κ₀ is based on the moderately coagulated thin ice mantle dust model of Ossenkopf & Henning (1994) at a frequency, ν₀. At a frequency of ∼93 GHz, we adopt a value of κ_ν ≈ 0.186 cm²g⁻¹ (extrapolating from κ₀ = 0.899 cm²g⁻¹ at ν₀ = 230 GHz with β = 1.75; e.g. Battersby et al. 2011). From Draine (2011, table 23.1), the hydrogen-to-(refractory-component-)dust-mass ratio is R_{gd, H} ∼ 100. Therefore we adopt a value of R_gd = 141 for the total (gas plus dust)-to-(refractory-component-)dust-mass ratio (assuming a typical interstellar composition of H, He, and metals). As discussed in Section 4.1, there are currently no gas or dust temperature estimates for G035.39–00.33 at a resolution equivalent to those studied here. For the leaf analysis, we assume that T_c < T_f, that T_d = T_c, and T_d = 13 K (at the lower end of the range derived by Nguyen Luong et al. 2011). The corresponding column density sensitivity derived from our 3σ_rms flux level of 0.21 mJy beam⁻¹ is N_H ∼ 1.3 × 10²³ cm⁻².

The derived column densities range from 3.6 × 10²³ cm⁻² < N_{H, c} < 8.0 × 10²³ cm⁻² (⁠|$1.5\times 10^{23}\,{\rm cm^{-2}}<N_{\rm H,c}^{\rm b}<5.9\times 10^{23}\,{\rm cm^{-2}}$|⁠), with a mean value of N_{H, c} ∼ 4.7 × 10²³ cm⁻² (⁠|$N^{\rm b}_{\rm H,c}\sim 2.5\times 10^{23}\,{\rm cm^{-2}}$|⁠). Given the uncertainties in the flux calibration (∼10 per cent), dust opacity per unit mass (∼30 per cent; accounting for different degrees of coagulation in the Ossenkopf & Henning 1994 models), total (gas plus dust)-to-dust mass ratio (∼30 per cent), and temperature (±3 K), the uncertainty in the derived beam-averaged column density is ∼50 per cent (after summing in quadrature). These values, as well as those estimated below, can be found in Table 2.

For comparison, we also estimate the mass from the mass surface density map of Kainulainen et al. (2013). We extract this mass estimate, M_MIREX, from within the boundary defining the maximum (projected) physical extent of each leaf (see Fig. 3 and those in Appendix A). We find 4.8 M_⊙ < M_MIREX < 26.8 M_⊙. Comparing the masses of directly, we find 0.4 < M_MIREX/M_c < 1.1, with an average value of 〈M_MIREX/M_c〉 ∼ 0.84. Due to the lack of short spacings in our interferometric map, we caution against drawing firm conclusions from this comparison. However, the fact that the estimates agree (within the uncertainties) gives us confidence in the reliability of our continuum derived masses. Finally, comparing the total mass of (all) the leaves with the mass of the inner filament estimated by Hernandez et al. (2012), we find that the leaves make up ∼10 per cent of the total mass of the region.

Dynamical properties of the dendrogram leaves

Using the physical properties of the dendrogram leaves derived in Section 4.2, we can now assess whether the leaves themselves are susceptible to gravitational collapse (and potentially further fragmentation). In the following sections, we evaluate the support provided by different mechanisms.

Thermal support

Thermal + turbulent support

Deviations from spherical symmetry are accounted for in a_θ. Bertoldi & McKee (1992) consider a triaxial ellipsoid with equatorial radius, R, and an extent in the third dimension, 2Z, such that the aspect ratio is y = Z/R. They show that for log₁₀(Z/R) < |1|, a_θ ≈ 1.0 ± 0.2. Consequently, we ignore the effect of clump elongation in the following analysis.

Fig. 6 displays the radial flux density profiles for leaves #5 and #9 (as in Fig. 3). The left-hand panel displays the radial flux density profile for leaf #5. For reasons discussed in Section 3.2, this leaf is not included in the analysis. However, in addition to the fact that the northern and southern portions of the leaf can be attributed to different kinematic components, this image supports our decision to reject this leaf. As can be seen from Fig. 6, the flux density profile exhibits additional peaks, possibly signifying the presence of underlying substructure (evidence, in this example, for the superposition of fragments associated with different sub-filaments). The right-hand panel displays the radial flux density profile for leaf #9. In contrast to the profile of leaf #5, the flux density decreases smoothly as a function of radius. The implication is that the leaf is monolithic (at the spatial resolution of our PdBI observations). The remaining radial flux density profiles can be found in Appendix A. The light and heavy lines in Fig. 6 signify the synthetic radial flux density profiles before and after beam correction, respectively. The solid orange line represents our best-fitting model solution to the radial flux density profile for leaf #9 (κ_ρ = 2.14). On average, we find κ_ρ = 1.9.

Figure 6.

Flux density as a function of equivalent radius for leaves #5 (left-hand panel) and #9 (right-hand panel). Black lines are model radial flux density profiles before (light) and after (heavy) beam correction (see Section 4.3.2). Each line corresponds to a density profile of the form |$F_{\rm \nu }\propto r^{-\kappa _{\rho }}$|⁠. The corresponding values of κ_ρ are included in the legend. The dark grey block indicates the extent of the geometric mean radius of the synthesized PdBI beam (1.83 arcsec or 0.025 pc at a distance of 2900 pc). The light grey block is twice this value. The solid orange line in the right-hand panel indicates the closest-matching model solution to the observations. The lack of a corresponding model solution in the left-hand panel reflects the fact that leaf #5 has been rejected from our analysis. This figure supports this decision, since the additional peaks in the radial flux density profile may signify the presence of underlying substructure. Radial flux density profiles for the additional leaves can be found in Appendix A.

Open in new tabDownload slide

The above considerations enable us to define a critical value, α_{vir, cr} ≡ a, which serves as a gauge to assess the stability of a cloud fragment. In this simplistic formalism, cloud fragments with α_vir < α_{vir, cr} are susceptible to gravitational collapse in the absence of additional internal support and cloud fragments with α_vir > α_{vir, cr} may expand in the absence of pressure confinement (more detailed stability analysis, accounting for the effects of surface pressure, finds α_{vir, cr} ≈ 2 for isothermal, non-magnetized cloud fragments in equilibrium; McKee & Holliman 1999; Kauffmann, Pillai & Goldsmith 2013). Normalization of α_vir by a gives α_{vir, cr} = 1, allowing for ease of comparison between virial parameters of a cloud fragment estimated assuming different density profiles. For simplicity, we do not carry out this analysis on the background-subtracted leaves, because of complications in estimating background-subtracted velocity dispersions. For leaf masses without background-subtraction, we find 0.3 < α_vir < 2.1 (with a mean value, 〈α_vir〉 = 1.0) and 0.2 < α_vir < 1.2 (〈α_vir〉 = 0.6) for κ_ρ = 0 and κ_ρ = 2.0, respectively. Incorporating our density profile analysis returns virial parameters spanning the range 0.2 < α_vir < 1.3 (with a mean value, 〈α_vir〉 = 0.7). Table 3 lists the virial parameters estimated for each of the dendrogram leaves.

This analysis indicates that, when taking σ_NT as an upper limit on the level of turbulent support, all dendrogram leaves are consistent with being either sub-virial or approximately virial (α_vir ≲ 1 to within a factor of ≳ 2 uncertainty for κ_ρ = 1.5–2). In the absence of additional support, dendrogram leaves that are strongly sub-virial should undergo fairly rapid (∼t_{ff, c}) global collapse. Leaves #7 and #9 are consistent with this picture (α_vir ∼ 0.4 and ∼0.2, respectively), and other regions of massive star formation where low virial parameters have been reported (e.g. Csengeri et al. 2011; Pillai et al. 2011; Li et al. 2013; Peretto et al. 2013; Battersby et al. 2014; Beuther et al. 2015b; Lu et al. 2015). However, low virial parameters such as these may instead be indicative of strong magnetic support (as suggested by e.g. Tan et al. 2013).

Thermal + turbulent + magnetic support

IMPLICATIONS FOR STAR FORMATION WITHIN G035.39–00.33

In Section 4.1, we find that the spacing between continuum sources throughout G035.39–00.33 is significantly (a factor of ∼8) smaller than that predicted by gravitational instabilities in purely hydrodynamical fluid cylinders. At face value this appears to suggest that magnetic fields may have a significant role to play in the reconciliation of the observed and predicted spatial distribution. However, complex line-of-sight structure and the presence of sub-filaments observed throughout G035.39–00.33, may make a significant contribution to the discrepancy.

The idea that the continuum sources may be associated with different sub-filaments is qualitatively supported by leaves that appear close to one another in projection, but show clear differences in their radial velocities. For example, leaves #9 and #10 are separated by a projected distance of ∼0.15 pc but their line centroid velocities differ by ∼0.8 km s⁻¹ (see Table 3). Conversely, leaves #9 and #11 have a similar spatial separation but their line centroids differ by ∼0.1 km s⁻¹. In this particular example, leaves #9 and #11 seem to be consistent with the mean velocity of filament F3 of Henshaw et al. (2014), (46.86 ± 0.04) km s⁻¹, whereas leaf #11 is consistent with the mean velocity of F2b, (46.00 ± 0.05) km s⁻¹ (filament F2a has a mean centroid velocity of [45.34 ± 0.04] km s⁻¹, for reference). If the continuum sources can indeed be attributed to different sub-filaments, then the assumption that G035.39–00.33 can be described simplistically, as a single cylindrical filament, is invalid. If confirmed, the observed spacing is most likely influenced by a combination of several important factors, including, the number of sub-filaments, differences in the individual sub-filament properties (e.g. density, inclination, velocity dispersion), and the strength and orientation of the magnetic field.

The origins of the complex physical and kinematic gas structure of G035.39–00.33 are currently unknown. Whether the observed sub-filaments are a result of the fragmentation process (cf. the ‘fray and fragment’ scenario proposed by Tafalla & Hacar 2015) or whether they first formed at the stagnation points of a turbulent velocity field and have been brought together by gravitational contraction on larger scales (cf. the ‘fray and gather’ scenario proposed by Smith et al. 2016), remains an open question. However, the presence of widespread emission from shocked gas tracers (e.g. SiO; Jiménez-Serra et al. 2010) throughout G035.39–00.33 (and other molecular clouds, e.g. Nguyen-Luong et al. 2013; Duarte-Cabral et al. 2014), may point towards a dynamical origin. The prevalence of sub-filaments in many observational studies (e.g. Hacar et al. 2013; Peretto et al. 2013, 2014; Alves de Oliveira et al. 2014; Fernández-López et al. 2014; Lee et al. 2014; Panopoulou et al. 2014; Dirienzo et al. 2015) emphasizes the importance of considering the underlying physical structure in any fragmentation analysis. Putting this another way, this highlights the danger of using simple geometric models, without prior consideration of the kinematic information.

The formation of sub-filaments, followed by the formation of cores native to those sub-filaments (and potentially, further fragmentation of those cores), may signify a multilayered fragmentation process within G035.39–00.33, similar to that proposed in other molecular clouds (e.g. Teixeira et al. 2006; Hacar et al. 2013; Kainulainen et al. 2013; Takahashi et al. 2013; Wang et al. 2014; Beuther et al. 2015a; Tafalla & Hacar 2015). We find that the majority of the cores within G035.39–00.33, including the two most massive objects (leaves #7 and #9), are located towards the H6 region ({Δα, Δδ} ∼ {3 arcsec, 20 arcsec}; Butler & Tan 2012). This is also the location at which several of the sub-filaments meet (Henshaw et al. 2013, 2014; Jiménez-Serra et al. 2014), which is reminiscent of studies highlighting the formation of star clusters at the junctions of complex filamentary systems (e.g. Myers 2009; Schneider et al. 2012; Peretto et al. 2014), and consistent with simulations (e.g. Dale & Bonnell 2011; Myers et al. 2013; Smith et al. 2013).

Interestingly, the analysis presented in Section 4.2 shows that there is only a factor of ∼3 difference between the highest mass (leaf #9; ∼24 M_⊙) and lowest mass (leaf #11; ∼8 M_⊙) cores identified in this study. This is in spite of the fact that our PdBI data are theoretically sensitive to masses that are a factor of ∼4 lower than this (∼2.0 M_⊙).⁶ The steep slope of the locally invariant stellar initial mass function implies that many low-mass stars form within clusters alongside high-mass stars (Bastian, Covey & Meyer 2010; Offner et al. 2014). Assuming that the mass distribution of pre-stellar cores, the core mass function, takes a form dN/d(log m) ∝ M^−Γ, where Γ = 1–1.5, for masses M > 0.5 M_⊙ (e.g. Motte, Andre & Neri 1998), it follows that for every 25 M_⊙ core (equivalent to the mass of the most massive leaf detected in this paper) one might expect to find 10 ± 3 cores in the range 2–8 M_⊙, i.e. the mass range covering our sensitivity and the lowest mass leaf detected in the present investigation.

In a recent study by Zhang et al. (2015), who present Atacama Large Millimeter Array observations of IRDC G28.34 + 0.06, an apparent dearth of low-mass dense cores was also noted. Zhang et al. (2015) explain that it would be counterintuitive for stars to form first in the lower density regions surrounding massive clumps within which high-mass stars are forming in G28.34 + 0.06, and instead favour the interpretation that low-mass cores and stars form at a later stage, after the formation of massive stars. However, this is in contrast to the work of Foster et al. (2014), who detect a population of low-mass protostars in the IRDC G34.43 + 00.24. This newly identified population of low-mass stars is situated in the interclump medium of the filamentary cloud. Their presence leads the authors to suggest that the population of low-mass stars may have formed before, or perhaps coevally with, the high-mass stars.

Close inspection of the continuum map presented in Fig. 2, and studying the radial flux density profiles of the dendrogram leaves (see Appendix A), indicates that several of the leaves exhibit substructure. Our ability to detect lower mass cores may therefore be limited by our angular resolution. To assess this further, in Section 4.3 we sought to establish whether the leaves which are kinematically coherent (i.e. those leaves which can be attributed to a single velocity structure but may harbour underlying substructure in the continuum), are susceptible to collapse, and potentially further fragmentation. This analysis reveals that in the absence of additional support, possibly from magnetic fields with strengths of the order 230 μG < B < 670 μG (determined by equating the leaf masses with a critical core mass that incorporates the effect of both gas and magnetic pressure in providing support to the fragment; Section 4.3.3), the leaves may collapse.

The above implies that our ascent of G035.39–00.33's structure tree is not yet complete. Future, high angular resolution (approaching the Jeans length of the individual leaves, λ_{J, c} ≲ 0.03 pc; Section 4.3) and high sensitivity continuum observations (Henshaw et al. 2016b), as well as observations of molecular lines with higher critical densities, are needed to probe further sub-fragmentation (akin to that observed in low-mass cores; e.g. Pineda et al. 2011, 2015). This will determine whether the lack of cores identified between 2 and 8 M_⊙ has a physical origin or if this can be explained by observational bias. Such observations will also aid in testing the predictions of hydrodynamical simulations of collapsing cloud cores, which show an increased level of fragmentation in cores with shallower (ρ ∝ r⁻¹) density profiles compared to those with steeper (ρ ∝ r⁻²) profiles (Girichidis et al. 2011). Specifically, this will help to establish the fate of leaves #7 and #9, which appear centrally concentrated (ρ ∝ r^−1.96 and r^−2.14, respectively; Section 4.3.2) and monolithic at the resolution of our PdBI observations. With steep density profiles, estimated masses of the order of ∼20–25 M_⊙ (Section 4.2), and dark at 8 and 24 μm (note that leaf #7 has a 70 μm counterpart; Nguyen Luong et al. 2011), these are currently the best candidates for progenitors of intermediate-to-high mass stars within the mapped region.

SUMMARY AND CONCLUSIONS

We use high angular resolution (∼4 arcsec; 0.05 pc) 3.2 mm PdBI continuum observations to perform a structural analysis of the filamentary IRDC G035.39–00.33. To date, these are the highest angular resolution continuum observations of G035.39–00.33, surpassing previous observations by factors of ∼2–3 (i.e. the 70 μm Herschel data presented by Nguyen Luong et al. 2011 and the 1.2 mm observations presented by Rathborne et al. 2006). Our analysis leads us to conclude the following.

• The continuum emission is highly structured. It is segmented into a series of 13 quasi-regularly spaced (λ_obs ∼ 0.18 pc) cores, identified as ‘leaves’ in the dendrogram analysis, that follow the major axis of the G035.39–00.33.

• Comparison between continuum and N₂H⁺ (1−0) observations suggests that some of the identified leaves may reflect a superposition of structures associated with different velocity components. Although the translation between position–position–velocity and true three-dimensional space can be uncertain, this result emphasises the importance of exercising caution when attempting to classify structure in two-dimensional maps.

• Some leaves which appear to be kinematically coherent (i.e. they can be attributed to a single velocity component) can also exhibit structured continuum emission, which is evident in their radial flux density profiles. However, the scales at which this substructure resides is beyond the angular resolution of the observations in this study and further investigation is required.

• There is a significant (a factor of ∼8) discrepancy between the spatial separation of the leaves and that predicted by theoretical work describing the fragmentation of purely hydrodynamic fluid cylinders. Consistent with the kinematic analysis of Henshaw et al. (2014), who find evidence for the presence of sub-filaments observed throughout G035.39–00.33, this result emphasizes the importance of considering the underlying physical structure (and potentially, dynamically important magnetic fields) in any fragmentation analysis.

• The leaves exhibit a range in column density (3.6 × 10²³ cm⁻² < N_{H, c} < 8.0 × 10²³ cm⁻²), mass (8.1 M_⊙ < M_c < 26.1 M_⊙), and number density (6.1 × 10⁵ cm⁻³ < n_{H, c, eq} < 14.7 × 10⁵ cm⁻³).

• We used the derived physical properties of the leaves to assess their dynamical state, and determine the likelihood that they will undergo gravitational collapse. All dendrogram leaves are consistent with being either sub-virial or approximately virial (α_vir ≲ 1, within the 60 per cent uncertainty, for κ_ρ = 1.5–2). Absolute values span a range 0.2 < α_vir < 1.3. In the absence of additional support, possibly from magnetic fields with strengths of the order of 230 μG < B < 670 μG, leaves that are strongly sub-virial are susceptible to gravitational collapse, and possibly further fragmentation. Leaves #7 and #9 are consistent with this picture (α_vir ∼ 0.4 and ∼0.2, respectively).

• The formation of sub-filaments, followed by the formation of cores native to those sub-filaments, and the possibility of further fragmentation may imply a multilayered fragmentation process within G035.39–00.33.

• Additional fragmentation may explain the presence of multiple peaks observed in the radial flux density profiles of several of the kinematically coherent leaves (i.e. those continuum sources that can be attributed to a single velocity component). In contrast however, leaves #7 and #9, dark in the mid-infrared, centrally concentrated (ρ ∝ r^−1.96 and r^−2.14, respectively), monolithic (with no discernible substructure at our PdBI resolution), and with estimated masses of the order of ∼20–25 M_⊙, are good candidates for progenitors of intermediate-to-high mass stars.

Looking towards future investigations, higher angular resolution dust continuum observations will assist in determining whether or not the structures identified in this work have fragmented further. Similarly, constraining the strength and orientation of the magnetic fields, as well as searching for infall motions through molecular line observations, will help to assess whether cores such as these deviate from virial equilibrium.

We would like to thank the anonymous referee for the constructive report which has helped to improve the paper. We would like to thank Michael Butler and Jouni Kainulainen for providing the mass surface density map used in this work. Based on observations carried out under project number V008 with the IRAM PdBI. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). JDH would like to thank Gary Fuller and Tom Hartquist for their constructive comments on an early version of this project, as well as Nate Bastian, Yanett Contreras, Luke Maud, Fumitaka Nakamura, and Andy Pon for helpful discussions. PC and JEP acknowledge support from European Research Council (ERC; project PALs 320620). IJ-S acknowledges the funding received from the STFC through an Ernest Rutherford Fellowship (proposal number ST/L004801/1). JCT acknowledges NASA grant ADAP10-0110. RJP acknowledges support from the Royal Astronomical Society in the form of a research fellowship.

REFERENCES