Strong stellar-driven outflows shape the evolution of galaxies at cosmic dawn

TL;DR

This study uses the GAEA galaxy evolution model to demonstrate that strong stellar-driven outflows significantly influence galaxy formation and evolution during cosmic dawn, successfully matching observed luminosity functions and star formation rates up to high redshifts.

Contribution

The paper introduces an improved galaxy evolution model with strong stellar feedback that accurately reproduces galaxy properties from redshift 4 to 10, extending previous work to earlier cosmic times.

Findings

GAEA reproduces UV and optical luminosity functions up to z~10.

Model matches observed galaxy stellar mass function up to z~7.

Dust attenuation levels are consistent with recent estimates.

Abstract

We study galaxy mass assembly and cosmic star formation rate (SFR) at high-redshift (z$\gt$4), by comparing data from multiwavelength surveys with predictions from the GAlaxy Evolution and Assembly (GAEA) model. GAEA implements a stellar feedback scheme partially based on cosmological hydrodynamical simulations, that features strong stellar driven outflows and mass-dependent timescale for the re-accretion of ejected gas. In previous work, we have shown that this scheme is able to correctly reproduce the evolution of the galaxy stellar mass function (GSMF) up to $z\sim3$. We contrast model predictions with both rest-frame Ultra-Violet (UV) and optical luminosity functions (LF), which are mostly sensible to the SFR and stellar mass, respectively. We show that GAEA is able to reproduce the shape and redshift evolution of both sets of LFs. We study the impact of dust on the predicted LFs and we find that the required level of dust attenuation is in qualitative agreement with recent estimates based on the UV continuum slope. The consistency between data and model predictions holds for the redshift evolution of the physical quantities well beyond the redshift range considered for the calibration of the original model. In particular, we show that GAEA is able to recover the evolution of the GSMF up to z$\sim$7 and the cosmic SFR density up to z$\sim$10.

Full text

1. Introduction

Since the introduction of the drop-out technique (see e.g. Steidel et al., 1996), the study of galaxy populations at increasingly higher redshift has provided fundamental contributions to our understanding of the first stages of structure formation in the Universe. The advent of space based observatories (like the Hubble Space Telescope - HST - and Spitzer) has allowed us to push these studies to the $3<z<6$ redshift range and beyond (e.g. Bouwens et al., 2016). Programs like the Great Observatories Origins Deep Survey (GOODS), the Cosmic Assembly Near-infrared Deep Extragalactic Legacy Survey (CANDELS) and the Hubble Ultra Deep Field (HUDF) provide excellent datasets to select high-z galaxy candidates (see e.g. Bouwens et al. 2015 and Finkelstein et al. 2015 and references herein). Early work in the field focus on the determination of the luminosity function (LF) in the rest-frame Ultra-Violet (UV), easily accessible through optical photometry (e.g. using the Advanced Camera for Surveys - ACS - onboard the HST). Since rest-frame UV bands provide information about the unobscured star formation rate (SFR) for detected sources, the resulting UV-LF can be used to study the evolution of the cosmic SFR and to assess the galaxy contribution to the reionization of the Universe (Robertson et al., 2013; Fontanot et al., 2014). Later instruments like the Infrared Array Camera (IRAC) on Spitzer and the Wide Field Camera 3 (WFC3) on HST have considerably widened the wavelength range available for high-z studies. Recently, several groups have used this multiwavelength information to estimate the LF at rest-frame optical wavelengths (Stefanon et al., 2016) and/or the galaxy stellar mass function (GSMF) at $z>4$ (González et al., 2011; Grazian et al., 2015). The rest-frame optical/near-infrared information complements the rest-frame UV data so that, when considered together, they represent a powerful tool to constrain the physical mechanisms shaping the early stages of galaxy evolution.

This wealth of data has been contrasted against theoretical models of galaxy formation (see e.g. Lo Faro et al., 2009; Lacey et al., 2011; Cai et al., 2014). In order to reproduce the observed LFs, all these studies favor a scenario where considerable dust extinction needs to be considered, up to the highest redshifts. Bouwens et al. (2014), however, argue that the evolution of the UV continuum slope of high-z galaxies selected from the CANDELS and HUDF fields requires that dust attenuation decreases with decreasing luminosities and increasing redshift (being considerably smaller at z$\sim$5-6 than at z$\sim$2-3). This implies that the UV-derived SFR density accounts for most of the cosmic SFR at z$>$4. Recent observations of 16 galaxies in the HUDF with the Atacama Large Millimeter Array (ALMA) by Dunlop et al. (2017) reinforce the conclusion that at z$\;\buildrel>\over{\sim}\;$4 most of the cosmic SFR is unobscured. Recent results from Rowan-Robinson et al. (2016) challenge, however, this scenario. They use Herschel 500$\mu$m counts and estimate SFR densities at 4$<$z$<$6 significantly higher than those expected from UV. It is worth stressing that in this redshift range SFR densities are still relatively uncertain, and are based on an handful of exceptional objects with individual SFR estimates $>10^{3}{\rm M}_{\odot}\ {\rm yr}^{-1}$. Indeed, work by Bourne et al. (2017) based on the SCUBA-2 Cosmology Legacy Survey, suggests a transition at z$\sim$4 from an (almost) unobscured early phase of galaxy formation to later epochs dominated by dust-obscured SFR, which the authors interpret as driven by the formation of the most massive galaxies.

The redshift evolution of galaxies below the knee of the GSMF has long been a problem for theoretical models of galaxy formation, which typically predict these objects to form too early (Fontanot et al., 2009; Weinmann et al., 2012; Hirschmann et al., 2012). A number of recent studies (Henriques et al., 2013; White et al., 2015; Hirschmann et al., 2016) point out that this problem can be alleviated by modifications of the adopted stellar feedback scheme. In detail, the most successful solutions invoke a combination of a strong ejective feedback (in the form of strong stellar-driven outflows) and a mass-dependent timescale for the re-accretion of the ejected gas onto dark matter haloes. Most, if not all, of these studies have focused on the $z<3$ GSMF, due to the more stringent constraints available. Lo Faro et al. (2009) show that the problem of the evolution of intermediate-to-low mass galaxies affects the predicted shape of the high-z UV-LFs.

In this Letter, we investigate the impact of strong stellar-driven outflows on galaxy properties beyond $z\sim 4$, taking advantage of new datasets and with the aim of presenting a coherent picture of galaxy evolution over the widest redshift range available.

2. Semi-analytic Model

In this letter, we consider predictions from the model for GAlaxy Evolution and Assembly (gaea - Hirschmann et al. 2016), which represents an evolution of the De Lucia & Blaizot (2007) code. The new model features significant improvements both in the treatment of chemical enrichment (the code accounts for the non instantaneous recycling of metals, gas and energy from asymptotic giant branch stars, Type Ia and Type II Supernovae - De Lucia et al. 2014) and in the modeling of stellar feedback. In particular, Hirschmann et al. (2016) compared gaea runs with different stellar feedback schemes to the observed evolution of the GSMF. In this Letter, we consider two of these schemes. The first one (“fiducial”), corresponds to the standard “energy-driven” scheme implemented in De Lucia et al. (2004) and De Lucia et al. (2014); the second one (H16F) corresponds to the “FIRE” stellar feedback implementation considered in Hirschmann et al. (2016). In the latter model, gas reheating is parametrized by using the fitting formulae discussed in Muratov et al. (2015), based on the “FIRE” set of hydrodynamical simulations (Hopkins et al., 2014). The same physical dependencies are assumed for the modeling of the rate of energy injection, while the ejected gas mass (outside the dark matter haloes) is estimated following energy conservation arguments as in Guo et al. (2011). Both the ejection and reheating efficiencies are treated as free parameters. Finally, we assume that the time-scale of gas re-incorporation scales with the halo mass as assumed in Henriques et al. (2013). Hirschmann et al. (2016) show that an improved modeling of both ejection and re-accretion is critical to reproduce the evolution of galaxies below the knee of the GSMF in the redshift range $0<z<3$, as well as the observed evolution of the galaxy mass- and gas-metallicity relations up to $z\sim 2$. The stellar feedback strength in the H16F prescriptions is assumed to increase with redshift so that outflows are stronger at higher redshifts. We stress that in this work, we analyse model predictions on redshift and stellar mass ranges that go well beyond both the original analysis by Muratov et al. (2015, z$\;\buildrel<\over{\sim}\;$4, 10^{10}\lower 2.15277pt\hbox{;\buildrel<\over{\sim};}M_{\star}/{\rm M}_{\odot}\lower 2.15277pt\hbox{;\buildrel<\over{\sim};}10^{12}), and the calibration in Hirschmann et al. (2016, limited at z$\;\buildrel<\over{\sim}\;$3). We note that the H16F prescription is not the only ejective feedback modeling able to reproduce the evolution of the GSMF in gaea, although being the only one able to reproduce simultaneously the observed evolution of the mass-metallicity relation. It is for this reason that we have elected the H16F prescription as our reference model in the following work, as well as in this study.

We couple gaea with dark matter halo merger trees extracted from the Millennium Simulation (MS Springel et al., 2005), a high resolution cosmological simulation of a $\Lambda$CDM concordance model, with parameters111Despite the values assumed for cosmological parameters are slightly different from the most recent determination (Planck Collaboration XVI, 2014), we do not expect this to affect significantly our conclusions (see e.g. Wang et al., 2008). assuming a WMAP1 cosmology (i.e. $\Omega_{\Lambda}=0.75$, $\Omega_{m}=0.25$, $\Omega_{b}=0.045$, $n=1$, $\sigma_{8}=0.9$, $H_{0}=73\,{\rm km/s/Mpc}$). In order to extend model predictions to lower masses and fainter luminosities we also consider runs based on the Millennium-II Simulation (MSII Boylan-Kolchin et al., 2009), which assumes the same cosmology, but a smaller cosmological volume and a 125 times better mass resolution. Throughout this Letter, we assume a Chabrier (2003) initial mass function and we use the the stellar population synthesis model from Bruzual & Charlot (2003).

3. High-z LFs

Fig. 1 presents the gaea predicted LFs in the rest-frame UV ($\sim 160{\rm nm}$) and optical ($\sim 900{\rm nm}$). Absolute UV magnitudes have been computed using a top-hat filter centered at 160 ${\rm nm}$ and 20 ${\rm nm}$ wide; optical magnitudes are computed in the $z_{ACS}$ filter. We compare our predictions with a variety of high-z data coming from the recent compilations of Bouwens et al. (2015), Finkelstein et al. (2015) and Stefanon et al. (2016). These data are all based on HST Legacy Fields, and cover a wide wavelength range from the optical to the Near-Infrared.

In order to estimate the effect of dust extinction in model galaxies, we use the same approach as in De Lucia & Blaizot (2007). Young stars in dense birth clouds suffer larger extinction than evolved stars in more diffuse cirrus. The age-dependent composite extinction curve is scaled with the column density of dust in the disc, assuming the dust mass is proportional to the metallicity. We further assume a “slab” geometry (Devriendt et al., 1999), to provide an estimate of the total dust attenuation in each source. This approach assumes a universal composite extinction curve for model galaxies, as well as a metallicity dependent normalization. Both are based on observations of $z=0$ galaxies. Therefore, large uncertainties linger at the redshifts considered in this study and might affect the comparison between theoretical predictions and observational data. In order to keep them under control, in Fig. 1 we consider both dust-extincted and unextincted LFs at the relevant wavelengths222The reference rest-frame wavelength varies slightly with redshift in data samples, due to the different filter sets used to select galaxies at different redshift.: the gray shaded (red hatched) area refers to predictions from the fiducial (H16F) model, with the lower (upper) envelope corresponding to the dust extincted (unextincted) LF. With respect to data, the unextincted fiducial model tends to overestimate the LFs at all redshifts, and a substantial amount of dust obscuration is needed to recover observations. The extincted H16F feedback model reproduces the number densities of UV and optical sources at z$\;\buildrel<\over{\sim}\;$5, but it underpredicts the bright-end of the LFs at higher redshifts. On the other hand, at $z>5$ the evolution of bright sources is better traced by the intrinsic LFs. In order to recover the overall evolution of the bright-end of the LF in the UV and optical bands, we thus have to assume a decreasing importance of dust attenuation, in qualitative agreement with Bouwens et al. (2009).

The magnitude range accessible with the MS is not wide enough to sample the faint end of the LFs. In order to study the shape of the LFs below the knee, we consider runs based on the MSII. We show predictions for the dust attenuated LFs from the Fiducial and H16F models (dotted black and dashed dark red lines, respectively) in Fig. 1. In both cases the convergence between the MSII predictions and those based on the MS is satisfactory. Dust attenuation is small for faint galaxies over the entire redshift range, in agreement with the analysis of Bouwens et al. (2009). In particular, the H16F run on the MSII reproduces the redshift evolution of the faint-end of the UV-LFs up to z$$\;\buildrel<\over{\sim}\;$$8. In contrast, the Fiducial model tends to over-predict the space density of faint UV galaxies. This result extends to the redshift range $4<z<7$ the evidence that strong stellar feedback represents a key ingredient to reproduce the observed evolution of the faint end of the LFs. At z$>$5 the available optical rest-frame data are not deep enough to firmly discriminate between the two schemes. Over this redshift range, the Fiducial model run on the MSII consistently predicts space densities for faint sources larger than the H16F feedback scheme.

4. Discussion

We compare in Fig. 2 the observed evolution of the GSMF and cosmic SFR density ($\rho_{\rm SFR}$), with model predictions. The latter have been convolved with an estimate of the observational errors, following Fontanot et al. (2009), i.e. a log-normal error distribution with amplitudes 0.25 and 0.3 for stellar masses and star formation rates, respectively. In the upper panel, we compare the predicted evolution of the GSMF at z$>$4 in the MS runs with observational determination based on stellar masses derived from spectral fitting techniques (González et al., 2011; Grazian et al., 2015) or from mass-to-light ratios in the optical (Stefanon et al., 2016). Results confirm and extend to higher-redshifts the conclusions from Hirschmann et al. (2016): models implementing the H16F feedback scheme are able to reproduce the shape and redshift evolution of the GSMF, while those based on the Fiducial scheme largely over-predict the number densities of galaxies below the knee of the mass function. The H16F feedback scheme is characterized by ejection rates that are larger than those assumed in the fiducial model (Hirschmann et al., 2016, see e.g. Fig 4). This implies that large amounts of reheated gas coupled with hot gas associated with dark matter haloes, are ejected in a reservoir that is is assumed to be unavailable for cooling. The observed evolution of the faint-end of the GSMF is then recovered by also assuming a dependence of gas re-accretion time-scale on halo mass. The agreement with data apparently worsens at the highest redshifts considered, but we stress that at z$\;\buildrel>\over{\sim}\;$6 model predictions are computed on the closest snapshots available for the MS and MSII (i.e. z=6.2 and z=7.2), that lie at a different redshift than the mean redshift of observed samples (z$\sim$6 and z$\sim$7). This mismatch could account for at least part of the disagreement.

In the bottom panel of Fig. 2, we consider the $\rho_{\rm SFR}$ evolution. In order to perform a meaningful comparison with high-z data from Bouwens et al. (2015) or Finkelstein et al. (2015) we only consider the contribution of model galaxies with $M_{\rm 160nm}<-17$, which roughly corresponds to $0.03L^{\star}_{\rm z=3}$. Given this faint integration limit, we consider in these panels predictions from the MSII runs. Once again, the run implementing the H16F feedback scheme provides an excellent agreement with observational determinations up to z$\sim$10, while the Fiducial model over-predicts the SFR density already at z$\sim$4. This confirms previous results suggesting that strong stellar feedback regulates the early evolution of $\rho_{\rm SFR}$ (see e.g. Vogelsberger et al., 2013, and references herein) We study the contribution of different galaxy populations (binned in stellar mass) to the cosmic SFR. Our results show that the main difference between the Fiducial and H16F feedback schemes is seen in the evolution of the smallest galaxies. The impact on $10^{11}<M_{\star}/{\rm M}_{\odot}<10^{10}$ galaxies is limited and mainly seen at z$\;\buildrel>\over{\sim}\;$4; more massive galaxies are those less affected by the different feedback schemes.

5. Summary

In this letter we contrast predictions from our semi-analytic model gaea with the latest constraints on the evolution of the high-z galaxies, using a combination of photometry (LFs in various bands) and derived physical properties. In particular, we consider predictions from our “fiducial” feedback scheme (that represents the standard “energy-driven” implementation from De Lucia et al. 2004) and from a feedback scheme based on the results of hydrodynamical simulations (H16F). Model predictions are compared to observational measurements out to the highest redshifts probed by state-of-the-art surveys, i.e. the edge of the epoch of reionization. Our results confirm and extend the conclusions in Hirschmann et al. (2016), i.e. they clearly show the need for strong stellar-driven outflows coupled with mass-dependent re-accretion timescales in order to correctly reproduce the evolution of the rest-frame UV and optical LFs over the redshift range 4$<$z$<$7. In addition, gaea runs implementing the H16F feedback agree well with the evolution of the GSMF up to the highest redshifts probed by the latest determinations (z$\sim$7) and with the cosmic SFR derived by Bouwens et al. (2015) (or Finkelstein et al., 2015) up to z$\sim$10.

We conclude that gaea is able to reproduce the overall evolution of the LFs, GSMF and cosmic SFR over the redshift range 0$<$z$<$10. It is worth stressing that the good level of agreement shown in this letter is obtained without any retuning of the feedback parameters (which have been calibrated by Hirschmann et al. 2016 against lower redshift observables). Our findings agree with recent results from the MUFASA simulation suite Davé et al. (2016), that implements a kinetic feedback scheme with scalings based on the parametrization provided in Muratov et al. (2015). They show that the MUFASA runs are able to reproduce the evolution of cosmic SFR and the GSMF at $z<4$ and the evolution of its low-mass end slope up to z$\sim$6 (where the simulated volume is too small to efficiently sample the high-mass end of the GSMF).

A critical point of the analysis presented here is the treatment of dust attenuation. In this letter, we consider both intrinsic and dust attenuated magnitudes. Despite the simplified dust model adopted, we show that the bright-ends of the corresponding LFs bracket the observed LFs. Our results are in qualitative agreement with the Bouwens et al. (2009) inferences based on the UV continuum of the sources, i.e. that dust attenuation should decrease at increasing redshifts, becoming negligible at $z>5$. We note that a significant level of dust attenuation is required by the Fiducial feedback scheme to match the shape of the observed LFs, in agreement with previous studies (Lo Faro et al., 2009; Lacey et al., 2011; Cai et al., 2014). Intrinsic and attenuated LFs converge below the knee of the LFs, again in qualitative agreement with the estimates by Bouwens et al. (2009) of a lower dust attenuation for fainter sources.

The predicted level of dust attenuation at high redshift thus represents a powerful discriminant between different feedback schemes. Therefore, further insight in our understanding of galaxy evolution at z$\;\buildrel>\over{\sim}\;$4 is tightly connected to a better description of this key aspect. Form a theoretical perspective, the implementation of a self-consistent dust treatment in theoretical models of galaxy evolution, able to follow dust production and destruction alongside with galaxy assembly (see e.g., the recent work by Popping et al. 2016) represents a promising avenue. In addition, dust provides an important channel for molecular hydrogen formation, and can therefore play a crucial role in regulating star formation at different cosmic epochs. This modeling is beyond the aims of this letter, and will be the subject of future work. Here, we just highlight that our conclusions will change dramatically if a substantial fraction of the SFR density is related to highly obscured objects, as some studies suggest at $3<z<6$ (Rowan-Robinson et al., 2016). In order to clearly assess the contribution of dusty sources to the SFR density at $z>4$, and to distinguish between different incarnations for stellar feedback, forthcoming or proposed facilities (like the James Webb Space Telescope - JWST

• or the SPace Infrared telescope for Cosmology and Astrophysics
• SPICA) will be of paramount importance.

Acknowledgments

Galaxy catalogs from our new gaea model (implementing the H16F feedback scheme) based on the Millennium merger trees will be publicly available at http://www.mpa-garching.mpg.de/millennium. FF and GDL acknowledge financial support from the MERAC foundation and from the PRIN INAF 2014 “Glittering kaleidoscopes in the sky: the multifaceted nature and role of Galaxy Clusters.” MH acknowledges financial support from the European Research Council via an Advanced Grant under grant agreement no. 321323 (NEOGAL).