Chemistry in \mbox{high-}z galaxies

The impact of chemistry on the structure of  galaxies

A. Pallottini, A. Ferrara, S. Bovino, L. Vallini, S. Gallerani,R. Maiolino, S. Salvadori.

Centro Fermi, Museo Storico della Fisica e Centro Studi e Ricerche “Enrico Fermi”, Piazza del Viminale 1, Roma, 00184, Italy
Cavendish Laboratory, University of Cambridge, 19 J. J. Thomson Ave., Cambridge CB3 0HE, UK
Kavli Institute for Cosmology, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK
Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126 Pisa, Italy
Kavli IPMU, The University of Tokyo, 5-1-5 Kashiwanoha, Kashiwa 277-8583, Japan
Hamburger Sternwarte, Universität Hamburg, Gojenbergsweg 112, D-21029 Hamburg, Germany
Nordita, KTH Royal Institute of Technology and Stockholm University, Roslagstullsbacken 23, SE-10691 Stockholm, Sweden
Dipartimento di Fisica e Astronomia, Universita’ di Firenze, Via G. Sansone 1, Sesto Fiorentino, Italy
INAF/Osservatorio Astrofisico di Arcetri, Largo E. Fermi 5, Firenze, Italy
GEPI, Observatoire de Paris, PSL Research University, CNRS, Place Jule Janssen 92190, Meudon, France
andrea.pallottini@centrofermi.it; ap926@mrao.cam.ac.uk
Abstract

To improve our understanding of  galaxies we study the impact of  chemistry on their evolution, morphology and observed properties. We compare two zoom-in high-resolution () simulations of prototypical galaxies at . The first, “Dahlia”, adopts an equilibrium model for  formation, while the second, “Althæa”, features an improved non-equilibrium chemistry network. The star formation rate (SFR) of the two galaxies is similar (within 50%), and increases with time reaching values close to at . They both have SFR-stellar mass relation consistent with observations, and a specific SFR of . The main differences arise in the gas properties. The non-equilibrium chemistry determines the H  transition to occur at densities , i.e. about 10 times larger than predicted by the equilibrium model used for Dahlia. As a result, Althæa features a more clumpy and fragmented morphology, in turn making SN feedback more effective. Also, because of the lower density and weaker feedback, Dahlia sits away from the Schmidt-Kennicutt relation; Althæa, instead nicely agrees with observations. The different gas properties result in widely different observables. Althæa outshines Dahlia by a factor of () in [C  ( ) line emission. Yet, Althæa is under-luminous with respect to the locally observed [C ]-SFR relation. Whether this relation does not apply at  or the line luminosity is reduced by CMB and metallicity effects remains as an open question.

keywords:
galaxies: high-redshift, formation, evolution, ISM – infrared: general – methods: numerical
pagerange: The impact of chemistry on the structure of  galaxiesBpubyear: 2017

1 Introduction

Understanding the properties of the interstellar medium (ISM) of primeval galaxies is a fundamental challenge of physical cosmology. The high sensitivity/spatial resolution allowed by current observations have dramatically improved our understanding of the ISM of local and moderate redshift () galaxies (Osterbrock, 1989; Stasińska, 2007; Pérez-Montero, 2017; Stanway, 2017). We now have a clearer picture of the gas phases and thermodynamics (Daddi et al., 2010a; Carilli & Walter, 2013), particularly for what concerns the molecular component, representing the stellar birth environment (Klessen & Glover, 2014; Krumholz, 2015).

For galaxies located in the Epoch of Reionization (EoR, ) optical/near infrared (IR) surveys have been very successful in their identification and characterization in terms of stellar mass and star formation rate (Dunlop, 2013; Madau & Dickinson, 2014; Bouwens et al., 2015). However, only recently we have started to probe the internal structure of such objects. With the advent of the Atacama Large Millimeter/Submillmeter Array (ALMA) it is now possible to access the far infrared (FIR) band at  with an unprecedented resolution and sensitivity. Excitingly, this enables for the first time studies of ISM energetics, structure and composition in such pristine objects.

Since C  is one of the major coolant of the ISM, [C ] ALMA detections (and upper limits) have so far mostly used this line for the above purposes (Maiolino et al., 2015; Willott et al., 2015; Capak et al., 2015) and to determine the sizes of early galaxies (Fujimoto et al., 2017). Line emission from different species (e.g. [O ]) have been used to derive the interstellar radiation field (ISRF) intensity (Inoue et al., 2016; Carniani et al., 2017), while continuum detections give us a measure of the dust content and properties (Watson et al., 2015; Laporte et al., 2017). Finally, some observations are beginning to resolve different ISM components and their dynamics by detecting spatial offsets and kinematic shifts between different emission lines, i.e. [C ] and optical-ultraviolet (UV) emission (Maiolino et al., 2015; Capak et al., 2015), [C ] and Ly (Pentericci et al., 2016; Bradac et al., 2017) and [C ] and [O (Carniani et al., 2017).

In spite of these progresses, several pressing questions remain unanswered. A partial list includes: (a) What is the chemical composition and thermodynamic state of the ISM in  galaxies? (b) How does the molecular gas turns into stars and regulate the evolution of these systems? (c) What are the optimal observational strategies to better constrain the properties of these primeval objects?

Theoretically, cosmological numerical simulations have been used to attack some of these problems. The key idea is to produce a coherent physical framework within which the observed properties can be understood. Such learning strategy is also of fundamental importance to devise efficient observations from current (e.g. HST/ALMA), planned (JWST) and proposed (SPICA) instruments. Before this strategy can be implemented, though, it is necessary to develop reliable numerical schemes catching all the relevant physical processes. While the overall performances of the most widely used schemes have been extensively benchmarked (Scannapieco et al., 2012; Kim et al., 2014, 2016), high-resolution simulations of galaxy formation introduce a new challenge: they are very sensitive to the implemented physical models, particularly those acting on small scales.

Among these, the role of feedback, i.e. how stars affect their own formation history via energy injection in the surrounding gas by supernova (SN) explosions, stellar winds and radiation, is far from being completely understood, despite considerable efforts to improve its modeling (Agertz & Kravtsov, 2015; Martizzi et al., 2015) and understand its consequences on  galaxy evolution (Ceverino et al., 2014; O’Shea et al., 2015; Barai et al., 2015; Pallottini et al., 2017; Fiacconi et al., 2017; Hopkins et al., 2017).

Additionally, we are still lacking a completely self-consistent treatment of radiation transfer. This is an area in which intensive work is ongoing in terms of faster numerical schemes (Wise et al., 2012; Rosdahl et al., 2015; Katz et al., 2016), or improved physical modelling (Petkova & Maio, 2012; Roskar et al., 2014; Maio et al., 2016).

A third aspect has received comparatively less attention so far in  galaxy formation studies, i.e. the implementation of adequate chemical networks. While various models have been proposed and tested (Krumholz et al., 2009; Bovino et al., 2016; Grassi et al., 2017), the galaxy-scale consequences of the different prescriptions are still largely unexplored (Tomassetti et al., 2015; Maio & Tescari, 2015; Smith et al., 2017). Besides, there is no clear consensus on a minimal set of physical ingredients required to produce reliable simulations.

The purpose of this paper is to analyze the impact of  chemistry on the internal structure of  galaxies. To this aim, we simulate two prototypical Lyman Break Galaxies (LBG) at , named “Dahlia” and “Althæa”, respectively. The two simulations differ for the  formation implementation, equilibrium vs. non-equilibrium. We show how chemistry has a strong impact on the observed properties of early galaxies.

The paper is organized as follows. In Sec. 2 we describe the two simulations highlighting common features (Sec. 2.1), separately discussing the different chemical models used for Dahlia (Sec. 2.2) and Althæa (Sec. 2.3). Results are presented as follow. First we perform a benchmark of the chemical models (Sec. 2.4), and compare the star formation and feedback history of the two galaxies (Sec. 3.1). Next, we characterize their differences in terms of morphology (Sec. 3.2), thermodynamical state of the ISM (Sec. 3.3), and predicted [C ] and  (Sec. 4) emission line properties. Our conclusions are summarized in Sec. 5.

2 Numerical simulations

To assess the impact of  chemistry on the internal structure of  galaxies, we compare two zoom-in simulations adopting different chemical models. Both simulations follow the evolution of a prototypical LBG galaxy hosted by a dark matter (DM) halo (virial radius ).

The first simulation has been presented in Pallottini et al. (2017, hereafter P17). The targeted galaxy (which includes also about 10 satellites) is called “Dahlia” (see also Gallerani et al., 2016, for analysis of its infall/outflow structure). In such previous work we showed that Dahlia’s specific SFR (sSFR) is in agreement with both analytical calculations (Behroozi et al., 2013), and with observations (González et al., 2010, see also Sec. 3.1).

In the second, new simulation we follow the evolution of “Althæa”, by using improved thermo-chemistry, but keeping everything else (initial conditions, resolution, star formation and feedback prescriptions) unchanged with respect to the Dahlia simulation. We describe the implementation of these common processes in the following Section. Next we describe separately the chemical model used for Dahlia (Sec. 2.2) and Althæa (Sec. 2.3).

2.1 Common physical models

Both simulations are performed with a customized version of the Adaptive Mesh Refinement (AMR) code ramses (Teyssier, 2002). Starting from cosmological IC111We assume cosmological parameters compatible with Planck results: CDM model with total matter, vacuum and baryonic densities in units of the critical density , , , Hubble constant with , spectral index , (Planck Collaboration et al., 2014). generated with music (Hahn & Abel, 2011), we zoom-in the DM halo hosting the targeted galaxy. The total simulation volume is that is evolved with a base grid with 8 levels (gas mass ); the zoom-in region has a volume of and is resolved with 3 additional level of refinement, thus yielding a gas mass resolution of . In such region, we allow for 6 additional level of refinement, that allow to follow the evolution of the gas down to scales of at , i.e. the refined cells have mass and size typical of Galactic molecular clouds (MC, e.g. Federrath & Klessen, 2013). The refinement is performed with a Lagrangian mass threshold-based criterion. , i.e. a cell is refined if its total (DM+baryonic) mass exceed the the mass resolution by a factor 8.

Metallicity () is followed as the sum of heavy elements, assumed to have solar abundance ratios (Asplund et al., 2009). We impose an initial metallicity floor since at our resolution is still insufficient to catch the metal enrichment by the first stars (e.g. O’Shea et al., 2015). Such floor is compatible with the metallicity found at  in cosmological simulations for diffuse enriched gas (Davé et al., 2011; Pallottini et al., 2014a; Maio & Tescari, 2015); it only marginally affects the gas cooling time.

Dust evolution is not explicitly tracked during simulations. However, we make the simple assumption that the dust-to-gas mass ratio scales with metallicity, i.e. , where for the Milky Way (MW, e.g. Hirashita & Ferrara 2002; Asano et al. 2013.

2.1.1 Star formation

Stars form according to a linearly -dependent Schmidt-Kennicutt relation (Schmidt, 1959; Kennicutt, 1998) i.e.

(1)

where is the local SF rate density, the SF efficiency, the  mass fraction, and is density of the gas of mean molecular weight . Eq. 1 is solved stochastically, by drawing the mass of the new star particles from a Poisson distribution (Rasera & Teyssier, 2006; Dubois & Teyssier, 2008; Pallottini et al., 2014a). In detail, in a star formation event we create a star particle with mass , with an integer drawn from

(2)

where the mean of the Poisson distribution is

(3)

with the simulation time step. For numerical stability, no more than half of the cell mass is allowed to turn into stars. Since we prevent formation of star particle with mass less then , cells with density less then (for ) are not allowed to form stars.

We set , in accordance with the average values inferred from MC observations (Murray, 2011, see also Agertz et al. 2013); depends on the adopted thermo-chemical model, as described later in Sec. 2.2 and Sec. 2.3.

2.1.2 Feedback

Similarly to Kim et al. (2014), we account for stellar energy inputs and chemical yields that depend both on time and stellar populations by using starburst99 (Leitherer et al., 1999). Stellar tracks are taken from the padova (Bertelli et al., 1994) library with stellar metallicities in the range , and we assume a Kroupa (2001) initial mass function. Stellar feedback includes SNs, winds from massive stars and radiation pressure (Agertz et al., 2013). We model the thermal and turbulent energy content of the gas according to the prescriptions by Agertz & Kravtsov (2015). The turbulent (or non-thermal) energy is dissipated as (Teyssier et al., 2013, see eq. 2), where, following Mac Low (1999), the dissipation time scale can be written as

(4)

where is the turbulent velocity dispersion. Adopting the SN blastwave models and OB/AGB stellar winds from Ostriker & McKee (1988) and Weaver et al. (1977), respectively, we account for the dissipation of energy in MCs as detailed in Sec. 2.4 and App. A of P17.

2.2 Dahlia: equilibrium thermo-chemistry

In the Dahlia simulation we compute by adopting the KTM09 analytical prescription (Krumholz et al., 2008, 2009; McKee & Krumholz, 2010). In KTM09, the  abundance is derived by modelling the radiative transfer on an idealized MC and by assuming equilibrium between  formation on dust grains and dissociation rates. For each gas cell, can then be written as a function of , and hydrogen column density (). By further assuming pressure equilibrium between CNM and WNM (Krumholz et al., 2009), turns out to be independent on the intensity of the ISRF, and can be written as

(5a)
with
(5b)
(5c)

and where is the Heaviside function; is the dust UV optical depth and it can be calculated by linearly rescaling the MW value,

(6)

In Dahlia cooling/heating rates are computed using grackle 2.1222https://grackle.readthedocs.org/ (Bryan et al., 2014). We use a H and He primordial network, and tabulated metal cooling/photo-heating rates from cloudy (Ferland et al., 2013). Inverse Compton cooling is also present, and we consider heating from a redshift-dependent ionizing UV background (UVB, Haardt & Madau, 2012). Since  is not explicitly included in the network, we do not include the corresponding cooling contribution.

2.3 Althæa: non-equilibrium thermo-chemistry

In Althæa we implement a non-equilibrium chemical network by using krome333https://bitbucket.org/tgrassi/krome (Grassi et al., 2014). Given a set of species and their reactions, krome can generate the code needed to solve the system of coupled ordinary differential equations that describe the gas thermo-chemical evolution.

2.3.1 Chemical network

Similarly to Bovino et al. (2016, hereafter B16), our network includes H, H, H, He, He, He, H, H and electrons. Metal species are not followed individually in the network, as for instance done in Model IV from B16; therefore, we use an equilibrium metal line cooling calculated via cloudy tables444As a caveat, we point out that there is a formal inconsistency in the modelling. Metal line cooling tables are usually calculated with cloudy by assuming a Haardt & Madau (2012) UV background, while the ISRF spectral energy density we adopt is MW-like. To remove such inconsistency one should explicitly track metal species, adopt a non-equilibrium metal line cooling and include radiative transfer. As noted in B16 (see their Fig. 16), using non-equilibrium metal line cooling can typically change the cooling function by a factor . This will be addressed in future work.. The adopted network contains a total of 37 reactions, including photo-chemistry (Sec. 2.3.2), dust processes (Sec. 2.3.3) and cosmic rays (CR, Sec. 2.3.4). The reactions, their rates, and corresponding references are listed in App. B of B16: specifically we use reactions from 1 to 31 (Tab. B.1 in B16), 53, 54, and from 58 to 61 (Tab. B.2 in B16).

2.3.2 Photo-chemistry

Photo-chemistry cross sections are taken from Verner & Ferland (1996) and by using the SWRI555http://phidrates.space.swri.edu. and Leiden666http://home.strw.leidenuniv.nl/~ewine/photo/. databases. In the present simulation, the ISRF is not evolved self-consistently and it is approximated as follows. For the spectral energy density (SED), we assume a MW like spectrum (Black, 1987; Draine, 1978), and we specify the SED using 10 energy bins from eV to eV. Beyond 13.6 eV the flux drops to zero, i.e. we do not include ionizing radiation.

We consider a spatially uniform ISRF whose intensity is rescaled with the SFR such that , where is the far UV (FUV) flux in the Habing band () normalized to the average MW value (Habing, 1968). Because of their sub-kpc sizes (Shibuya et al. 2015, Fujimoto et al. 2017) high values are expected in typical LBG at , as inferred also by Carniani et al. (2017). A similar situation is seen in some local dwarf galaxies (Cormier et al., 2015) that are generally considered as local counterparts of  galaxies. It is worth noting that the spatial variation of is very small in the MW, with an r.m.s. value (Habing, 1968; Wolfire et al., 2003). Nonetheless, spatial fluctuations of the ISRF, if present, might play some role in the evolution of  galaxies (e.g. Katz et al., 2016). We will analyze this effect in future work.

On top of the ISRF, we consider the cosmic microwave background (CMB), that effectively sets a temperature floor for the gas. Additionally, we neglect the cosmic UVB, since the typical ISM densities are sufficiently large to ensure an efficient self-shielding (e.g. Gnedin, 2010). For example, Rahmati et al. (2013) have shown that at the hydrogen ionization due to the UVB is negligible for , the typical density of diffuse ISM. The self-shielding of  to photo-dissociation is accounted by using the Richings et al. (2014) prescription777The self-shielding formulation by Richings et al. (2014) does not account for a directional dependence as done in more computationally costly models (Hartwig et al., 2015)., thus in each gas cell the shielding can be expressed as an analytical function of its  column density, temperature and turbulence (cfr. with Wolcott-Green et al., 2011).

2.3.3 Dust processes

As for Dahlia, the dust mass is proportional to the metal mass. Here we also specify the dust size distribution to be the one appropriate for  galaxies, the Small Magellanic Cloud one, following Weingartner & Draine (2001). Dust grains can affect the chemistry through cooling888Dust cooling is not included in the current model, as it gives only a minor contribution for , i.e. see Fig. 3 in B16. (Hollenbach & McKee, 1979), photoelectric-heating (Bakes & Tielens, 1994), and by mediating the formation of molecules (Cazaux & Spaans, 2009). In particular, the formation rate of  on dust grains is approximated following Jura (1975)

(7)

where is the hydrogen density. Note that for this dust channel is dominant with respect to gas-phase formation (e.g. reactions 6–7 and 9–10 B.1 in B16).

2.3.4 Cosmic rays

CR ionization can become important in regions shielded from radiation, like MC interiors. We assume a CR hydrogen ionization rate SFR (Valle et al., 2002) and normalized to the MW value (Webber, 1998):

(8)

The rate includes the flux of CR and secondary electrons (Richings et al., 2014). In the network, CR ionizations are proportional to and to coupling constants that depend on the specific ions; such couplings are taken from the kida database (Wakelam et al., 2012). Additionally we account for Coulomb heating, by assuming that every CR ionization releases an energy999For a more accurate treatment of Coulomb heating refer to Glassgold et al. (2012). of eV.

2.3.5 Initial abundances of the species

Finally, following Galli & Palla (1998), we calculate IC for the various species by accounting for the primordial chemistry101010For a possible implementation of the Galli & Palla (1998) chemical network see the “earlyUniverse” test contained in krome. at , for a density and temperature evolution corresponding to gas at the mean cosmic density.

2.4 Benchmark of  formation models

As a benchmark for our simulations, we compare the formation of  in different physical environments. For the Dahlia KTM09 model we compute from eq. 5 as a function of and . We choose an expression for resulting from the mass threshold-based AMR refinement criterion for which . We restate that the equilibrium KTM09 model is independent on and the gas temperature .

For the Althæa B16 model we use krome to perform single-zone tests varying and . In this case we assume an initial temperature111111The initial temperature corresponds to the virial temperature of the first star-forming halos present in the zoomed region. The results depend very weakly on this assumption. , and we let the gas patch evolve at constant density until thermo-chemical equilibrium is reached. This typically takes 100 Myr.

The comparison between the two models is shown in Fig. 1 as a function of for different metallicities. For and ,  formation is hindered in B16 with respect to KTM09, i.e. higher are needed to reach similar fractions. For and the two models are roughly in agreement: this is expected since KTM09 is calibrated on the MW environment. Finally, for and at (the metallicity floor in our simulation set) the  formation in the B16 model is strongly suppressed, e.g. for . Note that these fractions are comparable to the ones expected for  formation in a pristine environment where  formation proceeds via gas-phase reactions.

As noted in P17, Dahlia’s star formation (SF) model (eqs 1 and 5) is roughly equivalent to a density threshold criterion with metallicity-dependent critical density . Physically this corresponds to the density at which (see also Agertz et al. 2013). Thus, Fig. 1 quantifies the density threshold required to spawn stars in the simulation.

Dahlia forms stars in gas with and at a rate of about at . If Althæa has a similar SFR history (this is checked a posteriori, see Fig. 2), the resulting metallicity and ISRF intensity () should also be similar. Then, by inspecting Fig. 1 (middle-left panel) one can conclude that Althæa forms stars in much denser environments where for . As noted in Hopkins et al. (2013), although variations in the density threshold lead to similar total SFR, they might severely affect the galaxy morphology. We will return to this point in the next Section.

We remind that in both simulation we use the cell radius to calculate column density, that are used e.g. to calculate the gas self-shielding. This is done to mainly to ensure a fair comparison between the two simulations. In other simulations, e.g. MC illuminated by an external radiation field, the prescriptions adopted accounts for the contribution of column density from nearby cells, i.e. by using Jeans or Sobolev-like length (see e.g. Hartwig et al. 2015 for a comparison between different prescriptions). However in our simulation we expect stars to be very close or embedded in potential star forming regions. Using the contribution to the column density from the surrounding gas would then overestimate the self-shielding effect. Such modelling uncertainty would be solved by including radiative transfer in the simulation. However, we note that at z=6 the radius of our cells as a function of density can be approximated as , while the jeans length is . Thus for typical values found for the molecular gas in Althaea ( and , see later Fig. 8), the two prescriptions gives similar results, i.e. and .

Figure 1: Benchmark of the formation of  for the model used in Dahlia (KTM09, Sec. 2.2) and in Althæa (B16, Sec. 2.3). In each panel we plot the  mass fraction as a function of density (), with different panels showing the results for different metallicities (). In each panels the dashed grey line indicates the KTM09 model, while the B16 models are plotted with solid lines, with different colours indicating a different impinging ISRF flux (). In the upper axis we indicate the free-fall times () corresponding to .

3 Results

We now turn to a detailed analysis of the two zoomed galaxies, Dahlia and Althæa121212We refrain from the analysis of the satellite population of the two galaxies due to the oversimplifying assumption of a spatially uniform ISRF artificially suppressing star formation in environments with metallicity close to the floor value .. We start by studying the star formation and the build-up of the stellar mass from to (Sec. 3.1). We then specialize at to inspect the galaxy morphology (Sec. 3.2), the ISM multiphase structure (Sec. 3.3) and the predicted observable properties (Sec. 4). An overall summary of the properties of the two galaxies is given in Tab. 1.

3.1 Star formation history

Figure 2: Star formation rate (SFR) as a function of galaxy age () for Dahlia (black line and hatched region) and Althæa (orange line and transparent region). Also shown (grey dashed line) is an analytical approximation (within a factor 2 for both galaxies) to the average SFR trend. The redshift () corresponding to is plotted on the upper axis, and note that corresponds to the first stellar formation event in Dahlia, and the plotted SFRs are averaged over .
Figure 3: SFR vs stellar mass () for Althæa (circles) and Dahlia (squares), with symbols coloured accordingly to the age . With crosses we overplot SFR and inferred from 27 galaxies observed at by J16. Following J16 analysis14, galaxies identified as young and old are plotted in blue and red, respectively. To guide the eye, the linear correlation between the data sets are also shown with a dashed lines. See the text for more details.

In Fig. 2 we plot the SFR history as a function of “galaxy age” () for Dahlia and Althæa; marks the first star formation event in Dahlia131313Note that even with the same modelling and IC, differences in the SFR may arise as a result of stochasticity in the star formation prescription eq. 1. Such differences vanish once the SFR is averaged on timescales longer than the typical free-fall time of the star forming gas.. For both Dahlia and Althæa the SFR has an increasing trend which can be approximated with good accuracy (within a factor of 2) as . However, on average the SFR in Dahlia is larger by a factor when averaged over the entire SFR history (). Thus, in spite of very different chemical prescriptions, the SFR in the two galaxies shows very little variation. Stated differently, the higher critical density for star formation arising from non-equilibrium chemistry does not alter significantly the rate at which stars form, as already noticed in Sec. 2.4. This also entails a comparable metallicity, and we note that in both galaxy most of the metal mass is locked in stars (see Tab. 1), as they are typically formed from the most enriched regions.

It is interesting to check the evolutionary paths of Dahlia and Althæa (Fig. 3) in the standard SFR vs. stellar mass () diagram, and compare them with data141414SFR and have been derived by assuming an exponentially increasing SFR, consistent with the history of both our simulated galaxies (Fig. 2). inferred from observations of 27 Lyman Alpha Emitters (LAE) and LBGs (Jiang et al., 2016, hereafter J16). By using multi-band data, precise redshift determinations, and an estimate of nebular emission from Ly, J16 were able to distinguish between a young () and an old () subsample. Each subsample exhibits a linear correlation in , albeit with a different normalization: the young (old) subsample has a .

The SFR vs stellar mass of our simulated galaxies for () is fairly consistent with the young subsample relation (keeping in mind stochasticity effects at low stellar masses). At later evolutionary stages ( or ), Dahlia and Althæa nicely shift to the lower sSFR values characterizing the old J16 subsample data. This shift must be understood as a result of increasing stellar feedback: as galaxies grow, the larger energy input from the accumulated stellar populations hinders subsequent SFR events. Note that at late times ( Myr), when , the sSFRs of Dahlia and Althæa are in agreement with analytical results by Behroozi et al. (2013), and with observations by González et al. (2010).

Figure 4: Ratio of mechanical () and radiative () energy deposition rates by stars as a function of galaxy age () for Dahlia (black line/hatched area) and Althæa (orange/transparent). Dashed lines indicate the time-averaged mean of the ratios for each galaxy. To guide the eye we plot the unity value (dotted grey line). Similar to Fig. 2 the ratios are averaged over . The upper horizontal axis indicates redshift.
Property Symbol Dahlia Althæa [units]
Star formation rate
Specific SFR
Stellar mass
Metal mass in stars
Gas mass
H mass
Metal mass
Disk radius
Disk scale height
Gas density
H density
Metallicity
Gas surface density
Star formation surface density
Luminosity [C 
Luminosity  
Table 1: Physical properties of Dahlia and Althæa at . The values refer to gas and stars within from the galaxy center (similar to the field of view in Fig.s 5 and 11). The effective radius, , and gas scale height, , are calculated from the principal component analysis of the density field. Values for , , , , and represent mass-weighted averages.

As feedback clearly plays a major role in the overall evolution of early galaxies, we turn to a more in-depth analysis of its energetics. This can be quantified in terms of the stellar energy deposition rates in mechanical (SN explosions + OB/AGB winds151515On average OB/AGB winds account only for of the SN power., ) and radiative () forms. These are shown as a function of time in Fig. 4. The ratio shows short-term () fluctuations corresponding to coherent burst of star formation/SN activity.

Barring this time modulation, on average the mechanical/radiative energy ratio increases up to , when it suddenly drops and reaches an equilibrium value. This implies that radiation pressure dominates the energy input; consequently it represents the major factor in quenching star formation. While this is true throughout the evolution, it becomes even more evident after , when the first stellar populations with enter the AGB phase. At that time, winds from AGBs enrich their surroundings with metals and dust. As dust produced by AGBs remains more confined than SN dust around the production region, it provides a higher opacity, thus boosting radiation pressure via a more efficient dust-gas coupling (see also P17).

For Dahlia the radiative energy input rate is about 20 times larger than the mechanical one, while for Althæa such ratio is on average times higher, although larger fluctuations are present. The latter are caused by the occurrence of more frequent and powerful bursts of SN events in Althæa. Why does this happen?

The answer has to do with the different gas morphology. As already noted discussing Fig. 1, the higher critical density for star formation imposed by non-equilibrium chemistry has a number of consequences: (a) each formation event can produce a star cluster with an higher mass; (b) star formation is more likely hosted in isolated high density clumps (see later, particularly Fig. 6); (c) in a clumpier disk, SN explosions can easily break into more diffuse regions. The combination of (a) and (b) increases the probability of spatially coherent explosions having a stronger impact on the surrounding gas; due to (c), the blastwaves suffer highly reduced radiative losses (Gatto et al., 2015), and affect larger volumes. Similar effects have been also noted in the context of single giant MCs (), where unless the SNs explode coherently, their energy is quickly radiated away because of the very high gas densities (Rey-Raposo et al., 2017). For the reminder of the work we focus on , when the galaxies have an age of .

3.2 Galaxy morphology

Figure 5: (Caption next page.)
Figure 5: (Previous page.) Face-on maps16 of Dahlia (left panels) and Althæa (right) at age (). Shown are line-of-sigh mass weighted average of the gas density (upper panels), temperature (middle), and  density (lower) fields with amplitude given by the colorbar. The maps are on a side.
Figure 6: Morphological comparison of the molecular gas at . In the four panels we plot the Minkowsky functionals () of the  density field (). Functionals are plotted with black line and hatched regions for Dahlia, with orange line and transparent region for Althæa. Note that Minkowsky functionals are indicated in comoving units. For detail on the calculation of the Minkowsky functional see App. A (in particular see Fig. 13).

Dahlia and Althæa sit at the centre of a cosmic web knot and accrete mass from the intergalactic medium mainly via 3 filaments of length . In both simulations, the large scale structure is similar, and we refer the reader to Sec. 3.1 of P17 for its analysis. Differences between the simulation are expected to arise on the ISM scale, whose structure is visible on scales. In Fig. 5 we show the gas density, temperature, and  density () fields for Dahlia and Althæa. The map161616The maps of this work are rendered by using a customized version of pymses (Labadens et al., 2012), a visualization software that implements optimized techniques for the AMR of ramses. centers coincide with Dahlia’s stellar center-of-mass.

3.2.1 Overview

Qualitatively, both galaxies show a clearly defined, albeit somewhat perturbed, spiral disk of radius , embedded in a lower density () medium. However the mean disk gas density for Dahlia is , while for Althæa (see Tab. 1). The temperature structure shows fewer differences, i.e. the inner disk is slightly hotter for Dahlia () than for Althæa (), which features instead slightly more abundant and extended pockets of shock-heated gas (). Such high- regions are produced by both accretion shocks and SN explosions. In both cases the typical  density is the same, i.e. , however, with respect to Dahlia, Althæa shows a slightly smaller disk, that also seems more clumpy.

To summarize, the galaxies differ by an order of magnitude in atomic density, but have the same molecular density. In spite of this difference, the SFR are roughly similar. This can be explained as follows. To first order, in our model , where is the galaxy volume (Tab. 1). It follows that the larger density is largely compensated by the smaller Althæa volume.

3.2.2 In-depth analysis

Fig. 5 visually illustrates the morphological differences between the two galaxies. The gas in Althæa appears clumpier than in Dahlia. To quantify this statement we start by introducing the  clumping factor on the smoothing scale , which is defined as171717To calculate the clumping factor, first we construct the 3D unigrid cube of the  mass field, then we smooth it with a Gaussian kernel of scale and finally we calculate the mass-weighted average and variance of the smoothed  density field.

(9)

For Dahlia decreases from to 10 going from 30 pc to 1 kpc, while for Althæa is times larger on all scales.

A more in-depth analysis can be performed using the Minkowsky functionals (Schmalzing & Gorski, 1998; Gleser et al., 2006; Yoshiura et al., 2017, App. A) which can give a complete description of the molecular gas morphological structure. For a 3-dimensional field, 4 independent Minkowsky functionals can be defined. Each of the functionals, characterizes a different morphological property of the excursion set with  density : gives the volume filling factor, measures the total enclosed surface, is the mean curvature, quantifying the sphericity/concavity of the set, and estimates the Euler characteristic (i.e. multiple isolated components vs. a single connected one. Appendix A gives more rigorous definitions with an illustrative application (Fig. 13).

In Fig. 6 we plot the Minkowsky functionals () calculated for the  density field for Dahlia and Althæa. The functional analysis shows that Althæa is more compact, i.e. for each value Dahlia’s excursion set volume is larger and it plummets rapidly at large densities. On the other hand, the set surface of Althæa is larger by about a factor of 5, implying that this galaxy is fragmented into multiple, disconnected components. This is confirmed also by Althæa’s larger () Euler characteristic measure, , an indication of the prevalence of isolated structures. This feature becomes more evident towards larger densities, as expected if  is concentrated in molecular clouds181818 values at in Dahlia are mainly due to the presence of the 3 satellites/clumps outside the disk..

Further, in Dahlia most of the molecular gas resides in connected () disk regions, with a concave shape (). For Althæa there is a transition: for the gas has a concave (), disjointed (), filamentary structure, while for the galaxy is composed by spherical clumps ().

3.3 ISM thermodynamics

Figure 7: Equation of state (EOS) of the gas within for Dahlia (left panel) and Althæa (right panel) at (). EOS are shown as mass-weighted probability distribution function (PDF) in the density-temperature () plane, as specified by the colorbar. For both galaxies, the EOS projection on the () axis is additionally shown as an horizontal (vertical) inset. The 2D EOS are normalized such that the integral on the plane is unity; the projected EOS are normalized such that the sum of the bins is equal to .
Figure 8: EOS of the molecular () gas for Dahlia (left panel) and Althæa (right panel) i.e. the  mass-weighted PDF in the plane. Notation is similar to Fig. 7, albeit a different region of plane is shown.

The thermodynamical state of the ISM can be analyzed by studying the probability distribution function (PDF) of the gas in the density-temperature plane, i.e. the equation of state (EOS). In Fig. 7 we plot the mass-weighted EOS for Dahlia and Althæa at . We include gas within , or , from the galaxy center.

From the EOS we can see that in both galaxies of the gas in a photoionized state (), that in Dahlia is induced by the Haardt & Madau (2012) UVB, while in Althæa is mainly due to photo-electric heating on dust grains illuminated from the uniform ISRF of intensity . Only of the gas is in a hot K component produced by accretion shocks and SN explosions. A relatively minor difference descends from Althæa’s more effective mechanical feedback, already noted when discussing Fig. 4: small pockets of freshly produced very hot () and diffuse () gas are twice more abundant in Althæa, as it can be appreciated from a visual inspection of the temperature maps in Fig. 5.

Fig. 7 (in particular compare the upper horizontal panels) shows that the density PDF is remarkably different in the two galaxies. In Dahlia the distribution peaks at ; Althæa instead features a bi-modal PDF with a second, similar amplitude peak at . This entails the fact that the dense gas is about 2 times more abundant in the latter system. In addition, the very dense gas (), only present in Althæa, can cool to temperatures of 30 K, not too far from the CMB one.

The high-density part of the PDF is worth some more insight as it describes the gas that ultimately regulates star formation. This gas is largely in molecular form, and accounts (see Tab. 1) for () of the total gas mass in Dahlia (Althæa). Its  density-weighted distribution in the plane is reported in Fig. 8. On average, the  gas in Dahlia is 10 times less dense than in Althæa as a result of the new non-equilibrium prescription requiring higher gas densities to reach the same fraction; at the same time the warm ()  fraction drops from (Dahlia) to an almost negligible value. Clearly, the warm component was a spurious result as (a)  cooling was not included, and (b) was considered to be independent of gas temperature (see eq.s 5). Note that in Althæa traces of warm  are only found at large densities, in virtually metal-free gas in which  production must proceed via much less efficient gas-phase reactions rather than on dust surfaces. This tiny fraction of molecular gas can survive only if densities large enough to provide a sufficient  self shielding against photodissociation are present.

Finally, the sharp EOS cutoff at in Dahlia is caused by the density-threshold behavior mimicked by the enforced chemical equilibrium: above (Sec. 2.2) the gas is rapidly turned into stars. This spurious effect disappears in Althæa, implementing a full time-dependent chemical network.

4 Observational properties

As we already mentioned, the strongest impact of different chemistry implementations is on the gas properties, and consequently on ISM-related observables. In the following, we highlight the most important among these aspects.

Figure 9: Comparison of the observed and simulated Schmidt-Kennicutt relation expressed in terms of - . Observations are taken from single MCs (Heiderman et al., 2010; Lada et al., 2010), local unresolved galaxies (Kennicutt, 1998), and moderate redshift unresolved galaxies (Bouché et al., 2007; Daddi et al., 2010a, b; Tacconi et al., 2010; Genzel et al., 2010); the correlation (dispersion) for the observation found by Krumholz et al. (2012, see the text for details) is plotted with a black dashed line (grey shaded region). Dahlia and Althæa averaged value are plotted with black and orange stars, respectively (see Fig. 10 for the complete distribution in the simulated galaxies).
Figure 10: Schmidt-Kennicutt relation in Dahlia (left panel) and Althæa (right panel) at (). The relation is plotted using the  mass weighted PDF of the instantaneous SFR surface density versus the total gas surface density (). On both panels with dashed grey lines we overplot the relation observed from Kennicutt & Evans (2012), i.e. , for several normalizations that written inline. Otherwise notation is similar to Fig.s 8 and 7.
Figure 11: Synthetic emission maps16 of the simulated galaxies Dahlia (left panels) and Althæa (right panels) at age (). Integrated surface brightness of [C ] () and  () are shown in the upper and lower panels, respectively. The field of view is the same as in Fig. 5.
Figure 12: The [C ]-SFR relation. Shown are Althæa (orange star) and Dahlia (black) at 700 Myr or ; the errors refers to r.m.s. variation in the last . Lines refer to results from the V15 model: constant metallicity models with (solid black), (solid orange), (pink dashed), and a model with mean + density-metallicity relation extracted from cosmological simulations (Pallottini et al., 2014b, blue dot-dashed). Data for local dwarf galaxies (De Looze et al., 2014) are plotted with little circles and the grey hatched region gives the mean and r.m.s. variation in the sample. For high galaxies, detections (upper-limits) are plotted with filled (empty) symbols, according to the inset legend. The high sample include individual galaxies as BDF-3299 (Maiolino et al., 2015; Carniani et al., 2017), HCM6A (Kanekar et al., 2013), Himiko (Ouchi et al., 2013; Ota et al., 2014), IOK-1 (Ota et al., 2014), and data from Capak et al. (2015, ), Willott et al. (2015, ), Schaerer et al. (2015, ), (Pentericci et al., 2016, ), González-López et al. (2014, ), and lensed galaxies from Knudsen et al. (2016) and Bradac et al. (2017).

4.1 Schmidt-Kennicutt relation

We start by analyzing the classical Schmidt-Kennicutt (SK) relation. This comparison should be interpreted as a consistency check of the balance between SF and feedback, since in the model we assume a SFR law that mimics a SK relation (eq. 1).

The SK relation, in its most modern (Krumholz et al., 2012) formulation, links the SFR () and total gas () surface density per unit free-fall time, . The proportionality constant, often referred to as the efficiency per free-fall time following eq. 1, is simply . Experimentally, Krumholz et al. (2012) find (see Krumholz 2015 for a complete review on the subject). This result is supported also by a larger set of observations including single MCs (Heiderman et al., 2010; Lada et al., 2010), local unresolved galaxies (Kennicutt, 1998), and moderate redshift, unresolved galaxies (Bouché et al., 2007; Daddi et al., 2010a, b; Tacconi et al., 2010; Genzel et al., 2010). The SK relation is shown in Fig. 9, along with the location of Dahlia and Althæa at .

Dahlia appears to be over-forming stars with respect to its gas mass, and therefore it is located about above the KS relation. As Althæa needs about 10 times higher density to sustain the same SFR, its location is closer to expectations from the SK. We have checked that the agreement is even better if we use only data relative to MC complexes (e.g. Heiderman et al., 2010; Murray, 2011).

Dahlia’s is similar to the analog values found by Semenov et al. (2016), who compute such efficiency using a turbulent eddy approach (Padoan et al., 2012), with no notion of molecular hydrogen fraction. The difference is that Dahlia misses the high density gas. Althæa instead matches both the and the amount of high density gas found by Semenov et al. (2016). Also, its relation is consistent with Torrey et al. (2017), who use a star formation recipe involving self-gravitating gas with a local SK  dependent relation.

From our simulations it is also possible to perform a cell-by-cell analysis of the SK relation (Fig. 10). As expected, the results show the presence of a consistent spread in the local efficiency values which, however, has a different origin for Dahlia and Althæa. While in the former the variation is mostly due to a different enrichment level affecting  abundance (eq. 5), for Althæa the spread is larger because it results also from the individual evolutionary histories of the cells.

As noted by Rosdahl et al. (2017), for galaxy simulations with a SF model based on SK-like relation (eq. 1), the resulting depends on how the feedback is implemented. However, here we show that Althæa has a lower in spite of the fact that it implements exactly the same feedback prescription as Dahlia. The latter is qualitatively similar to a delayed cooling scheme used by Rosdahl et al. (2017) and others (Stinson et al., 2006; Teyssier et al., 2013). The lower efficiency is a consequence of chemistry. As under non-equilibrium conditions the gas must be denser to form , the ISM becomes more clumpy (Fig. 6). These clumps can form massive clusters of OB stars which, acting coherently, yield stronger feedback and may disrupt completely the star forming site.

4.2 Far and mid infrared emission

A meaningful way to compare the two galaxies is to predict their C  and  line emission, that can be observable at  with ALMA, and possibly with SPICA (Spinoglio & et al., 2017; Egami & et al., 2017, in preparation), respectively. Similarly to P17, we use a modified version of the [C ] emission model from Vallini et al. (2015, hereafter V15). Such model is based on temperature, density and metallicity grids built using cloudy (Ferland et al., 2013), as detailed in App. B. In Fig. 11 we plot the [C  and   surface brightness maps (); the field of view is the same as in Fig. 5.

4.2.1 Far infrared emission

Let us analyze first the C  emission. Dahlia has a [C ] luminosity of which is about 7 times smaller than Althæa, i.e. . Fig. 11 shows that the surface brightness morphology in the two galaxies is similar. Dahlia’s emission is concentrated in the disk, featuring and average surface brightness of with peaks up to along the spiral arms. The analogous values for Althæa are and , respectively.

This can be explained as follows. FIR emission from the warm (), low density () component of the ISM is suppressed at  by the CMB (Gong et al. 2012; da Cunha et al. 2013; Pallottini et al. 2015; V15; App. B), as the upper levels of the [C ] transition cannot be efficiently populated through collisions and the spin temperature of the transition approaches the CMB one (see Pallottini et al., 2015, for possibility of [C ] detection from low density gas via CMB distortions). Thus, of the [C ] emission comes from dense (, cold (), mostly molecular disk gas. As noted in V15 (see in particular their Fig. 4) even when the CMB effect is neglected, the diffuse gas () account only for of the emission for galaxies with and , while it can be important in smaller objects (Olsen et al., 2015). The emissivity (in ) of such gas can be written as in P17 (in eq. 8, see also Vallini et al., 2013, 2017; Goicoechea et al., 2015):

(10)

for , i.e. the critical density for [C ] emission191919As the suppression of the CMB affects only the diffuse component (), no significant difference is expected in the emissivity from the disks of the two galaxies (eq. 10), that is composed of much higher () density material.. As the metallicity in the disk of the two galaxies is roughly similar (, see Tab. 1), difference in the luminosities is entirely explained by the larger density in Althæa. We stress once again that such density variation is a result of a more precise, non-equilibrium chemical network requiring to reach much higher densities before the gas is converted to stars. It is precisely that dense gas that accounts for a larger FIR line emissivity from PDRs.

We can also compare the calculated synthetic [C ] emission vs. SFR with observations (Fig. 12) obtained for dwarf galaxies (De Looze et al., 2014), and available  detections or upper-limits. The [C ] emission from Dahlia is lower than expected based on the local [C ]-SFR relation; its luminosity is also well below all upper limits for  galaxies. Although Althæa is times more luminous, even this object lies below the local relation, albeit only by . We believe that the reduced luminosity is caused by the combined effects of the CMB suppression and relatively lower . Note, however, that the predicted luminosity exceeds the upper limits derived for LAEs (e.g. Ouchi et al., 2013; Ota et al., 2014), but is broadly consistent with that of the handful of LBGs so far detected, like e.g. the four galaxies in the Pentericci et al. (2016). In general, observations are still rather sparse, with few [C ] detections with SFR comparable to Althæa (e.g. Capak et al., 2015). Also unclear is the amplitude of the scattering of the relation for  objects compared with local ones. Improvements in the understanding of the ISM structure are expected from deeper observations and/or other ions (e.g. [O ] Inoue et al., 2016; Carniani et al., 2017). Also helpful would be a larger catalogue of simulated galaxies (cfr. Ceverino et al., 2017), to control environmental effects.

4.2.2 Mid infrared emission

By inspecting the lower panel of Fig. 11 showing the predicted   line emission, we come to conclusions similar to those for the [C ]. Althæa outshines Dahlia by by delivering a total line luminosity of . Differently from the [C ] case, also the deviations from the mean are much more marked in Althæa, as appreciated from the Figure.

Note that the   line emissivity is enhanced in high density, high temperature regions. Indeed,  emission mostly arises from shocked-heated molecular gas, for which and for (see App. B).

For Dahlia, the disk density is relatively low, ; in addition only of the gas is warm enough to allow some () emission. In practice, such emission predominantly occurs along the outer spiral arms of the galaxy where these conditions are met due to the heating produced by SN explosions. In the denser Althæa disk, the gas emissivity can attain already at moderate . The brightness peaks are associated to a few () pockets of thousand-degree gas; they can be clearly identified in the map. This is particularly interesting because a galaxy like Althæa might be detectable at very  with SPICA, as suggested by Egami & et al. (2017, in preparation).

5 Conclusions

To improve our understanding of  galaxies we have studied the impact of  chemistry on their evolution, morphology and observed properties. To this end, we compare two zoom-in galaxy simulations implementing different chemical modelling. Both simulations start from the cosmological same initial conditions, and follow the evolution of a prototypical galaxy at resolved at the scale of giant molecular clouds (30 pc). Stars are formed according to a  dependent Schmidt-Kennicutt relation. We also account for winds from massive stars, SN explosions and radiation pressure in a stellar age/metallicity dependent fashion (see Sec. 2.1). The first galaxy is named Dahlia and  formation is computed from the Krumholz et al. (2009) equilibrium model; Althæa instead implements a non-equilibrium chemistry network, following Bovino et al. (2016). The key results can be summarized as follows:

  • The star formation rate of the two galaxies is similar, and increases with time reaching values close to at (see Fig. 2). However, Dahlia forms stars at a rate that is on average times larger than Althæa; it also shows a less prominent burst structure.

  • Both galaxies at have a SFR-stellar mass relation compatible with Jiang et al. (2016) observations (Fig. 3). Moreover, they both show a continuous time evolution from specific SFR of to . This is understood as an effect of the progressively increasing impact of stellar feedback hindering subsequent star formation events.

  • The non-equilibrium chemical model implemented in Althæa determines the atomic to molecular hydrogen transition to occur at densities exceeding 300 , i.e. about 10 times larges that predicted by equilibrium model used for Dahlia (Fig. 1). As a result, Althæa features a more clumpy and fragmented morphology (Fig. 6). This configuration makes SN feedback more effective, as noted in point (a) above (Fig. 4).

  • Because of the lower density and weaker feedback, Dahlia sits away from the Schmidt-Kennicutt relation; Althæa, instead nicely agrees with observations (Fig. 9). Note that although the SF efficiency is similar in the two galaxies and consistent with other simulations (Semenov et al., 2016), Dahlia is off the relation because of insufficient molecular gas content (Fig. 8).

  • We confirm that most of the emission from the C  and  is due to the dense gas forming the disk of the two galaxies. Because of Dahlia’s lower average density, Althæa outshines Dahlia by a factor of () in [C  ( ) line emission (Fig. 11). Yet, Althæa has a 10 times lower [C ] luminosity than expected from the locally observed [C ]-SFR relation (Fig. 12). Whether this relation does not apply at  or the line luminosity is reduced by CMB and metallicity effects remains as an open questions which can be investigated with future deeper observations.

To conclude, both Dahlia and Althæa follow the observed  SFR- relation. However, many other observed properties (Schmidt-Kennicutt relation, C  and  emission) are very different. This shows the importance of accurate, non-equilibrium implementation of chemical networks in early galaxy numerical studies.

Acknowledgments

We are grateful to the participants of The Cold Universe program held in 2016 at the KITP, UCSB, for discussions during the workshop. We thank P. Capelo, D. Celoria, E. Egami, D. Galli, T. Grassi, L. Mayer, S. Riolo, J. Wise for interesting and stimulating discussions. We thank the authors and the community of ramses and pymses for their work. AP acknowledges support from Centro Fermi via the project CORTES, “Cosmological Radiative Transfer in Early Structures”. AF acknowledges support from the ERC Advanced Grant INTERSTELLAR H2020/740120. RM acknowledge support from the ERC Advanced Grant 695671 “QUENCH” and from the Science and Technology Facilities Council (STFC). SS acknowledges support from the European Commission through a Marie Skłodowska-Curie Fellowship, program PRIMORDIAL, Grant No. 700907. This research was supported in part by the National Science Foundation under Grant No. NSF PHY11-25915.

References

  • Agertz & Kravtsov (2015) Agertz O., Kravtsov A. V., 2015, ApJ, 804, 18
  • Agertz et al. (2013) Agertz O., Kravtsov A. V., Leitner S. N., Gnedin N. Y., 2013, ApJ, 770, 25
  • Asano et al. (2013) Asano R. S., Takeuchi T. T., Hirashita H., Inoue A. K., 2013, Earth, Planets, and Space, 65, 213
  • Asplund et al. (2009) Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ARA&A, 47, 481
  • Bakes & Tielens (1994) Bakes E. L. O., Tielens A. G. G. M., 1994, ApJ, 427, 822
  • Barai et al. (2015) Barai P., Monaco P., Murante G., Ragagnin A., Viel M., 2015, MNRAS, 447, 266
  • Behroozi et al. (2013) Behroozi P. S., Wechsler R. H., Conroy C., 2013, ApJ, 770, 57
  • Bertelli et al. (1994) Bertelli G., Bressan A., Chiosi C., Fagotto F., Nasi E., 1994, A&A Supp., 106, 275
  • Black (1987) Black J. H., 1987, in Hollenbach D. J., Thronson Jr. H. A., eds, Astrophysics and Space Science Library Vol. 134, Interstellar Processes. pp 731–744
  • Black & Dalgarno (1976) Black J. H., Dalgarno A., 1976, ApJ, 203, 132
  • Bouché et al. (2007) Bouché N., et al., 2007, ApJ, 671, 303
  • Bouwens et al. (2015) Bouwens R. J., et al., 2015, ApJ, 803, 34
  • Bovino et al. (2016) Bovino S., Grassi T., Capelo P. R., Schleicher D. R. G., Banerjee R., 2016, A&A, 590, A15
  • Bradac et al. (2017) Bradac M., et al., 2017, ApJL, 836, L2
  • Bryan et al. (2014) Bryan G. L., et al., 2014, ApJS, 211, 19
  • Capak et al. (2015) Capak P. L., et al., 2015, Nature, 522, 455
  • Carilli & Walter (2013) Carilli C. L., Walter F., 2013, ARA&A, 51, 105
  • Carniani et al. (2017) Carniani S., et al., 2017, preprint, (arXiv:1701.03468)
  • Cazaux & Spaans (2009) Cazaux S., Spaans M., 2009, A&A, 496, 365
  • Ceverino et al. (2014) Ceverino D., Klypin A., Klimek E. S., Trujillo-Gomez S., Churchill C. W., Primack J., Dekel A., 2014, MNRAS, 442, 1545
  • Ceverino et al. (2017) Ceverino D., Glover S., Klessen R., 2017, preprint, (arXiv:1703.02913)
  • Choudhury et al. (2001) Choudhury T. R., Padmanabhan T., Srianand R., 2001, MNRAS, 322, 561
  • Ciardi & Ferrara (2001) Ciardi B., Ferrara A., 2001, MNRAS, 324, 648
  • Coles & Jones (1991) Coles P., Jones B., 1991, MNRAS, 248, 1
  • Cormier et al. (2015) Cormier D., et al., 2015, A&A, 578, A53
  • Daddi et al. (2010a) Daddi E., et al., 2010a, ApJ, 713, 686
  • Daddi et al. (2010b) Daddi E., et al., 2010b, ApJL, 714, L118
  • Davé et al. (2011) Davé R., Finlator K., Oppenheimer B. D., 2011, MNRAS, 416, 1354
  • De Looze et al. (2014) De Looze I., et al., 2014, A&A, 568, A62
  • Draine (1978) Draine B. T., 1978, ApJS, 36, 595
  • Dubois & Teyssier (2008) Dubois Y., Teyssier R., 2008, A&A, 477, 79
  • Dunlop (2013) Dunlop J. S., 2013, in Wiklind T., Mobasher B., Bromm V., eds, Astrophysics and Space Science Library Vol. 396, Astrophysics and Space Science Library. p. 223 (arXiv:1205.1543), doi:10.1007/978-3-642-32362-1_5
  • Egami & et al. (2017) Egami E., et al. 2017, in preparation, 0, 0
  • Federrath & Klessen (2013) Federrath C., Klessen R. S., 2013, ApJ, 763, 51
  • Ferland et al. (2013) Ferland G. J., et al., 2013, Revista Mexicana de Astronomia y Astrofisica, 49, 137
  • Fiacconi et al. (2017) Fiacconi D., Mayer L., Madau P., Lupi A., Dotti M., Haardt F., 2017, MNRAS, 467, 4080
  • Fujimoto et al. (2017) Fujimoto S., Ouchi M., Shibuya T., Nagai H., 2017, preprint, (arXiv:1703.02138)
  • Gallerani et al. (2016) Gallerani S., Pallottini A., Feruglio C., Ferrara A., Maiolino R., Vallini L., Riechers D. A., 2016, preprint, (arXiv:1604.05714)
  • Galli & Palla (1998) Galli D., Palla F., 1998, A&A, 335, 403
  • Gatto et al. (2015) Gatto A., et al., 2015, MNRAS, 449, 1057
  • Genzel et al. (2010) Genzel R., et al., 2010, MNRAS, 407, 2091
  • Glassgold et al. (2012) Glassgold A. E., Galli D., Padovani M., 2012, ApJ, 756, 157
  • Gleser et al. (2006) Gleser L., Nusser A., Ciardi B., Desjacques V., 2006, MNRAS, 370, 1329
  • Gnedin (2010) Gnedin N. Y., 2010, ApJL, 721, L79
  • Goicoechea et al. (2015) Goicoechea J. R., et al., 2015, ApJ, 812, 75
  • Gong et al. (2012) Gong Y., Cooray A., Silva M., Santos M. G., Bock J., Bradford C. M., Zemcov M., 2012, ApJ, 745, 49
  • González-López et al. (2014) González-López J., et al., 2014, ApJ, 784, 99
  • González et al. (2010) González V., Labbé I., Bouwens R. J., Illingworth G., Franx M., Kriek M., Brammer G. B., 2010, ApJ, 713, 115
  • Grassi et al. (2014) Grassi T., Bovino S., Schleicher D. R. G., Prieto J., Seifried D., Simoncini E., Gianturco F. A., 2014, MNRAS, 439, 2386
  • Grassi et al. (2017) Grassi T., Bovino S., Haugbølle T., Schleicher D. R. G., 2017, MNRAS, 466, 1259
  • Haardt & Madau (2012) Haardt F., Madau P., 2012, ApJ, 746, 125
  • Habing (1968) Habing H. J., 1968, Bull. Astron. Inst. Netherlands, 19, 421
  • Hahn & Abel (2011) Hahn O., Abel T., 2011, MNRAS, 415, 2101
  • Hartwig et al. (2015) Hartwig T., Clark P. C., Glover S. C. O., Klessen R. S., Sasaki M., 2015, ApJ, 799, 114
  • Heiderman et al. (2010) Heiderman A., Evans II N. J., Allen L. E., Huard T., Heyer M., 2010, ApJ, 723, 1019
  • Hirashita & Ferrara (2002) Hirashita H., Ferrara A., 2002, MNRAS, 337, 921
  • Hollenbach & McKee (1979) Hollenbach D., McKee C. F., 1979, ApJS, 41, 555
  • Hopkins et al. (2013) Hopkins P. F., Narayanan D., Murray N., 2013, MNRAS, 432, 2647
  • Hopkins et al. (2017) Hopkins P. F., et al., 2017, preprint, (arXiv:1702.06148)
  • Inoue et al. (2016) Inoue A. K., et al., 2016, Science, 352, 1559
  • Jiang et al. (2016) Jiang L., et al., 2016, ApJ, 816, 16
  • Jura (1975) Jura M., 1975, ApJ, 197, 575
  • Kanekar et al. (2013) Kanekar N., Wagg J., Ram Chary R., Carilli C. L., 2013, ApJL, 771, L20
  • Katz et al. (2016) Katz H., Kimm T., Sijacki D., Haehnelt M., 2016, preprint, (arXiv:1612.01786)
  • Kennicutt (1998) Kennicutt Jr. R. C., 1998, ApJ, 498, 541
  • Kennicutt & Evans (2012) Kennicutt R. C., Evans N. J., 2012, ARA&A, 50, 531
  • Kim et al. (2014) Kim J.-h., et al., 2014, ApJS, 210, 14
  • Kim et al. (2016) Kim J.-h., et al., 2016, ApJ, 833, 202
  • Klessen & Glover (2014) Klessen R. S., Glover S. C. O., 2014, preprint, (arXiv:1412.5182)
  • Knudsen et al. (2016) Knudsen K. K., Richard J., Kneib J.-P., Jauzac M., Clément B., Drouart G., Egami E., Lindroos L., 2016, MNRAS, 462, L6
  • Kroupa (2001) Kroupa P., 2001, MNRAS, 322, 231
  • Krumholz (2015) Krumholz M. R., 2015, preprint, (arXiv:1511.03457)
  • Krumholz et al. (2008) Krumholz M. R., McKee C. F., Tumlinson J., 2008, ApJ, 689, 865
  • Krumholz et al. (2009) Krumholz M. R., McKee C. F., Tumlinson J., 2009, ApJ, 693, 216
  • Krumholz et al. (2012) Krumholz M. R., Dekel A., McKee C. F., 2012, ApJ, 745, 69
  • Labadens et al. (2012) Labadens M., Chapon D., Pomaréde D., Teyssier R., 2012, in Ballester P., Egret D., Lorente N. P. F., eds, Astronomical Society of the Pacific Conference Series Vol. 461, Astronomical Data Analysis Software and Systems XXI. p. 837
  • Lada et al. (2010) Lada C. J., Lombardi M., Alves J. F., 2010, ApJ, 724, 687
  • Laporte et al. (2017) Laporte N., et al., 2017, ApJL, 837, L21
  • Leitherer et al. (1999) Leitherer C., et al., 1999, ApJS, 123, 3
  • Mac Low (1999) Mac Low M.-M., 1999, ApJ, 524, 169
  • Madau & Dickinson (2014) Madau P., Dickinson M., 2014, ARA&A, 52, 415
  • Maio & Tescari (2015) Maio U., Tescari E., 2015, MNRAS, 453, 3798
  • Maio et al. (2016) Maio U., Petkova M., De Lucia G., Borgani S., 2016, MNRAS, 460, 3733
  • Maiolino et al. (2015) Maiolino R., et al., 2015, MNRAS, 452, 54
  • Martizzi et al. (2015) Martizzi D., Faucher-Giguère C.-A., Quataert E., 2015, MNRAS, 450, 504
  • McKee & Krumholz (2010) McKee C. F., Krumholz M. R., 2010, ApJ, 709, 308
  • Mecke et al. (1994) Mecke K. R., Buchert T., Wagner H., 1994, A&A, 288, 697
  • Murray (2011) Murray N., 2011, ApJ, 729, 133
  • Novikov et al. (2000) Novikov D., Schmalzing J., Mukhanov V. F., 2000, A&A, 364, 17
  • O’Shea et al. (2015) O’Shea B. W., Wise J. H., Xu H., Norman M. L., 2015, ApJL, 807, L12
  • Olsen et al. (2015) Olsen K. P., Greve T. R., Narayanan D., Thompson R., Toft S., Brinch C., 2015, ApJ, 814, 76
  • Osterbrock (1989) Osterbrock D. E., 1989, Astrophysics of gaseous nebulae and active galactic nuclei. Mill Valley, CA, University Science Books
  • Ostriker & McKee (1988) Ostriker J. P., McKee C. F., 1988, Reviews of Modern Physics, 60, 1
  • Ota et al. (2014) Ota K., et al., 2014, ApJ, 792, 34
  • Ouchi et al. (2013) Ouchi M., et al., 2013, ApJ, 778, 102
  • Padoan et al. (2012) Padoan P., Haugbølle T., Nordlund A., 2012, ApJL, 759, L27
  • Pallottini et al. (2014a) Pallottini A., Ferrara A., Gallerani S., Salvadori S., D’Odorico V., 2014a, MNRAS, 440, 2498
  • Pallottini et al. (2014b) Pallottini A., Gallerani S., Ferrara A., 2014b, MNRAS, 444, L105
  • Pallottini et al. (2015) Pallottini A., Gallerani S., Ferrara A., Yue B., Vallini L., Maiolino R., Feruglio C., 2015, MNRAS, 453, 1898
  • Pallottini et al. (2017) Pallottini A., Ferrara A., Gallerani S., Vallini L., Maiolino R., Salvadori S., 2017, MNRAS, 465, 2540
  • Pentericci et al. (2016) Pentericci L., et al., 2016, ApJL, 829, L11
  • Pérez-Montero (2017) Pérez-Montero E., 2017, Publ. Astr. Soc. Pac., 129, 043001
  • Petkova & Maio (2012) Petkova M., Maio U., 2012, MNRAS, 422, 3067
  • Planck Collaboration et al. (2014) Planck Collaboration et al., 2014, A&A, 571, A16
  • Rahmati et al. (2013) Rahmati A., Pawlik A. H., Raicevic M., Schaye J., 2013, MNRAS, 430, 2427
  • Rasera & Teyssier (2006) Rasera Y., Teyssier R., 2006, A&A, 445, 1
  • Rey-Raposo et al. (2017) Rey-Raposo R., Dobbs C., Agertz O., Alig C., 2017, MNRAS, 464, 3536
  • Richings et al. (2014) Richings A. J., Schaye J., Oppenheimer B. D., 2014, MNRAS, 442, 2780
  • Rosdahl et al. (2015) Rosdahl J., Schaye J., Teyssier R., Agertz O., 2015, MNRAS, 451, 34
  • Rosdahl et al. (2017) Rosdahl J., Schaye J., Dubois Y., Kimm T., Teyssier R., 2017, MNRAS, 466, 11
  • Roskar et al. (2014) Roskar R., Teyssier R., Agertz O., Wetzstein M., Moore B., 2014, MNRAS, 444, 2837
  • Roussel et al. (2007) Roussel H., et al., 2007, ApJ, 669, 959
  • Scannapieco et al. (2012) Scannapieco C., et al., 2012, MNRAS, 423, 1726
  • Schaerer et al. (2015) Schaerer D., Boone F., Zamojski M., Staguhn J., Dessauges-Zavadsky M., Finkelstein S., Combes F., 2015, A&A, 574, A19
  • Schmalzing & Buchert (1997) Schmalzing J., Buchert T., 1997, ApJL, 482, L1
  • Schmalzing & Gorski (1998) Schmalzing J., Gorski K. M., 1998, MNRAS, 297, 355
  • Schmidt (1959) Schmidt M., 1959, ApJ, 129, 243
  • Semenov et al. (2016) Semenov V. A., Kravtsov A. V., Gnedin N. Y., 2016, ApJ, 826, 200
  • Shibuya et al. (2015) Shibuya T., Ouchi M., Harikane Y., 2015, ApJS, 219, 15
  • Smith et al. (2017) Smith B. D., et al., 2017, MNRAS, 466, 2217
  • Spinoglio & et al. (2017) Spinoglio L., et al. 2017, in preparation, 0, 0
  • Stanway (2017) Stanway E. R., 2017, preprint, (arXiv:1702.07303)
  • Stasińska (2007) Stasińska G., 2007, preprint, (arXiv:0704.0348)
  • Stinson et al. (2006) Stinson G., Seth A., Katz N., Wadsley J., Governato F., Quinn T., 2006, MNRAS, 373, 1074
  • Tacconi et al. (2010) Tacconi L. J., et al., 2010, Nature, 463, 781
  • Teyssier (2002) Teyssier R., 2002, A&A, 385, 337
  • Teyssier et al. (2013) Teyssier R., Pontzen A., Dubois Y., Read J. I., 2013, MNRAS, 429, 3068
  • Timmermann et al. (1996) Timmermann R., Bertoldi F., Wright C. M., Drapatz S., Draine B. T., Haser L., Sternberg A., 1996, A&A, 315, L281
  • Tomassetti et al. (2015) Tomassetti M., Porciani C., Romano-Díaz E., Ludlow A. D., 2015, MNRAS, 446, 3330
  • Torrey et al. (2017) Torrey P., Hopkins P. F., Faucher-Giguère C.-A., Vogelsberger M., Quataert E., Keres D., Murray N., 2017, MNRAS, 467, 2301
  • Turner et al. (1977) Turner J., Kirby-Docken K., Dalgarno A., 1977, ApJS, 35, 281
  • Valle et al. (2002) Valle G., Ferrini F., Galli D., Shore S. N., 2002, ApJ, 566, 252
  • Vallini et al. (2013) Vallini L., Gallerani S., Ferrara A., Baek S., 2013, MNRAS, 433, 1567
  • Vallini et al. (2015) Vallini L., Gallerani S., Ferrara A., Pallottini A., Yue B., 2015, ApJ, 813, 36
  • Vallini et al. (2017) Vallini L., Ferrara A., Pallottini A., Gallerani S., 2017, MNRAS,
  • Verner & Ferland (1996) Verner D. A., Ferland G. J., 1996, ApJS, 103, 467
  • Wakelam et al. (2012) Wakelam V., et al., 2012, ApJS, 199, 21
  • Watson et al. (2015) Watson D., Christensen L., Knudsen K. K., Richard J., Gallazzi A., Michałowski M. J., 2015, Nature, 519, 327
  • Weaver et al. (1977) Weaver R., McCray R., Castor J., Shapiro P., Moore R., 1977, ApJ, 218, 377
  • Webber (1998) Webber W. R., 1998, ApJ, 506, 329
  • Weingartner & Draine (2001) Weingartner J. C., Draine B. T., 2001, ApJ, 563, 842
  • Willott et al. (2015) Willott C. J., Carilli C. L., Wagg J., Wang R., 2015, ApJ, 807, 180
  • Wise et al. (2012) Wise J. H., Abel T., Turk M. J., Norman M. L., Smith B. D., 2012, MNRAS, 427, 311
  • Wolcott-Green et al. (2011) Wolcott-Green J., Haiman Z., Bryan G. L., 2011, MNRAS, 418, 838
  • Wolfire et al. (2003) Wolfire M. G., McKee C. F., Hollenbach D., Tielens A. G. G. M., 2003, ApJ, 587, 278
  • Yoshiura et al. (2017) Yoshiura S., Shimabukuro H., Takahashi K., Matsubara T., 2017, MNRAS, 465, 394
  • da Cunha et al. (2013) da Cunha E., et al., 2013, ApJ, 766, 13

Appendix A Minkowsky functionals

In general, Minkowsky functionals are mathematical tools that give a complete characterization of the morphology of a field. In astrophysics they have been proposed as a mean to give a description of the large scale structure (e.g. Mecke et al., 1994), to study the topology of H   bubbles for reionization studies (e.g. Gleser et al., 2006; Yoshiura et al., 2017), and (for fields) analyze CMB anisotropies and non-gaussianity (e.g. Schmalzing & Gorski, 1998; Novikov et al., 2000).

A formal definition can be given following Schmalzing & Buchert (1997). Let denote a scalar field defined on a subset of with volume . Let us take such that it has zero mean () and variance . Then we can define the excursion set as the ensemble of regions in satisfying . Then, the Minkowsky functionals can be defined in terms of volume () and surface () integrals as a function of the threshold :

(11a)
(11b)
(11c)
(11d)

where is the Heaviside function, the surface beading the excursion set , and and the two principal curvatures of the surface. In practical terms, is a measure of the volume filling factor of the excursion set with threshold , of the surface, of the mean curvature (sphericity/concavity) and of the Euler characteristic (shape of components).

The curvatures on the surface can be expressed via the Koenderink invariant (Gleser et al., 2006, see appendix A and reference therein): by adopting the Einstein sum convention we can write

(12a)
(12b)
(12c)

where is the Levi-Civita symbol, is the Kronecker delta, and is the -th component of the partial derivative operator. Finally, eq.s 11 can be suitably expressed as integral over the volume by using the following relation

(13)

where is the Dirac delta.

As an illustrative example, we can compute the Minkowsky functionals for a zero-mean Gaussian random field with variance and variance of the tangent field . For such Gaussian field, the Minkowsky functionals can be expressed using the following analytical expression (Schmalzing & Buchert, 1997, see also Gleser et al. 2006)

(14a)
(14b)
(14c)
(14d)

where , , and are numerical constant with values and .

We numerically compute the Minkowsky functionals for th e field with on a unigrid box with volume , that thus result in a tangent field variance of . The resulting Minkowsky functionals are plotted in Fig. 13. We find a very good match with the analytical values.

Since the chosen Gaussian field is an approximation to the quasi-linear regime of the cosmic density field (e.g. Coles & Jones, 1991; Choudhury et al., 2001), it is intuitive to analyze the properties of its Minkowsky functionals. In Fig. 13 the filling factor gives the probability of finding regions with increasing overdensity ; note that at (mean density), , i.e. it is equiprobable to find voids () and overdense regions (). Both voids and overdense regions are isolated