What is controlling the fragmentation process in the Infrared Dark Cloud G14.225-0.506? Differet level of fragmentation in twin hubs

We present observations of the 1.3 mm continuum emission toward hub-N and hub-S of the infrared dark cloud G14.225-0.506 carried out with the Submillimeter Array, together with observations of the dust emission at 870 and 350 microns obtained with APEX and CSO telescopes. The large scale dust emission of both hubs consists of a single peaked clump elongated in the direction of the associated filament. At small scales, the SMA images reveal that both hubs fragment into several dust condensations. The fragmentation level was assessed under the same conditions and we found that hub-N presents 4 fragments while hub-S is more fragmented, with 13 fragments identified. We studied the density structure by means of a simultaneous fit of the radial intensity profile at 870 and 350 microns and the spectral energy distribution adopting a Plummer-like function to describe the density structure. The parameters inferred from the model are remarkably similar in both hubs, suggesting that density structure could not be responsible in determining the fragmentation level. We estimated several physical parameters such as the level of turbulence and the magnetic field strength, and we found no significant differences between these hubs. The Jeans analysis indicates that the observed fragmentation is more consistent with thermal Jeans fragmentation compared with a scenario that turbulent support is included. The lower fragmentation level observed in hub-N could be explained in terms of stronger UV radiation effects from a nearby HII region, evolutionary effects, and/or stronger magnetic fields at small scales, a scenario that should be further investigated.


INTRODUCTION
One ubiquitous fact in the process of formation of intermediate/high-mass stars is that they are usually found associated with clusters of lower-mass stars (e. g., Pudritz 2002;Lada & Lada 2003). However, the process behind the formation of such clusters is unclear as it remains ambiguous how a cloud core fragments to finally form a cluster. A common approach is to assume that the fragmentation is controlled by gravitational instability, where the velocity dispersion is, at best, accounted for only as an extra source of pressure support. In massive star-forming regions several studies show that thermal fragmentation alone does not account for the observed masses and/or core (fragment) separation (e. g., Zhang et al. 2009;Bontemps et al. 2010;Pillai et al. 2011;Wang et al. 2011;Naranjo-Romero et al. 2012;van Kempen et al. 2012;Wang et al. 2014;Zhang et al. 2015;Liu et al. 2015). Most of these works only address the formation of massive dense cores, which are far more massive than the thermal Jeans mass, and thus do not follow thermal Jeans fragmentation. However, in stellar clusters, most of stars are low-mass stars and their masses are related to thermal Jeans mass, as found by Takahashi et al. (2013), Palau et al. (2015), and Teixeira et al. (2015), raising the question of what is controlling the fragmentation process? The observational results suggest that cloud fragmentation leading to the formation of pre-stellar cores must be controlled by a complex interaction of gravitational instability, turbulence, magnetic fields, cloud rotation, and stellar feedback (e. g., Padoan & Nordlund 2002;Hosking & Whitworth 2004;Machida et al. 2005;Girart et al. 2013).
The infrared dark cloud (IRDC) G14.225-0.506 (hereafter G14.2), also known as M17 SWex (Povich & Whitney 2010), is part of an extended (77 pc×16 pc) and massive (>10 5 M ⊙ ) molecular cloud discovered by Elmegreen & Lada (1976) and located southwest of the Galactic HII region M17. The distance to M17 has been recently determined through trigonometric parallaxes of CH 3 OH masers, to be 1.98 +0.14 −0.12 kpc (Xu et al. 2011;Wu et al. 2014). Similar local standard of rest velocity (Elmegreen et al. 1979;Busquet et al. 2013) suggests that the bright HII region M17 and the IRDC are located at the same distance. The cloud is associated with star formation activity as revealed by the presence of several H 2 O masers (Jaffe et al. 1981;Palagi et al. 1993;Wang et al. 2006) and intermediate-mass YSOs (Povich & Whitney 2010). Through the analysis of near-/mid-infrared photometry Spitzer data, Povich & Whitney (2010) report an absence of early O stars and suggest a delay in the onset of massive star formation in the region, concluding that G14.2 represents a proto-OB association that has not yet formed its most massive stars.
Recent high angular resolution observations of the NH 3 dense gas, obtained by combining VLA and Effelsberg 100 m data (Busquet et al. 2013), reveal a network of filaments constituting two hub-filament systems. Hubs are associated with H 2 O maser emission (Wang et al. 2006) and mid-infrared sources (Povich & Whitney 2010). They appear more compact, are warmer (T rot ∼15 K) and show larger velocity dispersion and larger masses per unit length than filaments (Busquet et al. 2013), suggesting that they are the main sites of stellar activity within the cloud. In Figure 1 we provide a general perspective of the cloud, showing the dense filaments and hubs and their association with IR sources and/or water masers.
In this paper we present 1.3 mm observations performed with the Submillimeter Array (SMA 11 ), 870 µm observations carried out with the LABOCA bolometer at the APEX telescope, and 350 µm observations with the SHARC II bolometer at the Caltech Submillimeter Observatory (CSO), toward the two hubs, dubbed hub-N and hub-S (see Fig. 1). These observations are aimed at studying the cloud fragmentation process that leads to the formation of dense core fragments. Section 2 describes the observations and data reduction process. In Section 3 we present the results of the dust continuum emission observations, and their analysis is presented in Section 4. We discuss the interplay between density, turbulence, magnetic fields, UV radiation feedback, and evolutionary effects in the fragmentation process in Section 5. Finally, in Section 6 we list the main conclusions.

OBSERVATIONS AND DATA REDUCTION 2.1. Submillimeter Array Observations
The observations were carried out with the SMA (Ho et al. 2004) between 2008 and 2015 in different array configurations. Table 1 lists the observation dates, the configuration used, the number of antennas, and the projected baselines for each observing date. We performed a mosaic of 3 pointings toward hub-S and of 2 pointings toward hub-N, covering an area of ∼ 2 ′ × 1 ′ and ∼ 1 ′ × 1. ′ 5, respectively (see blue contour in Fig. 1).
For hub-S, the SMA correlator covered a 2 GHz bandwidth in each of the two sidebands, which are separated by 10 GHz. Each sideband is divided into 6 blocks, and each block consists of 4 chunks, yielding 24 chunks of 104 MHz width per sideband. For both compact and extended configurations the correlator was set to the standard mode with 256 channels per chunk, providing uniform channel spacing of 0.406 MHz (or 0.52 km s −1 ) per channel across the full bandwidth of 2 GHz. On the other hand, for hub-N the observations were conducted as filler projects in different runs using the 4 GHz receivers. In this case, each sideband was divided into 12 blocks with 4 chunks of 104 MHz each. The configuration of the correlator was different for each observing run. Table 1 summarizes the details of the SMA correlator, listing the local oscillator (LO) center frequency and the total bandwidth used to build the continuum.
The zenith opacity (τ 225 GHz ), measured with the National Radio Astronomy Observatory (NRAO) tipping radiometer 11 The Submillimeter Array is a joint project between the Smithsonian Astrophysical Observatory and the Academia Sinica Institute of Astronomy and Astrophysics, and is funded by the Smithsonian Institution and the Academia Sinica. Fig. 1.-Overview of the large scale structure of the G14.2 complex. White contours represent the 3σ contour level (27 mJy beam −1 km s −1 ) of the NH 3 (1,1) integrated intensity map (Busquet et al. 2013). The NH 3 synthesized beam, 8. ′′ 0×7. ′′ 2 (P. A.= −15 • ) is shown in the bottom left corner. Color scale is the 8 µm Spitzer image. Red stars indicate IRAS sources in the field, and pink crosses mark the position of H 2 O masers (Wang et al. 2006). Color stars depict the positions of the YSOs identified by Povich & Whitney (2010) where colors and symbols indicate the YSO evolutionary stage; three-point red star: Stage 0/I (dominated by infalling envelope); four-point yellow star: Stage II (optically thick circumstellar disk); five-point green star: Stage III (optically thin disk); six-point cyan star: ambiguous (see Povich et al. 2009;Povich & Whitney 2010 for further details on the YSOs classification). The SMA field of view of the two regions mosaiced, hub-N and hub-S, are indicated with the blue contour. located at the Caltech Submillimeter Observatory, was stable during the observations in the compact configurations, with values ∼ 0.1-0.2, and around 0.1 during the observations in the extended and very extended configurations.
The visibility data were calibrated using the standard procedures of MIRIAD (Sault et al. 1995). Information on the calibrators used is given in Table 1. The flux calibration is estimated to be accurate within ∼ 20 %.
The continuum was constructed from the line-free channels in the visibility domain. To produce the CLEANed single images of the two hubs, we used the tasks "invert" (option mosaic), "mossdi", and "restor" as recommended in the MIRIAD users guide 12 (Sault et al. 1995). The final combined images were obtained by weighting the visibilities in inverse proportion to the noise variance. We summarize the basic parameters (synthesized beam and rms noise level) of the different maps in Table 2. In this work we focus on the dust continuum emission, leaving the analysis of the molecular line emission for a   (Palau et al. 2010). b Root-mean-square noise level of the continuum data. c FWHM: full-width at half-maximum; P. A.: position angle. forthcoming paper.

APEX Observations
Continuum observations at 870 µm were carried out using LABOCA bolometer array, installed on the Atacama Pathfinder EXperiment (APEX 13 ) telescope. The array consists of 259 channels, which are arranged in 9 concentric hexagons around the central channel. The field of view of the array is 11. ′ 4, and the angular resolution of each beam is 18. ′′ 6 ± 1 ′′ .
The data were acquired on 2008 August 24 and 31 during the ESO program 081.C-0880A, under excellent weather conditions (zenith opacity values ranged from 0.15 to 0.24 at 870 µm). Observations were performed using a spiral raster mapping. This observing mode consists of a set of spirals with radii between 2 ′ and 3 ′ at a combination of 9 and 4 raster positions separated by 60 ′′ in azimuth and elevation, with an integration time of 40 seconds per spiral. This mode provides a fully sampled and homogeneously covered map in an area of 15 ′ ×15 ′ . The final map consisted of a mosaic of two points centered at α(J2000)=18 h 18 m 17 s .5, δ(J2000)=16 • 44 ′ 00. ′′ 0, and α(J2000)=18 h 18 m 17 s .5, δ(J2000)=16 • 57 ′ 00. ′′ 0, respectively. The total on-source integration time was ∼ 2 hours per position. Calibration was performed using observations of Mars as well as secondary calibrators. The absolute flux calibration uncertainty is estimated to be ∼ 8 %. The telescope pointing was checked every hour, finding a rms pointing accuracy of 2 ′′ . Focus settings were checked once per night and during the sunset.
We reduced the data using MiniCRUSH software package (see Kovács 2008). The pre-processing steps consisted of flagging dead or cross-talk channels frames with too high telescope accelerations and with unsuitable mapping speed, as well as temperature drift correction using two blind bolometers. Data reduction process included flat-fielding, opacity correction, calibration, correlated noise removal (atmospheric fluctuations seen by the whole array, as well as electronic noise originated in groups of detector channels), and despiking. Every scan was visually inspected to identify and discard corrupted data. The final map was smoothed to a final angular resolution of ∼ 22 ′′ , and the rms noise level achieved was 25 mJy beam −1 . In this paper we focus our attention on the dust emission associated with the two hubs, hub-N and hub-S, identified using high angular resolution NH 3 data (Busquet et al. 2013). The analysis of the entire cloud will be the subject of a forthcoming paper.

CSO SHARC II 350 µm
High angular resolution continuum observations at 350 µm were carried out using the SHARC II bolometer array, installed on the Caltech Submillimeter Observatory (CSO 14 ). The array consists of 12×32 pixels (approximately 85 % of these pixels work well). The simultaneous field of view provided by this array is 2. ′ 59 × 0. ′ 97, and the diffraction limited beam size is ∼ 9 ′′ .
The data were acquired on 2014 March 26 (τ 225 GHz ≃0.08). The telescope pointing and focusing were checked every 1.5-2.5 hours. Mars was observed for the absolute flux calibration, with a flux calibration uncertainty of ∼ 20 %. We used the standard 10 ′ × 10 ′ on-thefly (OTF) box scanning pattern, centered on the two positions α(J2000)=18 h 18 m 13 s .99; δ(J2000)=−16 • 51 ′ 00. ′′ 40, and α(J2000)=18 h 18 m 13 s .997; δ(J2000)=−16 • 59 ′ 00. ′′ 40. The final map covered an area of 20 ′ × 10 ′ . The total on-source time was 30 minutes for each of these two pointings. Data calibration was performed using the CRUSH software package (Kovács 2008). The final map was smoothed to an angular resolution of 9. ′′ 6 and the rms noise level achieved was ∼ 85 mJy beam −1 .   Figure 2 (left panel) and Fig. 3 (top panel) show the dust emission at 870 µm (blue dashed contours) and 350 µm (grey scale) obtained with the APEX and the CSO telescopes toward hub-N and hub-S, respectively. The smaller scale structure of the dust continuum emission follows the larger scale filaments seen in extinction, i. e., hub-N appears elongated in the north-south direction. Similarly, the large scale dust continuum emission of hub-S appears elongated along the eastwest direction. The morphology of the dust emission resembles that of the dense gas traced by NH 3 (Busquet et al. 2013). Table 3 summarizes the main physical parameters of the dusty envelope of both hubs, listing their peak position, peak intensity, and flux density. The deconvolved size was obtained by fitting a 2-dimensional Gaussian function. The two hubs have the same size, 0.6 × 0.3 pc at 870 µm and 0.3 × 0.2 pc at 350 µm. Note that the peak positions at 870 and 350 µm are offset from each other by ∼ 6 ′′ in hub-N and ∼ 3 ′′ in hub-S. This is most likely due to pointing errors. Further details on the physical properties of hubs are presented in Section 4.1 where we model the radial intensity profile and the spectral energy distribution of each hub.

Dust emission at 0.03 pc scale
In this section we present first the results obtained with the SMA in the compact configuration with an angular resolution of ∼ 3 ′′ − 4 ′′ , corresponding to ∼ 0.03 pc at the distance of the source. Then, in Sect. 3.3 we show the images of each hub obtained by combining all configurations. This allows us to obtain high angular resolution maps (∼ 1. ′′ 5) at the best possible sensitivity. From these maps, we used a 6σ threshold to identify the millimeter condensations in each hub, with σ being the rms noise level and requesting that the 6σ contour level is closed. Figure 2 (left panel) and Fig. 3 (top panel) present the SMA 1.3 mm continuum emission (in white contours) obtained with the compact configuration toward hub-N and hub-S, respectively. At a similar angular resolution and similar sensitivity, 2-3 mJy beam −1 , the dust envelopes of both hubs appear fragmented into dust cores. We detected 3 and 7 compact continuum condensations above a 6σ contour level toward hub-N and hub-S, respectively, in the compact configuration.
In hub-N, most of the 1.3 mm emission arise from a bright source, MM1, located at the peak of the single-dish submillimeter dust emission. This appears to be elongated along the north-south direction with a position angle of 7 • , similar to the large scale 870 µm emission (see Table 3). There are two faint and more flattened millimeter condensations associated with hub-N: MM2 located about 46 ′′ south of MM1, and MM3 located ∼ 12 ′′ south-east of MM1. All millimeter condensations in hub-N are deeply embedded within the dust envelope and the NH 3 (1,1) dense gas emission. , which were obtained using all available SMA configurations. The sensitivity of the final SMA images is the same: 1 mJy beam −1 . At a 6σ level, this sensitivity corresponds to a mass sensitivity of 0.7 M ⊙ , assuming a dust temperature of 17 K (see below), a gas-to-dust mass ratio of 100, and a dust mass opacity coefficient at 1.3 mm per unit mass density of dust and unit length of 0.899 cm −2 g −1 , which corresponds to coagulated grains with thin ice mantles in cores of densities ∼ 10 6 cm −3 (Ossenkopf & Henning 1994). We then identified sources having peak fluxes above the 6σ rms noise level. The number of millimeter continuum condensations in hub-N is 9 while in hub-S the number of millimeter condensations is 17. In Table 4 we report the identified condensations for each hub, listing their position, peak flux, total flux density, and deconvolved size obtained for a 2-dimensional Gaussian fitting, when possible. For those sources that cannot be fitted with a Gaussian function we estimated the flux density by integrating the emission above the 4σ level. We additionally report the association of a dust millimeter source with a near-IR source, and the presence of infrared excess at 4.5 µm according to the catalog of Povich & Whitney (2010), as well as its association with H 2 O maser emission (Wang et al. 2006).
As can be seen in the close-up image of hub-N ( Fig. 2-right panel) the high angular resolution map shows that in addition to MM2 and MM3, the strongest millimeter continuum source MM1, which has a peak intensity of 56 mJy beam −1 , breaks up in 6 additional sources. The strongest source MM1a is associated with a H 2 O maser spot (Wang et al. 2006). There are two sources located north of MM1a, labelled as MM1e and MM1f, and three faint sources south of MM1a, labelled as MM1b, MM1c, and MM1d. MM1f is associated with a YSO, identified in the Spitzer image by Povich & Whitney (2010), in Stage 0/I of evolution, which means it is dominated by an infalling envelope. Moreover, this source presents a 4.5 µm excess (Povich & Whitney 2010). Sources with this feature are typically called "Extended Green Objects" (EGOs) and is due to protostellar outflow activity, indicating that MM1f could be a massive YSO candidate (Cyganowski et al. 2008). We detected an additional source, called MM1g, which lies about 5. ′′ 6 south east of MM1a whose emission is very com-  Table 1. Blue dotted contour corresponds to the −3σ contour level. The dashed contours at 870 µm start at 3σ rms level and increase in steps of 20σ, where σ is 25 mJy beam −1 . The CSO 350 µm grey-scale image start at 2 % of the maximum, with levels increasing in 10 % of the maximum value. The synthesized beam of the SMA image is shown in the bottom left corner of the image. Right: Close-up of the MM1 and MM3 millimeter sources that show the SMA 1.3 mm continuum emission obtained by combining all configurations. The red thick contour level denotes the 6 σ level used to identify the continuum sources. Contours range from 3 to 21σ in steps of 3σ, and from 21 to 51σ in steps of 10σ, where σ is the rms of the map listed in Table 1. Dashed contours represent the NH 3 (1,1) integrated intensity emission (Busquet et al. 2013)  pact and faint. While MM3 is still detected in the high-angular resolution map, the emission of MM2, detected only at a 4σ level, is filtered out due to a combination of a faint and flattened structure.
Regarding the high angular resolution 1.3 mm image of hub-S, Fig. 3 (bottom panel) shows that the snake-shape structure located at the center of the image, which consisted in 4 fragments in the compact configuration, breaks up into 7 condensations. The strongest and compact source MM5a is associated with a H 2 O maser (Wang et al. 2006) and MM5b has its peak position very close to a YSO in the Stage 0/I of evolution that presents an excess in the 4.5 µm Spitzer band. MM1a and MM1b lie at the border of the dust envelope, about 15 ′′ to the west of MM5a. The brightest source, MM1a, is associated with an infrared source classified by Povich & Whitney (2010) as a YSO in the Stage 0/I of evolution (i. e., dominated by an infalling envelope). MM6 splits in two sources separated 2. ′′ 2, and MM7 breaks up in three sources, MM7a is the strongest one, and MM7b and MM7c are 3 ′′ and 6 ′′ to the south, respectively. The higher sensitivity of the combined map allowed to detect two additional sources, MM8 and MM9, with a high-enough signal-to-noise ratio.
Assuming that the dust emission at 1.3 mm is optically thin and that the temperature distribution is uniform, we estimated the mass of each millimeter condensation adopting a gas-to-dust ratio of 100 and using the dust mass opacity coefficient at 1.3 mm per unit mass density of dust and unit length of 0.899 cm 2 g −1 (Ossenkopf & Henning 1994). In both hubs, the rotational temperature obtained from NH 3 observations is 15 K (Busquet et al. 2013). We converted T rot to gas kinetic temperature using the expression provided by Tafalla et al. (2004), and obtained T kin ≃17 K. Assuming T kin =T d the masses of the continuum sources are estimated to be in the range 1 − 13 M ⊙ in hub-N and 0.7 − 18 M ⊙ in hub-S. We report the mass of each millimeter condensation in Table 4. Due to the uncertainty in the dust mass opacity coefficient, the values of the derived masses are good to within a factor of 2. While in hub-N only MM1a has a mass larger than 10 M ⊙ , in hub-S there are four sources with M > 10 M ⊙ . The sum of all mass estimated from the 1.3 mm dust continuum emission is 34.2 M ⊙ and 95.1 M ⊙ in hub-N and hub-S, respectively. It is worth noting that the mass of all millime-  Table 1. Bottom: Contour map of the SMA dust emission at 1.3 mm toward G14.2-hub-S obtained by combining the compact and extended configuration observations (angular resolution ≃ 1. ′′ 5) overlaid on the dust emission at 350 µm from the CSO SHARC II bolometer (grey scale) and the NH 3 (1,1) integrated intensity map (dashed contours) from Busquet et al. (2013). Contours at 1.3 mm range from 3 to 30σ increasing in steps of 3σ, and from 30 to 60σ in steps of 10σ, where σ is the rms of the map listed in Table 1. Blue dotted contour corresponds to the −3σ contour level. The red thick contour level denotes the 6σ level used to identify the continuum sources. The SMA synthesized beam is shown in the bottom left corner of each image. Symbols are the same as in Fig. 1. ter condensations, except for MM2 in hub-N, have been measured using the combined image.

Fragmentation level
In the previous sections we have shown the SMA maps at different resolutions. Now, the aim is to estimate the fragmentation level in each hub. To do so, we need first to obtain images with the same angular resolution to be easily comparable. We then performed the imaging of each hub using a common uv-range of visibilities (8-160 kλ). The observations of hub-N have more visibilities at longer uv-distances than hub-S, which results in a slightly smaller beam, especially in one direction (roughly the north-south direction), so we convolved the map to the beam obtained for hub-S (i. e., 1. ′′ 53 × 1. ′′ 41). In addition, since our SMA observations consist of a small mosaic instead of a single pointing, and the total field of view of each hub is different, we need to define a region in which we can evaluate in a consistent way the fragmentation level in each hub. The typical radii of compact and embedded clusters and subclusters lie in the range 0.1 − 0.2 pc (e. g., Testi et al. 1999;Alexander & Kobulnicky 2012;Kuhn et al. 2014). Then, we defined a radius of 0.15 pc, which corresponds to a region of ∼ 30 ′′ at the distance of the source. In Fig. 4 we present the 1.3 mm continuum emission in a field of view of 30 ′′ for hub-N (left panel) and hub-S (right panel). Adopting this field of view and the 6σ detection  (Wang et al. 2006). d Properties of MM2 determined using the SMA image at 4 ′′ resolution. e A 2D gaussian cannot be fitted and the peak intensity and total flux density were obtained by integrating the emission within the 4σ contour level. threshold, we found 4 fragments in hub-N and 13 fragments in hub-S. These values are listed in Table 6. Therefore, using the proper comparison we still clearly see that hub-S has a higher fragmentation level than hub-N.

ANALYSIS 4.1. Radial Intensity Profiles
Recently, Palau et al. (2014) studied the relation between the fragmentation level and the density structure in a sample of 19 massive dense cores selected to be in a similar evolutionary stage, i. e., all regions harbor intermediate/high-mass protostars deeply embedded in massive dense cores that have not yet developed an ultra-compact HII region, hence comparable to our targets hub-N and hub-S. They find a weak (inverse) trend of fragmentation level and density power-law index, with steeper density profiles tending to show lower fragmentation and vice-versa. One of the main results of this work is that, within a given radius, the fragmentation level increases with average density as a combination of flat density profile and high central density. A comparison with magnetohydrodynamic simulations (Commerçon et al. 2011) suggests that the cores showing no fragmentation could be related with a strong magnetic field. In a follow-up work using the same sample of massive dense cores, Palau et al. (2015) analyze whether the fragmentation level can be explained by turbulent or thermal support, and conclude that the observed fragmentation level is more consistent with pure thermal Jeans fragmentation than fragmentation including turbulent support.
In order to investigate whether the differences of the fragmentation level in hub-N and hub-S could be explained in terms of the density structure of the envelope, we applied the model presented in Palau et al. (2014) to fit simultaneously the radial intensity profile at 870 µm and 350 µm and the spectral energy distribution (SED) assuming spherical symmetry. It is important to mention that both hubs are embedded in filamentary structures, and thus the radial profile along the filament's main axis has contributions from the filament and the hub. Since we are interested in modelling the structure of the hub, not the filament, we extracted the radial profiles in the direction perpendicular to the filament main axis to avoid contamination from the filament itself at large distances. The radial intensity profiles at 870 µm and 350 µm were obtained in rings of 10 ′′ and 6 ′′ width, respectively, as a function of the projected distance from the clump center. The radial intensity profiles of hub-N and hub-S at 870 µm and 350 µm are presented in Fig. 5.
The SED of each hub was built considering the flux density at 350 and 870 µm from this work, and the flux densities measured by the ESA Herschel Space Observatory (Pilbratt et al. In the present model, we describe the density structure using a Plummer-like function of the form where ρ c is the central density, r c is the radius of the flat inner region, and p is the asymptotic power index. The temperature was assumed to be a power-law of radius where q = 2/(4 + β). We assumed a dust opacity law as a power-law of frequency with index β, κ ν = κ 0 (ν/ν 0 ) β , where ν 0 is an arbitrary reference frequency. At ν 0 = 230 GHz, we used the dust mass opacity coefficient adopted in Sect. 3.3. We adopted a constant temperature of 8 K for radii larger than the radius where the temperature distribution of index −q drops to 8 K. We stress that a minimum temperature of 10 K does not alter the results. Concerning the density distribution, we adopted, as a maximum radius of the envelope, the radius for which the envelope density achieved a value comparable to that of the ambient gas of the intecore medium, taken to be 5 × 10 3 cm −3 , the same value adopted in Palau et al. (2014).
The model has five free parameters: the dust emissivity index β, the envelope temperature T 0 at the reference radius r 0 , the envelope density ρ 0 extrapolated to the reference radius r 0 , the radius of the inner region where density is flat (in units of r 0 ) r c /r 0 , and the density power law index p. The reference radius r 0 is arbitrary and has been taken as r 0 = 1000 AU. It is important to remark that ρ 0 is not the actual density at r 0 , but the value given by the asymptotic density power law, ρ 0 = ρ c (r 0 /r c ) −p . In fact, if r 0 < r c , the density ρ 0 is higher than the maximum density ρ c . Only in the case r 0 ≫ r c does the density ρ 0 correspond to the actual density at r = r 0 . We used this model to fit simultaneously the observed radial intensity profiles at 870 µm and 350 µm and the observed SED. We adopted an uncertainty in the radial intensity profiles at both 870 and 350 µm of 20 %. The uncertainty of the flux densities used in the SED was 20 % at all wavelengths, plus an additional 20 %, added quadratically, for the 70 and 160 µm values (the PACS data). Note that SHARC II and LABOCA instruments cannot detect structures larger than 1 ′ scale.
For each set of model parameters we computed the intensity map at 870 µm and 350 µm and convolved it with a Gaussian of the size of the beam of the telescope used in the observation. The intensity profile was obtained from the convolved map, and compared with the circularly averaged observed radial intensity profiles. On the other hand, in order to compare the model with the observed SED, we computed the model flux density by integration of the model intensity profile within the aperture used to estimate the flux for each data point.  Table 5). Note that SHARC II and LABOCA cannot detect structures > 1 ′ scale.
The initial search ranges for the 5 parameters were β = 1.5 ± 1.5, T 0 = 300 ± 300 K, ρ 0 = (4.0 ± 4.0) × 10 −16 g cm −3 , r c /r 0 = 20 ± 20, and p = 1.5 ± 1. The search was carried out for 13 loops, each loop consisting in 4000 samples of the parameter space, with a search range reduced a factor 0.75 around the best-fit value of the parameters found for the last loop. A second search was performed with the best values found in the first run as starting values of the fit parameters. For this second run the initial search ranges were 0.2, 20 K, 2.0 × 10 −16 g cm −3 , 8, and 0.15, for β, T 0 , ρ 0 , r c /r 0 , and p, respectively, and the search range was reduced a factor 0.8 for each loop. We refer to the work of Palau et al. (2014) for further details on the model and the fitting procedure used to find the best-fit model (see also Estalella et al. 2012;Sánchez-Monge et al. 2013, for details on the χ 2 minimization). Figure 5 shows the best fit to the 870 and 350 µm radial intensity profile and the SED for hub-N (top panel) and hub-S (bottom panel). The fitted parameters are reported in Table 5. Interestingly, despite the fragmentation level being significantly higher in hub-S than in hub-N, Table 5 shows that the parameters inferred from the model are remarkably similar. The dust emissivity index is ∼ 1.8-1.9, which corresponds to a temperature power law index of q ≃ 0.34, consistent with internal heating. The fitted temperature at 1000 AU is 51±2 K and 45±2 K in hub-N and hub-S, respectively. The central density, ρ c , is (1.0-1.3)×10 −18 g cm −3 and the radius of the flat inner region is r c ≃ 0.1 pc. In Table 5 we also show the radius where the density reaches 5 × 10 3 cm −3 and the radius where the temperature drops down to 8 K. In both hubs we obtained the same index p ≃ 2.2 of the asymptotic density power law.

General properties of the hubs
Once we obtained the best-fit model, we derived different physical quantities, which are reported in Table 6. We estimated the total mass of each hub by integrating the model intensity up to the radius where the density profile could be measured. Within the errors, the mass of both hubs is essentially the same, being slightly more massive hub-N than hub-S. Given that we evaluated the fragmentation level within a region of 0.15 pc in radius, we estimated the mass (M 0.15 pc ), average density (n 0.15 pc ), and average temperature (T 0.15 pc ) within a radius of 0.15 pc. As expected, the inferred values in both hubs are remarkably similar. Note that T 0.15 pc is 18 K and 16 K in hub-N and hub-S, respectively, which is consistent with the kinetic temperature, T kin =17 K, derived using the NH 3 data (Busquet et al. 2013).
Given that our aim is to investigate the interplay between the different physical agents acting during the fragmentation process, in the following we estimated some relevant physical quantities to evaluate the relative importance of the different forces within a region of 0.3 pc in diameter, the region where we assessed the fragmentation level in each hub.
Velocity dispersion. We extracted the averaged NH 3 (1,1) spectrum of both hubs and fitted the hyperfine structure using the CLASS package of the GILDAS software to measure the NH 3 (1,1) line width, ∆v obs , over a region of 0.3 pc of diameter (i. e., 30 ′′ ). Within this region, infall and rotation motions may contribute to the observed linewidth, which is 2.27 and 2.38 km s −1 in hub-N and hub-S, respectively. Thus, in order to obtain the contribution due to turbulent motions we estimated the velocity gradient within this region of 0.3 pc using the NH 3 data and found 3.8 km s −1 pc −1 , which can be attributed to rotation and/or infall motions (see the  TABLE 5 Best-fit parameters to the radial intensity profiles and SED (pc) (pc) hub-N 1.81 ± 0.08 51 ± 2 (1.3 ± 0.2) × 10 −15 21 ± 4 2.24 ± 0.04 0.69 (1.4 ± 0.6) × 10 −18 0.34 1.05 0.62 hub-S 1.89 ± 0.08 45 ± 2 (1.0 ± 0.2) × 10 −15 20 ± 3 2.24 ± 0.04 0.76 (1.2 ± 0.5) × 10 −18 0.34 0.79 0.57 a Free parameter fitted by the model: β is the dust emissivity index; T 0 and ρ 0 are the temperature and density at the reference radius r 0 = 1000 AU; r c /r 0 is the radius of the flat inner region in units of r 0 ; p is the asymptotic power index; χ r is the reduced chi-square, [χ 2 min /(n − 5)] 1/2 . b Parameters derived from the modeling: ρ c is the central density, ρ c = ρ 0 (r c /r 0 ) −p ; q is the temperature power-law index; r 8K is the radius of the clump where the temperature drops to 8 K; r max is the radius at the assumed ambient density of 5 × 10 3 g cm −3 .    Palau et al. 2015); M tot Jeans is the total Jeans mass using the total (thermal + non-thermal) velocity dispersion (Eq. 4 of Palau et al. 2015). e σ 1D,nth is the non-thermal turbulent component of the velocity dispersion estimated as σ 1D,nth = σ 2 obs − σ 2 th using the NH 3 (1,1) data.
f M is the Mach number defined as σ 3D,nth /c s where σ 3D,nth = √ 3 σ 1D,nth ; M A is the Alfvén Mach number, expressed as cr is the mass-to-magnetic-flux ratio; and β rot is the rotational-to-gravitational energy ratio as defined in Chen et al. (2012) and derived using NH 3 (1,1) data. derivation of the rotational-to-gravitational energy ratio for further details) and subtract, in quadrature, this contribution from the observed line width. The resulting values are 1.96 and 2.09 km s −1 in hub-N and hub-S, respectively. We subsequently obtained the observed velocity dispersion, σ 1D,obs , and separate the thermal (σ 1D,th ) from the non-thermal turbulent (σ 1D,nth ) contribution of the observed velocity dispersion following the procedure described in Palau et al. (2015) and using the average kinetic temperature within a radius of 0.15 pc, which has been obtained from the model (see Table 6). Once we derived the thermal contribution of the velocity dispersion, σ 1D,th = κ B T/(µ m H ) with µ = 17 for NH 3 , which is ∼ 0.1 km s −1 , we obtained the non-thermal turbulent component, which is listed in Table 6. The isothermal sound speed at the temperature of the two hubs is c s ≃ 0.24 km s −1 and the total velocity dispersion, defined as σ 1D,tot = c 2 s + σ 2 1D,nth , is ∼ 0.9 km s −1 .

Mach number.
We computed the Mach number M as in Palau et al. (2015) and obtained a value of 5.6-6.4, which means that hubs show supersonic non-thermal gas motions.
Jeans mass and Jeans number. We investigated whether the fragmentation can be explained in terms of pure thermal Jeans fragmentation and including both thermal and non-thermal support adopting Equations (6) and (4) of Palau et al. (2015), respectively. Using the average temperature and density within a radius of 0.15 pc (listed in Table 6), we derived a thermal Jeans masses of 1.3 M ⊙ and 1.2 M ⊙ for hub-N and hub-S, respectively. The Jeans mass increases up to values of 70 M ⊙ and 90 M ⊙ if we include the non-thermal support.
The number of expected fragments, N Jeans , is given by the ratio between the mass of each hub inside a region of 0.15 pc in radius, M 0.15 pc , and the Jeans mass, M Jeans , scaled by a core formation efficiency (CFE) to take into account that not all the mass of the hub will be accreted and converted into compact fragments (equation (7) of Palau et al. 2015). Thus, to estimate the number of fragments we need first to compute the CFE of each hub, which can be calculated from the ratio of the sum of the masses of all fragments measured with the SMA, listed in Table 4, and the total mass of the hub obtained from the model and listed in column (4) of Table 6. The values of the CFE are 3 % and 13 % for hub-N and hub-S, respectively, consistent with the values reported in other works (e. g., Bontemps et al. 2010;Palau et al. 2013Palau et al. , 2015Louvet et al. 2014). Using the values of CFE, we obtained that the number of expected fragments considering pure thermal Jeans fragmentation is 3 and 12 fragments in hub-N and hub-S, respectively, while if we add turbulent support the expected number of fragments is smaller than unity. Therefore, the observed number of fragments in each hub is in very good agreement with pure thermal Jeans fragmentation.
Alfvén Mach number. To evaluate the importance of turbulence against the magnetic field we calculated the Alfvén Mach number, M A , defined as M A = √ 3σ 1D,nth /v A , where v A = B tot / 4πρ is the Alfvén speed and B tot is the total magnetic field. Recently, Santos et al. (in preparation) has measured the angular dispersion factor of the magnetic field and the sky-projected magnetic field strength through near-infrared polarimetric data, which trace the magnetic field in the diffuse medium surrounding all filaments and hubs previously identified in Busquet et al. (2013). For the diffuse gas, the magnetic field is 0.3 mG in hub-N and 0.5 mG in hub-S. If we now assume that we are in the magnetic field controlled contraction regime, B ∝ n 1/2 H 2 , the values of the magnetic field extrapolated to the densities reported in Table 6 are 1.0 mG in hub-N and 1.5 mG in hub-S, which yields to M A of 0.4 and 0.3 in hub-N and hub-S, respectively (see Table 6). Clearly M A < 1, indicating sub-Alfvénic conditions, even if we take into account the uncertainty in the derived values of the sky-projected magnetic field, of the order of 50 %. Thus, magnetic energy dominates over turbulence. Similar values are found by Pillai et al. (2015) towards two massive IRDCs using submillimeter polarization data, and by Franco et al. (2010) toward the Pipe nebula using optical data.

Mass-to-magnetic-flux ratio.
The relevance of magnetic field with respect to the gravitational force can be assessed by the mass-to-magnetic-flux ratio (M/Φ). This ratio is expressed in terms of a critical value as λ = (M/Φ)/(M/Φ) cr = 7.6 × 10 −21 [N H 2 /cm −2 ][B tot /µG] −1 (Crutcher et al. 2004). For the column density, N(H 2 ), we used the mass obtained from the model within a region of 0.3 pc of diameter. Assuming spherical symmetry N H 2 = M/(2m H π R 2 ), we obtained N(H 2 )= 1.1 × 10 23 , and 0.9×10 23 cm −2 in hub-N and hub-S, respectively. The resulting mass-to-flux ratio are 0.9 and 0.5 in hub-N and hub-S, respectively (see Table 6). Taking into account the uncertainty associated with the sky-projected magnetic field strength, the mass-to-flux ratio is close to unity, as expected, since there is ongoing star formation activity in these hubs, and hence part of the mass has already been accreted onto the protostars. Although global gravitational collapse is certainly expected in these hubs, magnetic field cannot be ignored.
Rotational-to-gravitational energy ratio. Finally, we computed the rotational-to-gravitational energy ratio, β rot = E rot /E grav , following Equation (1) of Chen et al. (2012). We used the velocity field map of NH 3 (1,1) shown in Busquet et al. (2013) and estimated the velocity gradient within a region of 0.3 pc in diameter, finding similar values in both hubs, around 3.8 km s −1 pc −1 . Both rotation and infall motions contribute to the observed velocity gradient, and gravitational collapse is certainly expected in these two hubs as they contain deeply embedded protostellar objects. However, with the current data we cannot separate the contribution due to rotation from that due to infall motions. Assuming a mass infall rate of 10 −5 − 10 −3 M ⊙ /yr (Kirk et al. 2013;Peretto et al. 2013Peretto et al. , 2014 we derived an infall velocity, v infall =Ṁ/π R 2 ρ, lying in the range of 0.015-0.15 km s −1 . The lack of information from molecular line observations prevent us to estimate the velocity gradient due to infall or collapse motions, and hence, we derived β rot assuming that the observed velocity gradient is due only to rotation. The inferred values are the same in both hubs, β rot = 0.02, which is within the range of values observed in low-and high-mass star-forming regions (e. g., Goodman et al. 1993;Chen et al. 2012;Palau et al. 2013). This value should be regarded with caution as the velocity gradients used to estimate β rot was attributed as due to rotation motions only but it may be affected by the presence of several components and/or infall/outflow motions. Nevertheless, this value is slightly above the threshold value of 0.01 obtained in simulations of a rotating core that could fragment if the rotational energy is large enough compared to the gravitational energy (e. g., Boss 1999).

DISCUSSION
In the previous sections we showed that the two hubs associated with the IRDC G14.2 seem to be twin hubs. The model applied to obtain the underlying density and temperature structure of these hubs yield the same physical parameters, indicating a similar internal structure of the envelope where the fragments are embedded in, and all the derived quantities reported in Table 6 are notably similar as well. However, they present different level of fragmentation, hub-S being more fragmented than hub-N. With the aim of exploring whether these two similar hubs undergo similar fragmentation process, in this section we will discuss the interplay between density structure, turbulence, magnetic field, UV radiation feedback, and core evolution in the fragmentation process.
5.1. Density structure, turbulent fragmentation, and magnetic field One of the main results of Palau et al. (2014), who studied the relation between the fragmentation level and the density structure in a sample of 19 massive dense cores, is that there is a clear trend of fragmentation level increasing with average density within a given radius as a result of flat density profile and high central density (see Fig. 6 of Palau et al. 2014), indicating that density structure plays an important role in the fragmentation process. In Palau et al. (2014) the massive dense cores were observed with interferometers reaching spatial scales down to 1000 AU, picking up only the most compact fragments, and the region used to evaluate the fragmentation level was 0.05 pc in radius. Apparently, our results seem to be in contradiction with Palau et al. (2014) since hub-N and hub-S show the same density structure but they clearly present different fragmentation level, and according to Palau et al. (2014) we should expect higher central density and a steeper density profile in hub-N as it presents a lower fragmentation level. However, in the present work we assessed the fragmentation level in a region of 30 ′′ (i. e., 0.15 pc in radius) and the SMA combined data is sensitive to structures of 3000-10000 AU, and hence sensitive to flattened condensations. Therefore, to properly compare our result with the work of Palau et al. (2014) we evaluated the fragmentation level within the same field of view of Palau et al. (2014), i. e., 0.1 pc of diameter, and we obtained 3-4 fragments in each hub. From the results of the best-fit model obtained in Section 4.1 the density at 0.05 pc is 2.9×10 5 cm −3 and 2.4×10 5 cm −3 in hub-N and hub-S, respectively. Including these numbers in Fig. 6(c) of Palau et al. (2014) we can see that the observed number of fragments in hub-N and hub-S, N mm ≃ 3 − 4, is in good agreement with the trend observed by Palau et al. (2014), and the fragmentation level under these conditions would be the same. This result suggests that the fragmentation of most compact structures (which will probably form stars in the future) depends on the density structure while fragmentation including larger scales structures ( 5000 AU) might be determined by different processes. Thus, we conclude that the differences in the fragmentation level at the scales measured in this work cannot be produced by density and temperature structure of the hub, given that the model results are similar.
The Jeans analysis performed in Section 4.2 revealed that the fragmentation process can be explained without invoking turbulent support. In fact, as can be seen in Table 4, the most massive cores have masses > 10 M ⊙ which are closer to the turbulent Jeans mass rather than thermal Jeans mass. On the other hand, most cores have masses ∼ 1 M ⊙ and are close to thermal Jeans mass, which is consistent with the results of Zhang et al. (2015) and Palau et al. (2015) regarding the formation of massive dense cores and low-mass fragments, respectively. Both the Jeans mass of most cores and the expected number of fragments point towards thermal Jeans fragmentation since the inclusion of turbulence as an additional form of support implies the formation of 4-5 fragment at most and a CFE higher than 100 % (see Palau et al. 2015 for further discussion on CFE). Despite internal motions in both hubs belong to the supersonic regime, the level of turbulence at these scales (from 0.01 to 0.15 pc) does not seem to play an important role in the fragmentation of these hubs. In particular, the non-thermal velocity dispersion is ∼ 0.5 km s −1 in both hubs, suggesting that the effects of turbulence are the same. While turbulence at larger scales could be responsible for the cloud fragmentation into filamentary structures, at smaller scales gravity might be more relevant. In fact, Busquet et al. (2013) show that the fragmentation of filaments constituting the cloud complex IRDC G14.2 can be described by the fragmentation of an isothermal cylinder determined by thermal pressure with a (small) additional contribution from turbulent pressure. High angular resolution observations of molecular lines would provide information on the relative velocity between the different cores embedded in the hub, which would confirm the nature of the internal gas motions inside each core (Busquet et al. in prep.).
Another important physical agent that presumably determines the fragmentation level in a cloud core is the magnetic field. Observationally, Palau et al. (2013) investigate the effects of the magnetic field in a sample of dense cores and find that the low-fragmented regions are well reproduced in the magnetized core case, while the highly fragmented regions are consistent with cores where turbulence dominates over magnetic field as indicated by a comparison with magnetohydrodynamic simulations (Commerçon et al. 2011). Recent studies toward high-mass star-forming regions also suggest that magnetic field plays a crucial role during the collapse and fragmentation of massive clumps and the formation of dense cores (Zhang et al. 2015;Li et al. 2015). Regarding the IRDC G14.2, as shown in Busquet et al. (2013), the network of parallel filaments traced by the NH 3 dense gas could be the result of gravitational instability of a thin gas layer threaded by magnetic fields. In particular, magnetic fields seem to play an important role in the alignment of filaments as revealed by the preliminary near-infrared polarization measurements of a small part of the cloud. A follow up work of optical and near-infrared polarization data, covering the entire group of filaments, provides a measurement of the sky-projected field strength in each hub (Santos et al. in prep.). The inferred values are similar, ∼ 1 mG in hub-N compared to 1.5 mG in hub-S, and the derived Alfvén Mach numbers are smaller than 1, suggesting that during the fragmentation of the cloud core the magnetic field is an important ingredient that should not be ignored. Actually, Santos et al. show that magnetic fields are tightly perpendicular to the dense filaments and hubs but also to the cloud as a whole. However, polarization properties could be different when tracing different scales, i. e., different parts of the cloud (e. g., Alves et al. 2014), and hence the magnetic field deep into the dense hubs could be different from that inferred from the near-infrared polarimetric data. This might have the consequence that a stronger magnetic field pervades hub-N relative to hub-S resulting hub-N undergoing lesser fragmentation, a scenario that could be tested with ALMA.
Therefore, the physical properties investigated in this section, i. e., density structure, turbulence, and magnetic field, are remarkably similar in both hubs. We stress that we assumed that the magnetic field properties derived from optical and near-infrared data hold deep into the dense hubs, so any change in such properties, specially differences between the two hubs, could lead to different level of fragmentation in this region.

Evolutionary stage and UV radiation feedback effects
Another possible scenario that could explain the observed differences in the fragmentation level would be the effects of UV radiation from nearby high-mass stars and/or differences in the evolutionary stage.
As shown in Fig. 1, hub-N is located near the infrared source IRAS 18153-1651, which has a luminosity of ∼ 1.1 × 10 4 L ⊙ (Jaffe et al. 1982). The VLA 6 cm continuum emission reveal a cometary HII region ionized by a B1 star with the head of the cometary arch pointing toward hub-N (A. Sánchez-Monge 2012, private communication) that could compress and heat the gas. In fact, Busquet et al. (2013) find evidence for local heating in the western part of hub-N facing the HII region, and the position-velocity plot of the NH 3 dense gas along hub-N unveils an inverted C-like structure consistent with an expanding shell. Numerical simulations suggest that radiative heating from nearby high-mass protostars can strongly suppress fragmentation (e. g., Offner et al. 2009;Hansen et al. 2012;Myers et al. 2013). Radiative feedback from this IRAS source could thus explain why hub-N is less fragmented than hub-S since an increase in temperature would increase the Jeans mass and consequently the Jeans number would be smaller. However, both temperature and Jeans mass are extremely similar in these two hubs, suggesting that radiative feedback does not seem to play a dominant role in the fragmentation process of hub-N, and further observations are needed to investigate whether the UV radiation from the nearby IRAS source is affecting the properties of hub-N. Since molecular outflows are unknown in this region, we cannot account for their feedback effects.
Finally, the evolutionary stage of these hubs could also lead to the differences observed in the fragmentation level. We searched infrared sources classified as candidate to YSO by Povich & Whitney (2010) located within the SMA field of view. In Figs. 2 and 3, we show the position and the evolutionary stage (using different colors and symbols) of the YSOs identified in each hub. Hub-N contains 14 infrared sources whereas only 10 sources are associated with hub-S. The main difference arises in the number of infrared sources classified in the Stage 0/I (dominated by an infalling envelope): there are 10 YSOs in hub-N while hub-S contains only 6 YSOs in this stage of evolution. In both hubs, there are 3 YSOs in the more advanced evolutionary stage (Stage II, which have optically thick cicumstellar disk), and in each hub one source has been classified as ambiguous. Interestingly, the YSOs in hub-N appear to be more clustered around the bright millimeter condensation MM1a (see left panel in Fig. 4).
We report in Table 6 the ratio between the number of infrared sources without a millimeter counterpart and the number of millimeter sources within a region 0.3 pc of diameter. In this case, we do see differences between both hubs: N IR /N mm ≃ 1.7 in hub-N in front of 0.4 in hub-S. Thus, the relative number of infrared sources in front of the millimeter sources is larger in hub-N by a factor ∼ 4. Most of the objects associated with hub-N are more evolved YSOs than in hub-S, which still harbors deeply embedded objects in an earlier evolutionary stage. Therefore, the cluster of low-mass protostars detected in the infrared and associated with hub-N could effectively heat the gas and thus prevent fragmentation (Krumholz 2006;Krumholz et al. 2007).
Thus, the evolutionary stage of each hub could provide an explanation to the observed differences in the fragmentation level. Clearly, this scenario should be further investigated together with the information of the magnetic field at small scale.
6. CONCLUSIONS We present 1.3 mm continuum observations carried out with the SMA toward hub-N and hub-S of the IRDC G14.225-0.506 complemented with single-dish continuum observations at 870 and 350 µm conducted with the APEX and CSO telescopes. The main findings of this work are summarized as follows.
1. The two hubs presented in this work show different level of fragmentation in the millimeter range. Toward hub-N we find 4 millimeter condensations, while hub-S displays a higher degree of fragmentation, with 13 millimeter condensation detected within a field of view of 0.3 pc of diameter.
2. We fitted simultaneously the radial intensity profiles at 870 and 350 µm and the SED assuming spherical symmetry. In both hubs the radial density is well described by a Plummer-like function with a power-law index of 2.2, a central density of ∼ 1.3 × 10 −18 g cm −3 , and a radius of the flat inner region of 21000 AU, while temperature can be described as a power-law function with a power-law index of 0.34 and a dust emissivity index of 1.8. The differences in the fragmentation level cannot be explained by the density and temperature structure of these hubs.
3. The CFE is 3 % and 13 % in hub-N and hub-S, respectively. The derived masses and the observed number of fragments are consistent with pure thermal Jeans fragmentation, while adding non-thermal support implies that these hubs should not fragment. The level of fragmentation can thus be explained without invoking the presence of turbulent support.
4. All the derived properties such as the level of turbulence and magnetic field are remarkably similar in these hubs.
In particular, the non-thermal turbulent contribution to the velocity dispersion is ∼ 0.6 km s −1 in both hubs, and despite internal motions are slightly supersonic they do not seem to dominate at scales of 0.15 pc. On the other hand, both the value of the magnetic field as well as its orientation with respect to the main hub axis are remarkably similar in both hubs. Differences in the polarization properties at smaller scales, not traced by our H-band polarization data, could provide a reasonable explanation to the observed differences in the fragmentation level: a stronger magnetic field would suppress fragmentation, as it would be the case of hub-N. 5. Evolutionary effects and possibly feedback effects from the UV radiation of a nearby massive B1 star in hub-N could be an alternative explanation for the observed differences in the fragmentation level in these hubs.
In summary, our analysis of the fragmentation level in the twin hubs of the IRDC G14.2 reveals that the fragmentation process at scales 2000 AU is more consistent with pure thermal Jeans fragmentation than with fragmentation including turbulent support, and that agents such as the magnetic field, the evolutionary stage, and radiative feedback might be crucial and need to be considered in future work.