EDP Sciences
Free Access
Issue
A&A
Volume 545, September 2012
Article Number A98
Number of page(s) 12
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201118706
Published online 13 September 2012

© ESO, 2012

1. Introduction

It is well established that low-mass stars form from the collapse of dense prestellar cores. The star formation process can be divided both observationally and theoretically into an evolutionary sequence. From observations, protostars are classified according to the parametrized properties of the observed spectral energy distribution (SED) and appearance. For instance, a Class 0 source (André et al. 1993), which corresponds to the phase where the mass of the protostar is less than half the mass of the surrounding envelope, is defined according to its bolometric temperature Tbol (Chen et al. 1995) or to the ratio of bolometric to submillimeter luminosities Lbol/Lsmm (André et al. 1993). At later stages, the Class I, Class II, and Class III sources can be classified according to the slope αIR = dlog(λFλ)/dlog(λ) between 2.2 μm and 10–25 μm (Lada 1987). These classifications attempt to replicate the different evolutionary stages of a protostar. Observational ambiguities (e.g., projection effects, episodic accretion), however, imply that these classifications may not always reflect the true evolutionary stage of an object (e.g., Dunham et al. 2010).

The protostellar collapse phase is divided into two phases. First the prestellar dense core (gravitationally bound, no protostar, and no outflow) collapses isothermally and is able to freely radiate its gravitational energy into space, until it becomes dense (ρ ~ 10-13 g cm -3) and opaque enough for that radiation to be trapped. At this time, the gas begins to heat up and the first hydrostatic core (FHSC) is formed (Larson 1969). The FHSC is optically thick and its temperature increases adiabatically until it reaches ~2000 K, where H2 dissociation starts. This endothermic reaction allows the gas to undergo a second collapse phase until stellar densities are reached, i.e., the formation of the second hydrostatic core (SHSC), the protostar that corresponds to the protostellar Class 0 evolutionary stage. The FHSC is thus a transition phase from prestellar cores to Class 0 objects and has not been yet identified conclusively because its observational appearance is unclear. The predicted timescales of the FHSC phase range from a few hundred years in spherical models (Masunaga & Inutsuka 2000) to a few thousand years in very-low-mass non-magnetized rotating dense cores (Tomida et al. 2010a).

In recent years, significant progress has been achieved in both theory and observation thanks to increasing computing power and new observational capabilities provided by Spitzer and Herschel. For instance, the large Spitzer “Core to Disks” (c2d) Legacy Program survey (Evans et al. 2003) helped to increase the number of known protostars up to the point where first detections of FHSC can now be expected, given the short duration of the FHSC phase (Enoch et al. 2010; Dunham et al. 2011). Spitzer also helped to identify a new class of low-luminosity objects: the Very Low Luminosity Objects (VeLLOs, e.g., di Francesco et al. 2007). These are embedded in a dense envelope and have low intrinsic luminosity, L < 0.1 L. It is still a matter of debate whether VeLLOs are dense cores undergoing first collapse or low-mass protostars with a low-accretion rate. In a recent study, Kim et al. (2011) showed that at least in one case (the Bok globule CB130-1), the latter is the most likely explanation. Prestellar dense core populations have been identified over the past twenty years using space- and ground-based telescopes (e.g., Ward-Thompson et al. 1994; Motte et al. 1998; Johnstone et al. 2000; Nutter & Ward-Thompson 2007). More recently, Herschel observations have identified hundreds of prestellar cores and Class 0 sources (e.g., André et al. 2010; Bontemps et al. 2010; Henning et al. 2010). Combining different instrument data, Launhardt et al. (2010) presented the results of a comprehensive infrared, submillimeter, and millimeter continuum emission study of isolated low-mass star-forming cores, and identified nine Class 0 sources. Finally, a few FHSC candidate detections have been reported in recent years thanks to the increasing spectral and spatial scale coverage (Belloche et al. 2006; Chen et al. 2010; Enoch et al. 2010; Dunham et al. 2011; Pineda et al. 2011; Pezzuto et al. 2012). Their confirmation remains controversial, however, because the expected observable signatures of FHSC are still uncertain.

In parallel, enormous work has been put into developing numerical hydrodynamical models that include magnetic fields and radiative transfer (e.g., Price & Bate 2009; Commerçon et al. 2010; Tomida et al. 2010b) but relatively little has been done in producing synthetic observations that are directly comparable with actual observations. Most of the observables used to confront theory and observation remain statistical in nature, such as the initial mass function (IMF) or the binary distribution, whereas observations report mostly dust continuum and molecular line emission data. Therefore, it is vital that the interpretation of the data is sound, necessitating the use of non-symmetric models incorporating complex physics to derive predictions of dust continuum and molecular line emission observations.

For the early stages of star formation, Boss & Yorke (1995), Masunaga et al. (1998), Omukai (2007), Yamada et al. (2009), Tomida et al. (2010a), Tomisaka & Tomida (2011), and Saigo & Tomisaka (2011) made FHSC observational predictions of either dust emission or molecular line emission. Young & Evans (2005) and Dunham et al. (2010) also presented evolutionary signatures of the formation of low-mass stars, but used a crude approximation of the FHSC and the evolution of the envelope during collapse. To date, there has been no systematic study of observational predictions for FHSCs that explore physical conditions, and in particular on the magnetization level. Magnetic fields have indeed been found to control the first collapse and fragmentation phases (e.g., Hennebelle & Teyssier 2008; Commerçon et al. 2010) and strong magnetization seems to be a favored scenario for reproducing observational data of Class 0 multiplicity (Maury et al. 2010). The ability of a disk to fragment depends on the magnetization level, leading to potential changes in the SED.

In this paper, we derive SED evolutionary sequences and classical observational indicators such as Lbol and Tbol from state-of-the-art radiation-magneto-hydrodynamic (RMHD) 3D calculations of collapsing 1 M dense cores with different initial magnetic field strengths. We address the questions of whether or not SEDs can help in distinguishing physical conditions (e.g., magnetic field strengths) and how SEDs can help in targeting FHSC candidates. This paper is the first of a series and is associated with a companion paper in which we produce synthetic ALMA dust continuum emission maps (Commerçon et al. 2012, hereafter Paper II).

The paper is organized as follows: in Sect. 2, we present the numerical codes we used and our post-processing methodology. The RMHD calculations are qualitatively presented in Sect. 3. In Sect. 4, the RMHD calculations are post-processed to produce synthetic SEDs. In Sect. 5, we discuss our results and potential observational tests. Section 6 concludes our paper.

2. Method

2.1. The RAMSES code

For the dynamical RMHD calculations, we used the adaptive mesh refinement (AMR) code RAMSES (Teyssier 2002), which integrates the equations of ideal MHD (Fromang et al. 2006; Teyssier et al. 2006) and the equations of radiation hydrodynamics under the gray flux-limited diffusion approximation (Commerçon et al. 2011c). The RAMSES code is well designed for collapse calculations and has been extensively used in low-mass star formation studies (e.g., Hennebelle & Fromang 2008; Commerçon et al. 2008; Hennebelle & Ciardi 2009; Commerçon et al. 2010). It was also recently used in the high-mass star formation framework (Hennebelle et al. 2011; Commerçon et al. 2011b).

2.2. The RADMC-3D code

To post-process the RMHD calculations, we used the 3D radiative transfer code RADMC-3D1, developed by Dullemond. RADMC-3D can work with 3D AMR grids such as the RAMSES grid. RADMC-3D has the capability to perform dust continuum calculations with anisotropic scattering and gas/molecular line transfer calculations (LTE or non-LTE). The dust temperature, emission, and absorption are computed using a thermal Monte Carlo method. RADMC-3D is able to compute images and spectra from dust continuum radiative transfer. For non-LTE line radiative transfer calculations, RADMC-3D can estimate the level populations using approximate methods, such as the large velocity gradient (LVG) method (Shetty et al. 2011). In this paper, we will only use the dust continuum features of RADMC-3D and neglect scattering.

thumbnail Fig. 1

Top: scattering and absorption coefficient of Semenov et al. (2003) as a function of the wavelength for the homogeneous spheres model (with Fe/Fe + Mg = 0.3 “normal” silicate composition). Bottom: corresponding ratio between the scattering to absorption coefficients as a function of wavelength.

Open with DEXTER

2.3. Interface

First, we assumed that the dust and the gas are well coupled at high densities (ρ > 10-19 g cm-3 or n > 3 × 104 cm-3; Galli et al. 2002) to estimate the dust temperature. Second, we assumed that the correct temperature is the gas temperature computed in the dynamical RMHD calculations using the gray FLD approximation method. This second assumption is likely correct for FHSC formation and evolution, as shown recently by Vaytet et al. (2012) for 1D spherically symmetric gray and multi-frequency dense core collapse calculations. To obtain more accurate temperatures, we could have used the thermal Monte-Carlo capabilities of RADMC-3D. During the FHSC stage, however, the only sources of photons are from the compressional work and the accretion luminosity, which are not trivial to extract from the RMHD calculations.

Thanks to the AMR grid capabilities of RADMC-3D, we directly loaded the AMR structure of RAMSES within RADMC-3D. We thus did not lose any information between the RMHD and the post-processing radiative transfer calculations. We show in Appendix A the importance of using the AMR grid instead of remapping the data onto a uniform grid.

We used the low-temperature opacities of Semenov et al. (2003) for the homogeneous sphere model (with Fe/Fe + Mg = 0.3 “normal” silicate composition). For the RMHD calculations, we used a parametrization of the gray Rosseland and Planck means as functions of the gas temperature and density. For the post-processing radiative transfer calculations, we used the corresponding frequency-dependent opacities. Figure 1 shows the absorption and scattering opacities we used, as well as the ratio of scattering to absorption coefficients. The scattering coefficient is always less than the absorption coefficient, while their ratio is greater than 10% for visible to long near-infrared wavelengths (i.e., 0.2 μm < λ < 8 μm). Since we neglected scattering, we may have overestimated the amount of radiation that escapes at shorter wavelengths. Neglecting scattering remains a good approximation, however, since in this study we focus on the internal luminosity of collapsing dense cores at a stage earlier than the protostar+disk stage where SEDs peak at shorter wavelengths (e.g., Dunham et al. 2010).

thumbnail Fig. 2

Column density and temperature maps in the xy-plane (two upper rows) and xz-plane (two lower rows) for the three RMHD calculations: MU2 (left), MU10 (middle), and MU200 (right). Scales are logarithmic, and column density N is given in cm-2.

Open with DEXTER

3. Radiation-magneto-hydrodynamics models

3.1. Initial conditions

We used similar initial conditions as in Commerçon et al. (2010). We considered 1 M uniform density and temperature sphere in solid body rotation about the z-axis. The initial radius of the sphere is R0 = 4.93 × 1016 cm (~3300 AU). The initial angular velocity is set using a ratio of the rotational to gravitational energies β = 0.045. To favor fragmentation, we introduced an m = 2 azimuthal density perturbation with an amplitude of 10%. The corresponding free-fall time and orbital time are tff ~ 33 kyr and torb ~ 5.15 × 102 kyr, respectively. The initial density is ρ0 = 3.97 × 10-18 g cm-3 and the initial temperature is T0 = 11 K. The ratio of specific heats γ was set to 5/3. In this study, we did not include any additional microphysics such as H2 dissociation. Therefore, our models are only able to follow the dynamics up to a temperature of 2000 K, at which the second collapse starts.

The magnetic field is initially uniform and parallel to the rotation axis. The strength of the magnetic field is expressed in terms of the mass-to-flux to critical mass-to-flux ratio μ = (M0/Φ)/(M0/Φ)c. In this study, we explore three magnetization degrees: μ = 200, which corresponds to a quasi-hydrodynamical model (model MU200), μ = 10 for an intermediate magnetic field strength (model MU10), and μ = 2 for a strong magnetic field case (model MU2).

We used a coarse grid of 643 cells with ten levels of refinement, which corresponds to a maximum effective resolution of 0.2 AU. Calculations were performed using the HLLD Riemann solver (Miyoshi & Kusano 2005). Following up on former studies (Commerçon et al. 2008, 2010), we imposed at least 12 cells per Jeans length as a grid refinement criterion.

3.2. Qualitative description

In the following, we describe only qualitatively the results of the three RMHD calculations since the focus of this paper is the observational characteristics of the models and not their physical quantitative description. We synchronized the calculations at time t0 at which the maximum level of refinement is reached, which corresponds roughly to the formation of the FHSC. Column density and temperature maps for the three calculations are shown in Fig. 2. All fragments that are formed in the calculations are FHSCs.

3.2.1. The strongly magnetized case – μ = 2

In the strongly magnetized case, fragmentation is totally suppressed and no rotationally supported disk is formed because of the efficient magnetic braking (e.g., Joos et al. 2012). At time t0 + 0.78 kyr (with t0 = 50.4 kyr), there are 1.66 × 106 cells in the whole computational domain. An outflow is launched along the rotational axis with velocity ~3.7 km s-1 and it extends up to ~550 AU (not seen in the zoom in Fig. 2). Only one object is formed, surrounded by a thick pseudo-disk region. The mass of the FHSC is ~0.06 M and its mean temperature is ~375 K (Tmax ~ 630 K). The heated region around the central FHSC is almost spheroidal, which cannot be obtained in models using a barotropic approximation. In the latter case, the warm region would correspond to the dense region (the temperature is set by the density), i.e., the pseudo-disk (see Commerçon et al. 2010). These differences in the temperature distribution are fundamental to the SED estimates. We ran the MU2 calculations until time t0 + 1.2 kyr, where Tmax ~ 1160 K and the outflow extends to ~850 AU, i.e., just before the start of the second collapse. At temperature T > 1200 K, dust begins to evaporate, which considerably lowers the opacity by more than two orders of magnitude (e.g., Lenzuni et al. 1995). The high temperature region is then able to radiate efficiently and the temperature increases fast towards 2000 K, at which time H2 dissociation starts.

3.2.2. The intermediately magnetized case – μ = 10

In the intermediately magnetic field strength model, there is no fragmentation, but an outflow of extent ~300 AU and velocity ~1 km s-1 is launched and a rotationally supported disk of radius ~100 AU is formed. This result agrees well with previous studies of Hennebelle & Ciardi (2009) and Commerçon et al. (2010). At time t0 + 1.81 kyr (with t0 = 35.7 kyr), there are 1.3 × 106 cells in the whole computational domain. The mass of the FHSC is ~0.06 M and its mean temperature is ~490 K (Tmax ~ 810 K). From the temperature maps, we distinguish two shocked regions (i.e., where there is a weak spike in gas temperature) in addition to the accretion shock on the central FHSC itself. In the xy-plane, the shocked region corresponds to the border of the disk. In this optically thin region, it is then easy to see the spike in gas temperature that is characteristic of radiative shocks and has a length of a few photon mean-free paths (e.g., Mihalas & Mihalas 1984). In the xz-plane, we identify a shock at the outflow border, with a spike in the gas temperature of the same amplitude as at the disk border. In comparison with the MU2 model, the outflow is broader and denser, and the vertical extent of the pseudo-disk is strongly diminished. The calculations of the MU10 model were run until Tmax ~ 1200 K at time t0 + 2.94 kyr, just before the second collapse is expected to begin.

thumbnail Fig. 3

Spectral energy distributions for the three models at different times: MU2 (top row), MU10 (middle row), and MU200 (bottom row). The first panels in the left column show the SEDs obtained from the initial conditions. SEDs are computed at each epoch for three different inclination angles: θ = 0° (pole-on view, solid line), θ = 60° (dotted line), and θ = 90° (edge-on view, dashed line).

Open with DEXTER

3.2.3. The very weakly magnetized case – μ = 200

In the very weakly magnetized model, the collapsing core forms a prestellar disk that fragments into several objects (five at time t0 + 3.26 kyr), but no outflow is launched because of the weak initial magnetic field. The disk radius extends up to ~150 AU. At time t0 + 3.26 kyr (with t0 = 35.4 kyr), there are 1.95 × 106 cells in the whole computational domain. The mass of the fragments is about the same: M ~ 0.043 M for the central fragment and M ~ 0.036−0.039 M for the orbital fragments. The mean temperature is T ~ 270 K for the central fragment and T ~ 245−263 K for the orbital fragments. Heating of the gas around the central fragment is much less efficient than in the two previous simulations. In addition, we see from the temperature maps in the xz-plane that the radiative shock at the vertical borders of the disk is sharper and is located closer to the center compared to the MU10 run. The infall velocity and thus the shock strength are also higher. Before being accreted by a fragment, gas has already lost a part of its gravitational energy through that shock, and the accretion shock on the fragment is then weaker. The MU200 calculations were halted before the onset of the second collapse, the first core lifetime being relatively long for the hydrodynamical case (Tomida et al. 2010a).

From this qualitative description, we clearly see that the three models are different and relevant to the various conditions that may be found at the FHSC stage. In the following, we will see that these differences are also expressed in synthetic SEDs of these collapsing objects.

4. Synthetic observations

4.1. Definitions of classical physical stage indicators

The SEDs of the models were used to calculate classical observational signatures. We monitored for each model the temporal evolution of the bolometric luminosity Lbol, the bolometric temperature Tbol, and the ratio of bolometric to submillimeter luminosity Lbol/Lsmm. The bolometric luminosity Lbol is defined as (1)where d is the distance to the object, Sν is the flux density, and ν is the frequency. The bolometric temperature Tbol is defined following Myers & Ladd (1993), i.e., as the temperature of a blackbody with the same mean frequency as the observed SED, which can be approximated by (2)Following André et al. (1993), we also calculated the ratio Lbol/Lsmm, with the submillimeter luminosity estimated as (3)These three quantities are commonly used to indicate the physical stage of observed sources. Tbol and Lbol/Lsmm are both used to define class boundaries. For our purpose, only the boundary between Class 0 and Class I matters. André et al. (1993) defined Class 0 sources as objects with Lbol/Lsmm ≤ 200, while Chen et al. (1995) defined Class 0 objects with Tbol < 70 K. Both classifications are problematic for models because one has to define the region over which to integrate the luminosity. Finally, VeLLOs are defined as objects with an internal luminosity Lint < 0.1 L, where Lint is the luminosity of the object in excess of that supplied by the interstellar radiation field (di Francesco et al. 2007).

4.2. Spectral energy distribution time evolutionary sequence

In this subsection, we qualitatively describe the SEDs obtained for the three models. The fluxes of each SED shown below were computed and integrated over a lobe size of 3000 AU  ×  3000 AU (about half the initial size of the core) and observed at a distance of 150 pc (a typical distance to star-forming regions like Taurus). These choices give an aperture of ~20′′ for our synthetic observations, which compares well with current instrument resolutions. Our choice of fixing the same object size for all wavelengths is justified because observations often report flux over a given area whatever the beam size. In addition, there are many different instruments that report star-forming region observations and it is impossible to match all beam sizes in a single SED.

Figure 3 shows the SEDs at different epochs for the three models. For each epoch, SEDs are presented for three different inclination angles: θ = 0° (pole-on view, solid line), θ = 60° (dotted line), and θ = 90° (edge-on view, dashed line). We also computed SEDs for models at an inclination angle of θ = 45°, but these were indistinguishable from those obtained with an inclination of θ = 0°. As long as t < t0 in all models (see Sect. 3.2 for the definition of t0), SEDs are indistinguishable from the initial one and no emission is found at wavelengths shorter than 30 μm. This behavior is not surprising because the FHSC has not yet formed and the gas is still isothermal at T = 11 K.

Surprisingly, the SEDs of the MU2 and MU200 models have similar time evolutionary sequences compared to the MU10 model, despite their totally different physical conditions. In the MU2 and MU200 models, the departure from the initial conditions starts when t > t0 and the flux increases at all wavelengths as temperature increases within the FHSC. The peak of the SED, however, does not change after the FHSC formation. We also find that a significant amount of radiation is emitted between 20 μm and 100 μm in the late evolution of the MU2 and MU200 models, in contradiction to previous studies (e.g., Masunaga et al. 1998; Omukai 2007), which found no observable emission below 30−50 μm. This difference can be explained by geometrical effects. Masunaga et al. (1998) and Omukai (2007) assumed spherical symmetry, whereas we considered a full 3D structure. In our models, the envelope seen with inclination angle θ < 60° is much less dense than that of spherical models. In the same manner, the flux we obtain with θ = 90° is strongly reduced in this wavelength range because the envelope has reprocessed the radiation emitted by the FHSC. Interestingly, the MU10 model does not show any evolutionary sequence even at the latest stage before the second collapse, independent of the inclination angle. In the edge-on view (θ = 90°), this lack of variation is not surprising, since the MU10 model combines a disk and a pseudo-disk structure. The radiation emanating from the FHSC is thus fully reprocessed. But even looking pole-on through the outflow, no flux below 30 μm is seen. From Fig. 2, we can see that the outflow in the MU10 model is much denser than that in the MU2 model; it acts as an envelope in which radiation from the FHSC is reprocessed.

thumbnail Fig. 4

Radius at which the optical depth τλ = 1 as a function of the wavelength for the three models at the same time as in Fig. 2: MU2 (black), MU10 (red), and MU200 (blue). The solid line shows optical depth with an inclination of θ = 90°, and the dashed line that of θ = 0°. The dotted line at R = 10 AU roughly represents the FHSC’s outer boundary.

Open with DEXTER

Figure 4 shows the radius at which the optical depth τλ = 1 (the transition from the optically thin to thick regime) as a function of the wavelength for the three models at the same time as in Fig. 2 and for two inclinations: θ = 90° (solid line, edge-on) and θ = 0° (dashed line, pole-on). The dotted line at R = 10 AU roughly represents the FHSC outer boundary. The sharp jumps observed in the wavelength range between ~2 μm and ~50 μm are due to jumps in the opacity itself owing to ice or silicate for instance (see Fig. 1). Because the opacity increases when the wavelength decreases, the optically thick radius is larger at short wavelengths. For θ = 90°, the optical depth is about the same in the three models, since the pseudo-disk or the disk play the same role in obscuring radiation from the FHSCs. The central FHSC boundary can thus only be seen in an edge-on view at wavelengths >5 mm. The disk or pseudo-disk features alone (typical size of about 100 AU) are nevertheless observable at shorter millimeter wavelengths, although these features are not reliable diagnostics of the presence of FHSCs, since they will not disappear as the second collapse proceeds. For θ = 0°, the radius at which the optical depth reaches unity corresponds to the FHSC radius for wavelengths >100 μm, and even marginally for shorter wavelengths 3 μm < λ < 100 μm for the MU2 model. At inclination angle θ < 60° (since the θ = 0° SED closely resembles the θ = 60° one), the flux received at wavelengths >100 μm traces the FHSC outer boundary, and the accretion luminosity can escape. This result clearly shows the importance of properly resolving the accretion shock on the FHSC with a full RHD solver and can explain the differences found in the past literature. It also has some important implications for observability since an FHSC should then be observable as a point source at mid- to far-infrared wavelengths if the inclination is favorable.

thumbnail Fig. 5

Bolometric luminosity (top) and bolometric temperature (bottom) as a function of time for the three models: MU2 (left), MU10 (middle), and MU200 (right). Bolometric quantities are computed from models with three different inclination angles: θ = 0° (black), θ = 60° (red), and θ = 90° (blue).

Open with DEXTER

Last but not least, there is a clear FHSC lifetime dependance on the initial magnetization of the cloud. Whereas the FHSC has a lifetime >4 kyr in the MU200 model, it is shorter in the other models as magnetization increases, i.e., ~3 kyr for MU10 and ~1.2 kyr for MU2. The evolution is faster as magnetization increases, since the accretion rate increases with stronger magnetic braking that transports angular momentum outward (Commerçon et al. 2010; Joos et al. 2012). As already mentioned, Tomida et al. (2010a) reached a similar conclusion about the FHSC lifetime using RHD models of different initial mass rotating dense cores.

4.3. Lbol, Tbol, and Lbol/Lsmm

Figure 5 shows the evolution of the internal luminosity, i.e., the difference between the bolometric luminosity at time t and the initial bolometric luminosity (top panels) and the bolometric temperature (bottom panels) as a function of time for the three models and three inclination angles: θ = 0°, 60°, and 90°. Time t0 is indicated by the vertical dashed line. The dotted line in the bolometric luminosity plots indicates the VeLLO boundary (di Francesco et al. 2007), Lint < 0.1 L. The bolometric luminosity and temperature globally increase after the formation of the FHSC for inclination angles θ < 60°. Only the edge-on view gives rise to decreasing bolometric luminosity because of the disk or pseudo-disk structure. More interestingly, the luminosity in the θ < 60° cases increases relatively slowly and becomes greater than 0.1 L only at late times, long after the FHSC formation. FHSCs will thus appear as VeLLOs for most or all of their lifetime.

In the bolometric temperature plots, the horizontal dotted line shows the Class 0/Class I boundary of Chen et al. (1995), i.e., Tbol = 70 K. In terms of Tbol, FHSCs are indistinguishable from starless cores and remain at the lower end of the Class 0 Tbol range. The same behavior is obtained for the ratio Lbol/Lsmm. In summary, these classical definitions cannot help in distinguishing between an FHSC and a more evolved Class 0 source, or even between an FHSC and a starless core.

thumbnail Fig. 6

Spectral energy distribution evolution as a function of time and inclination for the MU2 model. The black line represents the SED when the FHSC has not yet been formed, whereas the red lines indicate SEDs after FHSC formation. The gray area around the initial SED reflects the uncertainties in the initial temperature of the prestellar core (9 K < T0 < 13 K). The blue cross represents the sensitivity at 24 μm (0.83 mJy at 3σ/24 s) for the Spitzer c2d Legacy program (Evans et al. 2003) and the horizontal blue line the sensitivity of the Herschel PACS instrument (~4.4 mJy in point source mode, 5σ/1 h; Poglitsch et al. 2010).

Open with DEXTER

4.4. FHSC signature

We have shown in the previous section that classical observational signatures are not reliable evolutionary indicators during the early collapse stage. Thanks to the increasing spectral coverage and sensitivity of new instruments, however, SEDs are increasingly more accurate, and new indicators can be defined to distinguish between a collapsing core with an FHSC and a prestellar (starless) core.

Figure 6 shows the SED evolution as a function of time and inclination for the MU2 model at stages before (black lines) and after (red lines) FHSC formation. The blue cross represents the sensitivity at 24 μm (0.83 mJy at 3σ/24 s) for the c2d Legacy program (Evans et al. 2003) and the horizontal blue line represents the inflight sensitivity of the HerschelPACS instrument at 70 μm, derived from science observations (~4.4 mJy in point source mode, 5σ/1 h, Poglitsch et al. 2010). Shortly after FHSC formation, the SED does not significantly depart from the initial one. After the FHSC begins to build up its mass, however, the accretion shock becomes stronger and departures from the initial SED are much clearer. To take into account the uncertainties on the initial temperature of the prestellar cores, we show SEDs of prestellar cores with initial temperatures ranging from 9 K to 13 K (about 20% upper limit uncertainty; Launhardt et al., in prep.). The SEDs with accreting FHSC are clearly distinct from the initial uncertainties on the temperature for inclination θ < 60° although the SED peak itself does not show evidence of an FHSC. According to the sensitivity of the Herschel PACS instrument, FHSCs can be identified with sufficiently deep observations. More interestingly, the sensitivity of Spitzer at 24 μm could also be sufficient for FHSC detection at late stages of their lifetime. The wavelength range over which the flux has increased the most after the FHSC formation is 20 μm < λ < 100 μm (at least a factor 10 at 100 μm for θ < 60°), whereas Lsmm increases only by a factor less than 2. Similarly, the ratio Lbol/Lsmm increases at maximum by only a factor ~2 after FHSC formation.

The estimate of Lbol from observations strongly depends on the number of points in the SED, i.e., the spectral coverage. Instead of the physical stage indicator Lbol/Lsmm, we chose here to express a related luminosity ratio using a particular wavelength. In this example, we chose the luminosity at 70 μm to show how it is possible with the Spitzer/MIPS or Herschel/PACS instruments to identify FHSC candidates in comparison to starless cores. We note that a similar analysis can be performed using the 100 μm, depending on the available observations. The flux at this wavelength is well constrained, and we just showed that Lsmm increases only by a factor 2 after the formation of the FHSC. Finally, we note that similar flux increases in the wavelength range 20 μm < λ < 100 μm are also observed in 10 M dense core collapse calculations since, even if the envelope is thicker, the accretion shock on the FHSC is also stronger (e.g., for massive dense core models, see Commerçon et al. 2011b).

thumbnail Fig. 7

Top: evolutionary diagram of the flux at at 24 μm (solid line), 160 μm (dashed line), and 250 μm (dotted line) as a function of the flux at 70 μm for the three models and for the the same three inclination angles as previously. The horizontal gray line indicates the aforementioned 24 μm sensitivity limit of Spitzer and the vertical dotted line the 70 μm sensitivity limit of the Herschel PACS instrument.Bottom: color–color plot representing log(λFλ160/λFλ70) as a function of log(λFλ70/λFλ24) for the three models with the same three inclination angles.

Open with DEXTER

thumbnail Fig. 8

Dust continuum emission maps for the MU2 models at four wavelengths, from left to right: 24 μm, 70 μm, 160 μm, and 1200 μm. Three different times are shown that trace the evolutionary stage of the source: t0 − 26.43 kyr (top row, i.e., starless core), t0 + 0.7 kyr (middle row, i.e., faint FHSC) and t0 + 1.13 kyr (bottom row, i.e., end of the FHSC lifetime). The fluxes are given in Jy/px (1 px = 0.08′′ × 0.08′′) for objects at a distance of 150 pc and observed pole-on.

Open with DEXTER

Figure 7 (top) shows the time evolution of the flux at 24 μm (solid line), 160 μm (dashed line), and 250 μm (dotted line) as a function of the flux at 70 μm and of the inclination angle for the three models. The horizontal gray line indicates the aforementioned 24 μm sensitivity limit of Spitzer and the vertical dotted line the 70 μm sensitivity limit of the Herschel PACS instrument. For θ < 60°, the flux F70 increases by about two orders of magnitude after FHSC formation for the MU2 and MU200 models, whereas the MU10 model does not show any increase, as expected, since its SED does not show any evolution. Note that the flux received at 70 μm is always above the sensitivity limits of the current instruments. The most dramatic increase is found in the flux at 24 μm, which is almost zero at FHSC formation and then increases to detectable values at late times. Similar to Lsmm, the fluxes at 160 μm and 250 μm remain almost constant throughout the FHSC lifetime. The fluxes at 24 μm and 70 μm can thus help in identifying FHSC candidates with respect to starless cores if the physical conditions are appropriate (either a strong or a very weak magnetic field), and if the objects are observed with a favorable inclination <60°. The evolutionary tracks of the FHSCs in the color–color plane of log(λFλ70/λFλ24) − log (λFλ160/λFλ70) are shown in Fig. 7 (bottom) for the three models at the three same inclination angles as discussed previously. There is a clear evolutionary track followed by the accreting FHSC that can be directly compared to observations.

5. Discussion

5.1. Comparison to previous work

Spectral energy distributions of FHSCs embedded in their parent cores have already been studied for instance by Boss & Yorke (1995), Masunaga et al. (1998), Omukai (2007), Tomida et al. (2010a), and Saigo & Tomisaka (2011). While they all found that the radiation emitted by the FHSC is essentially completely reprocessed by the dust in the surrounding core to 50–200 μm with no observable emission below 30 μm, we found that there is indeed observable emission at shorter far-infrared wavelengths down to 20 μm. These differences may be explained by the different physical models and post-processing of each study. While we directly post-processed 3D AMR RMHD calculations, Boss & Yorke (1995) used a spherical grid and integrated the RHD equations with a low numerical resolution. Their FHSC temperature is also relatively low (<200 K). Their results are therefore similar to ours shortly after FHSC formation. In addition, as mentioned previously, Masunaga et al. (1998) and Omukai (2007) used spherical models and their SEDs are then similar to ours with θ = 90°. Next, Saigo & Tomisaka (2011) post-processed 3D high-resolution collapse calculations with plane-parallel approximation but neglected magnetic fields and used a barotropic EOS. They observed an increase of the total energy in the SED but their SEDs were calculated over a very small area of 40 AU × 40 AU and therefore excluded the large-scale emission dominant at wavelengths >100 μm. Last, Tomida et al. (2010a) computed SEDs from 3D RHD calculations but their cores do not fragment (slow rotation). In addition, they used Bonnor-Ebert initial density profiles in which accretion rates are lower and therefore FHSCs are fainter.

5.2. What can SEDs tell about physical processes and evolutionary stages?

Figures 3 and 6 show that the evolution of an FHSC throughout its lifetime can be reflected in observable changes of the SED, but only under certain conditions. Only the intermediate MU10 model does not show any evolutionary sequence because of the thick envelope structure consisting of a system disk + pseudo-disk + outflow. The MU2 and MU200 models show an evolutionary sequence very different than that of the MU10 model, but very similar to each other. The MU2 model is strongly magnetized and exhibits a pseudo-disk + outflow system, whereas the MU200 model fragments and exhibits a rotationally supported disk. An FHSC with a clear evolutionary sequence can then be observed in simple geometrical configurations and with inclination angle θ < 60°. Nevertheless, it is impossible to distinguish between a magnetized and a non-magnetized scenario from a strict SED point of view. High angular resolution images of, at least, dust emission are required to distinguish between these scenarios. This kind of diagnostic can only be provided by interferometry at longer wavelengths (see Paper II).

We have shown than it can be relatively straightforward to distinguish a starless core from an FHSC thanks to flux received at mid- to far-infrared wavelengths, i.e., at 24 μm and 70 μm, but it remains hard to distinguish between an FHSC and a young stellar object (YSO) with a disk (Pineda et al. 2011), i.e., an accreting second hydrostatic core (SHSC). This difficulty is due to the thick envelope in which the FHSC and SHSC are each embedded throughout the Class 0 phase that reprocesses short-wavelength radiation to longer wavelengths. Here, models clearly fail to predict any detectable differences. Using only dust continuum emission is not the appropriate way to differentiate the two core stages, since almost all radiation emanating from both the FHSC and the SHSC will be reprocessed by the thick envelope. Another approach could consist of focusing on the outflow properties (collimation, velocity, momentum, etc.) and on the chemical properties of the sources, such as gas phase and ice molecules abundances, which reflect the local thermal history (e.g., Kim et al. 2011). We will study these properties in upcoming papers.

5.3. Example of methodology to target first core

We showed in Sect. 4.4 that the fluxes at 24 μm and 70 μm can be good indicators of FHSCs in comparison with starless cores. We also did not find any significant emission below 10 μm. We infer that the combination of a point source detection at 70 μm or even at 24 μm with no detection at wavelength <10 μm can be a good indicator of FHSCs from a strict SED point of view. Under the current telescope sensitivity limits, a detection at 70 μm alone is most likely, since the emission at 24 μm is faint, and reaches the Spitzersensitivity limit only in the MU2 model in the second half of the FHSC lifetime. Figure 8 shows dust continuum emission maps for the MU2 model at three different evolutionary stages and for four wavelengths (24 μm, 70 μm, 160 μm, and 1200 μm). It clearly shows that the emission at 70 μm and 24 μm should be point-like after the FHSC formation. Note that at all evolutionary stages, a flux at 70 μm can be detected but its emission would not look like a point source for a starless core. It will instead be much more extended (i.e., the size of the initial dense core) because of the ISRF reprocessing. Thanks to the recent Herschel and Spitzer observations, observers should be able to combine data from both telescopes and make comparisons with our predicted color–color plot in Fig. 7. For instance, di Francesco et al. (2010) already used this procedure and combined maps of Herschel data with Spitzer ones to distinguish between prestellar and protostellar cores. Once FHSC candidates are identified from these flux examinations, the next step will be to target them with high-resolution interferometry such as ALMA (see Paper II) to constrain the physical properties (e.g., strength of the magnetic field).

As mentioned earlier, a complementary approach would be to use chemical analysis to separate the FHSC and the SHSC stages. At temperatures above ~1000 K, grain evaporation starts and refractory chemical species such as C and Si are released into the gas phase, considerably changing the chemical network and the opacities (e.g., Lenzuni et al. 1995). Deuterium-related chemistry could also trace an SHSC because deuterium-bearing species should be under-abundant in high-temperature regions (Albertsson et al. 2011).

Another more direct and more promising way to stage differentiation would be to find evidence of a low-velocity outflow without its high-velocity counterparts. Recent theoretical studies have shown that a low-velocity outflow can be launched during the FHSC stage (e.g., Hennebelle & Fromang 2008; Commerçon et al. 2010; Tomida et al. 2010b) without the high-velocity outflow that is always present at the SHSC stage (Machida et al. 2008). Not only may the morphology and energetics of the outflow help in identifying the cores, but also the chemical composition of the gas in the outflow. Since temperatures are different during the FHSC and SHSC stages, we may expect that the low-velocity outflow would show evidence of low-temperature gas-grain chemistry and that the high-velocity outflow would show evidence of high-temperature chemistry. In addition, the low-velocity outflow is launched from regions with temperatures <1000 K where grains are still not fully evaporated. Grains should then be present in the FHSC outflow and since the velocity is relatively low (a few km s-1), they are not expected to be destroyed in the shock at the envelope interface. Therefore, we may expect FHSC outflows to be seen in thermal dust emission if they are dense enough. It is worth mentioning that some dust continuum observations may be also contaminated by CO emission from outflow gas, making them easier to see (e.g., Drabek et al. 2012).

Recent work by Price et al. (2012) reported the formation of a collimated jet driven from the FHSC, supporting recent observations of an FHSC candidate in Per-Bolo 58 (Dunham et al. 2011). The timescale needed to produce such a collimated jet (~2 kyr), however, is longer than the expected lifetime of FHSCs in a magnetized collapsing dense core. In addition, Price et al. (2012) used sink particles (scales smaller than 5 AU are not described) and a barotropic approximation with an adiabatic index of γ = 7/5, for which FHSC lifetimes are even shorter than those we predict in this study with γ = 5/3 (the FHSC contracts adiabatically faster with γ = 7/5). It is therefore very likely that FHSCs reported in Price et al. (2012) should have already undergone second collapse and formed an SHSC within ~2 kyr after their formation.

5.4. Limits of the model

In this study, we used the gray radiative transfer approximation for the RMHD calculations which may not be the most appropriate method at later stages of the evolutions. As mentioned earlier, however, Vaytet et al. (2012) showed that multi-frequency collapse calculations were well-reproduced by the gray FLD calculations of Commerçon et al. (2011a), at least for the first collapse and FHSC stages. In addition, we used a constant ratio of specific heats γ, set to 5/3, although it is known that γ should change to a value of 7/5 as temperature increases, and H2 rotational degrees of freedom become available (e.g., Machida et al. 2008). With γ = 7/5, FHSCs compress and cool more efficiently and thus reach the second collapse more quickly (e.g., Commerçon et al. 2010), which reduces the FHSC lifetime.

For the RADMC-3D post-processing, we used one specific opacity model of Semenov et al. (2003) for all temperatures and densities. This choice may not be appropriate since grain mantle compositions change as temperature increases. Nevertheless, the changes in opacity remain relatively minor at temperatures <1000 K and we do not expect them to have a significant impact on our results. Finally, we have neglected scattering although it has been found to be important for wavelengths <10 μm (Young & Evans 2005; Dunham et al. 2010). Since we did not find any emission from FHSCs at these wavelengths and we do not consider this wavelength range significant for our conclusions, scattering effects are not relevant for our purposes.

6. Summary and perspectives

We used state-of-the-art radiation-magneto-hydrodynamic models of the first collapse and first hydrostatic core formation in the context of low-mass star formation to predict spectral energy distributions and related observational indicators. We investigated the collapse and the fragmentation of a 1 M dense core with three different magnetization levels up to the start of the second collapse. We showed that the FHSC lifetime strongly depends on the magnetization of the collapsing core and that a strong magnetic field leads to a shorter lifetime because magnetic braking raises the accretion rate. We post-processed the RMHD calculations with a 3D radiative transfer code and produced synthetic spectral energy distributions. We showed that FHSCs can be identified thanks to the rapid increase in the flux received between 20 μm and 100 μm for inclination angle θ < 60° under different physical conditions (e.g., high or weak magnetization). We also found that FHSCs are globally consistent with VeLLOs. The classical stage evolutionary indicators Tbol, Lbol, and Lbol/Lsmm are not reliable for identifying FHSCs since they do not show any significant evolution through the FHSC lifetime. We found that a point source detection at wavelengths shorter than 200 μm is a good indicator of FHSC in combination with the absence of detection at wavelengths <10 μm. We currently cannot, however, provide any diagnostics to distinguish an FHSC from a low-luminosity SHSC.

Recently, Herschel has mapped hundreds of dense cores through, e.g., the EPOS and Gould Belt Survey key observation programs (e.g., Henning et al. 2010; André et al. 2010; Bontemps et al. 2010) with PACS data at 70 μm. The community has also gathered many Spitzer observations at 24 μm and 70 μm. Combining the observations for common sources of both datasets could identify FHSC candidates. Cataloging these candidates will provide the basis required for observation proposals aimed at high-resolution interferometers such as the Submillimeter Array (SMA) and ALMA, which will be indispensable to accurately measure the physical properties of the collapsing cores (see below). Future progress is also expected from the next generation of space telescopes such as SPICA with the SAFARI instrument (with 1-sigma rms sensitivities of ~2 μJy at 24 μm and ~50 μJy at 70 μm; Swinyard et al. 2009), which could then test the evolutionary tracks we predicted.

The SED information alone cannot help in distinguishing the important physical mechanisms at work at early Class 0 stages (large disk formation and fragmentation, outflow launching, etc.). Only interferometric observations would be capable of such a feat, and forthcoming ALMA observations will give a clearer, if not definitive, answer to the fragmentation properties of Class 0 sources. ALMA dust continuum emission synthetic maps will be presented in Paper II. Predictions other than dust continuum emission only are also needed to target FHSC candidates more efficiently. Chemistry post-processing and molecular line emission calculations are therefore a natural next step for FHSC observation predictions. Future investigations will also aim to follow the second collapse with the appropriate physics (mainly H2 dissociation).


Acknowledgments

We thank the anonymous referee for her/his comments and corrections that improved the quality of the manuscript. The calculations were performed at CEA on the DAPHPC cluster. The research of B.C. is supported by the postdoctoral fellowships from Max-Planck-Institut für Astronomie and from the CNES. BC also acknowledges support from the French ANR Retour Postdoc program. B.C. thanks N. Vaytet, F. Levrier and A. Stutz for useful comments and discussions.

References

Appendix A: AMR versus cartesian grid for postprocessing

In this appendix, we present SEDs computed with two different post-processing methods. The first is the one we used throughout this paper, i.e., the AMR grid of RAMSES is directly loaded in RADMC-3D. The second projects the AMR data of RAMSES into a uniform cube of extent ~3000 AU × 3000 AU. The second method is much simpler for interfacing, but a lot of information is lost when smoothing the finer levels onto the uniform grid. In our calculations, since the finer levels correspond to the high-temperature and high-density regions, smoothing to a lower resolution may affect the SEDs.

Figure A.1 shows SEDs for the MU10 models computed on three different grids: the AMR grid for the black lines, a uniform 2563 cube for the red lines, and a uniform 5123 cube for the blue lines. Two different inclination angles are shown: θ = 0° (solid line) and θ = 90° (dashed line). The three SEDs obtained with θ = 90° are similar with no emission below 30 μm independent of the grid used for post-processing, but they are different for θ = 0°. For the latter case, there is indeed emission between 3 μm and 30 μm and an increase of the flux around 100 μm. We clearly see that poor resolution may strongly affect the SEDs, which solely comes from the fact that high-density regions have been smoothed out by reprojection of the AMR grid onto a uniform one.

thumbnail Fig. A.1

Spectral energy distribution for the MU10 models computed on three different grids: the AMR grid for the black lines, a uniform 2563 cube for the red lines, and a uniform 5123 cube for the blue lines. Two different inclination angles are also shown: θ = 0° (solid line) and θ = 90° (dashed line, almost indistinguishable from one to another).

Open with DEXTER

All Figures

thumbnail Fig. 1

Top: scattering and absorption coefficient of Semenov et al. (2003) as a function of the wavelength for the homogeneous spheres model (with Fe/Fe + Mg = 0.3 “normal” silicate composition). Bottom: corresponding ratio between the scattering to absorption coefficients as a function of wavelength.

Open with DEXTER
In the text
thumbnail Fig. 2

Column density and temperature maps in the xy-plane (two upper rows) and xz-plane (two lower rows) for the three RMHD calculations: MU2 (left), MU10 (middle), and MU200 (right). Scales are logarithmic, and column density N is given in cm-2.

Open with DEXTER
In the text
thumbnail Fig. 3

Spectral energy distributions for the three models at different times: MU2 (top row), MU10 (middle row), and MU200 (bottom row). The first panels in the left column show the SEDs obtained from the initial conditions. SEDs are computed at each epoch for three different inclination angles: θ = 0° (pole-on view, solid line), θ = 60° (dotted line), and θ = 90° (edge-on view, dashed line).

Open with DEXTER
In the text
thumbnail Fig. 4

Radius at which the optical depth τλ = 1 as a function of the wavelength for the three models at the same time as in Fig. 2: MU2 (black), MU10 (red), and MU200 (blue). The solid line shows optical depth with an inclination of θ = 90°, and the dashed line that of θ = 0°. The dotted line at R = 10 AU roughly represents the FHSC’s outer boundary.

Open with DEXTER
In the text
thumbnail Fig. 5

Bolometric luminosity (top) and bolometric temperature (bottom) as a function of time for the three models: MU2 (left), MU10 (middle), and MU200 (right). Bolometric quantities are computed from models with three different inclination angles: θ = 0° (black), θ = 60° (red), and θ = 90° (blue).

Open with DEXTER
In the text
thumbnail Fig. 6

Spectral energy distribution evolution as a function of time and inclination for the MU2 model. The black line represents the SED when the FHSC has not yet been formed, whereas the red lines indicate SEDs after FHSC formation. The gray area around the initial SED reflects the uncertainties in the initial temperature of the prestellar core (9 K < T0 < 13 K). The blue cross represents the sensitivity at 24 μm (0.83 mJy at 3σ/24 s) for the Spitzer c2d Legacy program (Evans et al. 2003) and the horizontal blue line the sensitivity of the Herschel PACS instrument (~4.4 mJy in point source mode, 5σ/1 h; Poglitsch et al. 2010).

Open with DEXTER
In the text
thumbnail Fig. 7

Top: evolutionary diagram of the flux at at 24 μm (solid line), 160 μm (dashed line), and 250 μm (dotted line) as a function of the flux at 70 μm for the three models and for the the same three inclination angles as previously. The horizontal gray line indicates the aforementioned 24 μm sensitivity limit of Spitzer and the vertical dotted line the 70 μm sensitivity limit of the Herschel PACS instrument.Bottom: color–color plot representing log(λFλ160/λFλ70) as a function of log(λFλ70/λFλ24) for the three models with the same three inclination angles.

Open with DEXTER
In the text
thumbnail Fig. 8

Dust continuum emission maps for the MU2 models at four wavelengths, from left to right: 24 μm, 70 μm, 160 μm, and 1200 μm. Three different times are shown that trace the evolutionary stage of the source: t0 − 26.43 kyr (top row, i.e., starless core), t0 + 0.7 kyr (middle row, i.e., faint FHSC) and t0 + 1.13 kyr (bottom row, i.e., end of the FHSC lifetime). The fluxes are given in Jy/px (1 px = 0.08′′ × 0.08′′) for objects at a distance of 150 pc and observed pole-on.

Open with DEXTER
In the text
thumbnail Fig. A.1

Spectral energy distribution for the MU10 models computed on three different grids: the AMR grid for the black lines, a uniform 2563 cube for the red lines, and a uniform 5123 cube for the blue lines. Two different inclination angles are also shown: θ = 0° (solid line) and θ = 90° (dashed line, almost indistinguishable from one to another).

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.