Free Access
Issue
A&A
Volume 575, March 2015
Article Number A33
Number of page(s) 25
Section Numerical methods and codes
DOI https://doi.org/10.1051/0004-6361/201424462
Published online 19 February 2015

© ESO, 2015

1. Introduction

After more than twenty years, (hybrid)-semi analytical models (SAMs) are the best tool for the physical interpretation of large surveys. Under the hypothesis that baryonic processes cannot strongly influence the structuration of the dark-matter, the dark-matter evolution is decoupled from the computation of the baryonic component. Even if this decoupling is a strong assumption (see van Daalen et al. 2011), it allows a wide range of physical prescriptions to be explored in a realistic cosmological context and in a short computational time.

Originally proposed by White & Frenk (1991), SAMs are still being developed (e.g. Cole 1991; Cole et al. 2000; Hatton et al. 2003; Baugh 2006; Croton et al. 2006; Cattaneo et al. 2006; Somerville et al. 2008; Guo et al. 2011; Henriques et al. 2013). After a first version described in a dedicated paper, the updates of these models are often fragmented into a large number of publications. In this context, we found it necessary to clarify all the steps of the modelling and build an up-to-date model based on the latest results obtained by hydrodynamic simulations and the comparison of empirical models with observations. Inspired by the previous GalICS model (Hatton et al. 2003), we revisit all the standard mechanisms acting on the formation and the evolution of galaxies in our new SAM, step by step, and in a single paper.

We mainly focus on the gas cycle, which is the main actor in galaxy stellar mass assembly. It is thus described from the cosmological smooth accretion to the complex feedback mechanisms. Obviously, all physical processes act in the gas cycle with different time scales. Indeed, gas-accretion, star formation, and clump migration time scales are all different. For this reason, we propose an adaptive time-step scheme to solve the gas evolution equations.

Like its predecessor (Hatton et al. 2003), the new model is based on a dark-matter N-body simulation computed in the Λ cold dark-matter paradigm (Λ-CDM). Indeed, even if the growth of dark-matter haloes can be followed using the Press & Schechter (1974) formalism and their extensions (e.g. Cole 1991; Kauffmann et al. 1993; Cole et al. 1994; Somerville & Kolatt 1999), the evolution of the dark-matter component is now commonly extracted from huge dark-matter simulations (e.g. Kauffmann et al. 1999; Helly et al. 2003; Hatton et al. 2003). In such high-resolution simulations, the dark-matter halo population covers a wide range of mass (Mh ∈ [ 107:1013 ] M) and scale, from dwarf haloes to groups. To follow its time evolution, the hierarchical growth of dark-matter structures is commonly represented by merger trees in which each branch represents the growth of a given halo. The analysis of dark-matter simulations have set the current paradigm in which the dark-matter halo growth follows two modes (e.g. Lacey & Cole 1994; Tormen et al. 1997; Wechsler et al. 2002; van den Bosch 2002; Zhao et al. 2009; Fakhouri & Ma 2010; Fakhouri et al. 2010; Genel et al. 2010):

  • A smooth accretion that is difficult to define with accuracy. Indeed, it is closely linked to the mass resolution. For example, a 106 M dark-matter particle can be seen (or defined) as an unresolved dark-matter halo. However Genel et al. (2010) show that the smooth accretion process converges to non-zero values when resolution increases. More precisely, they demonstrate that the smooth accretion can account for up to 40% of the total mass of a dark-matter halo. These results prove that the naive expectation that all accretion would come through mergers when the resolution increases is wrong. The smooth dark-matter accretion is therefore a major vector of the growth of dark-matter structures. It is subject of a detailed description in this paper.

  • The second mode is obviously linked to the merging of pre-existing structures.

Applied to the evolution of the dark-matter component, the baryonic post-treatment is articulated around three new prescriptions:

  • A two-phase smooth baryonic accretion, with a cold and a hot mode that are built on the smooth dark-matter accretion.

  • A complete monitoring of the hot phase. The evolution of the mean temperature and density profile parameters are followed using explicit energy conservation laws, under the hydrostatic equilibrium hypothesis.

  • A new approach for describing of the disc gravitational instabilities. Based on the formation, migration and disruption of giant-clumps, this new prescription adapted for SAMs relies on the work by Dekel et al. (2009b).

Beyond these three new prescriptions, the model examines the standard mechanisms of star formation regulation. For the low-mass structures two processes are in play:

The galaxy/dark-matter potential well is too deep for high-mass structures (Mh> 1012 M). In this case, SNe and photoionisation processes are not powerful enough to eject some gas or to prevent the infall of baryons. SNe still produce winds but ejecta velocities are not sufficient to expel the gas. It is thus rapidly re-accreted in the disc. A popular interpretation that can explain the low star formation efficiency in Mh> 1012 M haloes is based on the activity of a super massive black hole (SMBH) hosted by the galaxy (e.g. Croton et al. 2006; Cattaneo et al. 2006; Bower et al. 2006; Malbon et al. 2007; Somerville et al. 2008; Ostriker et al. 2010; Bower et al. 2012). If a fraction of the energy produced by the accretion on the SMBH is used to heat the hot gas phase, the efficiency of the cooling process is decreased. The accretion of the gas in the galaxy, and thus the star formation, are then reduced.

With these standard recipes some profound disagreements between model predictions and observations are still present, even if observed local (z = 0) galaxy properties are generally reproduced well. For example, one of the most important problems is the strong excess of low-mass galaxies predicted at high z. Indeed the stellar mass functions measured at 1 <z< 4, for example by Ilbert et al. (2010, 2013), or Caputi et al. (2011), demonstrate that SAMs predict a number of low-mass galaxies roughly ten times larger than what observed. Models cannot also easily explain the cosmic star formation history in detail and the bulk of star formation activity observed at 1 <z< 3. In addition, recent studies show that the link between dark-matter haloes and their host galaxies is not well understood (e.g. Leauthaud et al. 2012). Some empirical models (e.g. Béthermin et al. 2012; Behroozi et al. 2013a,b) suggest a non-monotonically growth of the stellar mass with the dark-matter halo mass. Indeed, the star formation activity derived from observations seems to be very inefficient both in low-mass (Mh< 1010 M) and high-mass haloes (Mh> 1013 M). The star formation activity seems to peak for galaxies evolving at intermediate halo mass (Mh ≃ 1012 M) (Eke et al. 2004, 2005; Guo et al. 2010).

While the current paper describes all the standard prescriptions used in the new model in detail, the reasons for the tensions between models and observations are discussed in a companion paper where we compare the predictions with fundamental observations, such as the stellar mass function, Mh versus M relation, IR luminosity function, and specific star formation rate (Cousin et al. 2015). The companion paper describes, in particular, the impact of SNe feedback and photoionisation process on the stellar mass assembly. It emphasises the difficulty matching observations by using only the standard star formation regulation processes and proposes, in this context, an ad hoc modification of the standard star formation cycle.

The paper is organised as followed. In Sect. 2 we describe the dark-matter simulation used to process the dark-matter evolution. We then present an explicit description of the smooth dark-matter process based on background particles. In Sect. 3, we describe the bimodal baryonic accretion prescription adopted in the model. It is based on a photoionisation model and on two gas reservoirs dedicated to cold and hot phases. Section 4 presents a complete modelling of the monitoring of the hot gas phase. Ejecta processes coming from SNe and/or SMBH activity are also presented, together with the prescription for computating the gas escape fraction and the cooling rate. In Sect. 5, we present all the properties and evolution laws of discs and bulges. The new clumpy component that derives from disc instabilities is also fully described in this section. In Sect. 8 we describe the effect of merger events onto galaxy evolution. Section 9 closes the paper by giving an overview and some research perspectives. In addition, we have added in Appendix A.4, a description of the adaptive time-step scheme used to follow the different dynamical time scales of the gas cycle.

Table 1

Dark-matter symbols and their definition.

2. dark-matter

In this section, we first describe the dark-matter merger-tree structure and then define the smooth accretion process. We also give the useful dark-matter halo properties. All symbols and their definition are listed in Table 1.

2.1. Dark-matter evolution

Like its predecessor, the new version eGalICS uses a pure N-body simulation to follow the dark-matter evolution. We use an N-body cosmological simulation based on the WMAP-3yr cosmology (Ωm = 0.24, ΩΛ = 0.76, fb = 0.16, h = 0.73). In a volume of (100h-1)3 ≃ 150 Mpc3, 10243 particles evolve with an elementary mass of mp = 8.536 × 107 M. Dark-matter merger trees (see Fig. 1) are extracted with the halo-finder program AdaptaHOP proposed by Aubert et al. (2004) and extended by Tweed et al. (2009). The dark-matter assembly is therefore pre-computed for the different time steps. Its evolution cannot be modified by baryonic processes.

The halo-finder version used here (Tweed et al. 2009) allows us to extract and separate different levels of structures: main haloes (M-H) and subhaloes (S-H or satellites). In the two cases, we only consider dark-matter structures containing at least 20 dark-matter particles. This limit gives a minimum dark-matter halo mass Mhmin=1.707×109M\hbox{$M_{\rm h}^{\rm min} = 1.707\times10^9~\Msun$}. Subhaloes are extracted and identified as overdensities in a pre-existing structure that can be a M-H or another S-H. Currently, it is possible to extract up to six different levels of S-H. The halo finder generates all the links between all haloes in time and in structuration level.

thumbnail Fig. 1

Merger tree showing the evolution of the dark-matter smooth accretion. We only show the fifty most massive branches. The junctions with less massive branches are indicated with small arrowheads. The main branch (on the left) is built up through a smooth accretion process acc (Eq. (1)) and multiple mergers with other haloes. The evolution of the smooth dark-matter accretion through the merger tree is shown with colour. A halo without any dark-matter accretion is shown with a black symbol. Circles are used for main haloes (M-H), while squares are used for subhaloes (S-H). The size of the symbols is proportional to the dark-matter halo virial mass. It is clearly visible that before a merger event with a larger structure, a halo is systematically identified as a substructure of this more massive halo. During this subhalo phase, the dark-matter accretion rate strongly decreases and may even stop.

2.2. Dark-matter background accretion

To follow the dark-matter halo growth and, later on, the baryonic accretion process, we identify at a given time step n, all particles that are detected in a halo , and that have never been identified in an other halo ’ before. These background particles (mp\hbox{$m_{\rm p}^{\star}$}) play an important role. Indeed, we consider them as the only source of new baryonic mass. In this framework, the accretion rate on a given structure at a given time step is defined as acc=mptntn1·\begin{eqnarray} \dot{M}_{\rm acc} = \dfrac{\sum_{\Halo} m_{\rm p}^{\star}}{t^{n}-t^{n-1}}\cdot \label{dm_acc_rate} \end{eqnarray}(1)Figure 1 shows the evolution of the dark-matter accretion rate acc. We can see that before merging with a larger structure (M-H), a halo is systematically identified as a S-H of this more massive halo. During this pre-merger period, it is clear that the dark-matter accretion rate (acc) decreases strongly and may even stop. In these conditions, the diffuse accretion is then supported by the most massive structure. The decrease in background accretion is also visible during fly-by events, as identified in the 18th branch of the merger tree shown in Fig. 1.

By following with time the background accretion process for a given structure i, we can define an integrated accretion mass: Macc,in=Macc,in1􏽼0􏽻􏽺0􏽽I+mp􏽼0􏽻􏽺0􏽽II+jiMacc,jn1􏽼0􏽻􏽺0􏽽III.\begin{eqnarray} M_{{\rm acc},i}^{n} = \underbrace{{M}_{{\rm acc},i}^{n-1}}_{\rm I} + \underbrace{\sum_{\Halo} m_{\rm p}^{\star}}_{\rm II} + \underbrace{\sum_{j \neq i} M_{{\rm acc},j}^{n-1}}_{\rm III}. \label{dm_acc_mass} \end{eqnarray}(2)In this equation:

  • I is the integrated background mass accreted by the halo i up to the previous time step n − 1.

  • II is the background mass accreted during the time step [n − 1, n]. This mass is computed using all particles identified in the halo i at the time step n and that have never been identified in an other halo before.

  • III corresponds to the integrated background mass accreted by all the progenitors j of the halo i.

2.3. Remark on the dark-matter background accretion

For M-H, the integrated accretion mass (Macc) is closed to the well known halo finder mass Mfof, which is defined as the mass of all particles instantaneously detected in a structure. In contrast, the integrated background mass Macc for S-H is often higher than Mfof. Indeed, S-H are sensitive to the stripping process, and a fraction of the halo mass (Mfof) can be ripped by gravitational interactions acting before a merger event. In the context of our dark-matter background accretion process, one of the key points of the model resides in the fact that the stripped particles that have already been identified in a structure cannot be considered in any other structure.

The mass estimator Macc summarises the accretion history of a given halo. When the instantaneous mass is needed, we use the virial mass Mvir, as in other SAMs. This mass is computed using all particles that follow the virial theorem for a given halo.

2.4. Dark-matter halo properties

2.4.1. Density profile

The spatial distribution of a dark-matter halo mass is described by a Navarro Frenk and White (NFW) density profile (Navarro et al. 1995, 1996, 1997) defined as: ρ(r)=ρdmrrdm(1+rrdm)2\begin{eqnarray} \rho(r) = \dfrac{\rho_{\rm dm}}{\frac{r}{r_{\rm dm}}\left(1 + \frac{r}{r_{\rm dm}}\right)^2} \label{rho_dm} \end{eqnarray}(3)where rdm is the characteristic scale radius of the dark-matter halo and ρdm/ 4 is the density at r = rdm. With this density profile, the dark-matter halo mass Mh( <r) contained in the radius r is given by (e.g. Suto et al. 1998; Makino et al. 1998; Komatsu & Seljak 2001) Mh(<r)=4πρdmrdm3φ(x)\begin{eqnarray} M_{\rm h}({<}r) = 4\pi\rho_{\rm dm}r_{\rm dm}^3\phi\left(x\right) \label{dm_mass_r} \end{eqnarray}(4)where φ(x) is a geometrical function defined as φ(x)=ln(1+x)x1+xwithx=rrdm·\begin{eqnarray} \phi(x) = \ln(1 + x) - \dfrac{x}{1 + x}~~\mbox{with}~~x = \frac{r}{r_{\rm dm}}\cdot \label{phix} \end{eqnarray}(5)Obviously, if rvir is the virial radius of the dark-matter halo extracted from the dark-matter simulation (see Sect. 2.3 of Hatton et al. 2003 for more information): Mh(<rvir)=Mvir=4πρdmrdm3φ(cdm)\begin{eqnarray} M_{\rm h}({<}r_{\rm vir}) = M_{\rm vir} = 4\pi\rho_{\rm dm}r_{\rm dm}^3\phi\left(c_{\rm dm}\right) \label{dm_mass_rvir} \end{eqnarray}(6)where cdm is the dark-matter concentration parameter defined as cdm = rvir/rdm. With the dark-matter virial mass and the virial radius, it is possible to define the dynamical time of the dark-matter structure: tdyn=rvir3GMvir·\begin{eqnarray} t_{\rm dyn} = \sqrt{\dfrac{r_{\rm vir}^3}{GM_{\rm vir}}}\cdot \label{t_dyn} \end{eqnarray}(7)This time will be used in the next sections to define the baryonic accretion rate.

2.4.2. Virial temperature and escape velocity

In the hot atmosphere evolution monitoring, we use the virial temperature Tvir linked to the dark-matter halo: Tvir=35.9GMvirrvir[K].\begin{eqnarray} T_{\rm vir} = 35.9\dfrac{GM_{\rm vir}}{r_{\rm vir}}~~[{\rm K}]. \label{T_vir} \end{eqnarray}(8)We also use the escape velocity of the dark-matter halo defined as Vesc,dm=8πGρdmrdm2ln(1+x)xwithx=rrdm\begin{eqnarray} V_{\rm esc,dm} = \sqrt{8\pi G \rho_{\rm dm} r_{\rm dm}^2\dfrac{\ln\left(1 + x\right)}{x}}~~\mbox{with}~~x = \frac{r}{r_{\rm dm}} \label{Vesc_dm} \end{eqnarray}(9)and the circular velocity of the dark-matter halo defined as Vcirc,dm=4πGρdmrdm2φ(x)xwithx=rrdm·\begin{eqnarray} V_{\rm circ,dm} = \sqrt{4\pi G \rho_{\rm dm} r_{\rm dm}^2\dfrac{\phi(x)}{x}}~~\mbox{with}~~x = \frac{r}{r_{\rm dm}}\cdot \label{Vcirc_dm} \end{eqnarray}(10)

3. Baryonic accretion

The evolution of the dark-matter component is therefore known using a large set of merger trees. The dark-matter halo properties are listed and defined in the previous section. We now have to add the baryonic content on the hierarchical structure of the dark-matter. This section describes the link between the dark-matter smooth accretion and the baryons that will feed the galaxy. All parameters used in this section and their definition are listed in Table 2.

As explained previously, we consider that the dark-matter background-accreted mass (acc) is the only vector of baryonic accretion. The link between the dark-matter and baryonic accretion is defined through the effective accreted baryonic fraction, fb (Eq. (11)).

Table 2

Baryonic accretion symbols and their definition.

3.1. Photoionisation

The value taken by this parameter fb follows two distinct regimes. Before reionisation (zreion), we consider that the baryons feed the dark-matter structures following the universal baryonic fraction fb=ΩbΩm\hbox{$\left<f_{\rm b}\right>=\frac{\Omega_{\rm b}}{\Omega_{\rm m}}$}. After reionisation, the heating of the gas generated by the first generation of stars limits the baryonic gas accretion in the smallest structures (Kauffmann et al. 1993). Originally proposed by Doroshkevich et al. (1967), the photoionisation mechanism has been developed in the CDM paradigm by Couchman & Rees (1986), Ikeuchi (1986), and Rees (1986). This process is widely used in the context of galaxy evolution (Efstathiou 1992; Babul & Rees 1992; Shapiro et al. 1994; Quinn et al. 1996; Thoul & Weinberg 1996; Bullock et al. 2000; Gnedin 2000; Benson et al. 2002; Somerville 2002; Somerville et al. 2008, 2012; Croton et al. 2006; Hoeft et al. 2006; Okamoto et al. 2008; Guo et al. 2011; Henriques et al. 2013).

The effective baryonic fraction fb therefore depends on both the redshift z and dark-matter halo mass Mvir. The most commonly used formulation is that proposed by Gnedin (2000) and Kravtsov et al. (2004): fb=fb{\begin{eqnarray} \begin{small} f_{\rm b} = \left<f_{\rm b}\right>\left\{ \begin{array}{ll} \left[1 + (2^{\alpha/3}-1)\left(\dfrac{M_{\rm vir}}{M_{\rm c}(z)}\right)^{-\alpha}\right]^{-3/\alpha} & \! \!: {\rm if}\quad z<z_{\rm reion} \\ & \\ 1 & \!\!\mbox{: otherwise}. \end{array}\right. \end{small} \label{baryonic_fraction} \end{eqnarray}(11)Following Okamoto et al. (2008), we use α = 2. The characteristic halo mass follows1: Mc(z)=8.22×109exp(0.7z)[M].\begin{eqnarray} M_{\rm c}(z) = 8.22\times10^{9}\exp(-0.7z)~~[\Msun ]. \label{filtering_mass} \end{eqnarray}(12)According to this filtering mass, haloes with Mvir = Mc really collect only half of the baryons coupled with the accreted dark-matter.

Figure 2 shows the evolution of the effective accreted baryonic fraction as a function of both the redshift and the dark-matter halo mass. The mass resolution limit is 20 dm-particules. At this mass scale, the effective baryonic fraction is strongly reduced (fb< ⟨fb⟩ / 2) only in the smallest structures evolving at z< 2. As demonstrated in Cousin et al. (2015), the impact of photoionisation is limited. Even if the amount of accreted gas is reduced, it is not sufficient to significantly limit the star formation activity.

thumbnail Fig. 2

Evolution of the effective baryonic fraction (fb: Eq. (11)) with dark-matter halo mass and redshift. We use the Okamoto et al. (2008) prescription with a reionisation redshift zreion = 9 (see text for more details). The vertical grey line indicates our mass resolution limit (20 dm-particules). At this mass resolution and for a redshift of about 2, the effective baryonic fraction fb applied to the accretion flux is ΩbΩm\hbox{$\frac{\Omega_{\rm b}}{\Omega_{\rm m}}$}.

3.2. Bimodal accretion

Following the definitions summarised by Eqs. (1) and (11), the baryonic accretion rate on a dark-matter halo is given by bar=fb(z,Mvir)acc.\begin{eqnarray} \dot{M}_{\rm bar} = f_{\rm b}(z,M_{\rm vir})\dot{M}_{\rm acc}. \label{bar_acc_rate} \end{eqnarray}(13)As proposed by Khochfar & Silk (2009) we adopt an explicit bimodal accretion composed of a cold and a hot gas flow. These two parallel accretion modes have been discussed in Kereš et al. (2005); Dekel et al. (2009a,b); Khochfar & Silk (2009); van de Voort et al. (2011), and Faucher-Giguère et al. (2011).

thumbnail Fig. 3

Large-scale exchanges. On a large scale, dark-matter and baryons are coupled. During the accretion process and while the dark-matter smooth accretion (acc) feeds the dark-matter halo (Mh), the associated baryons (bar) feed a baryonic reservoir (Mbar). In the context of bimodal accretion, this baryonic reservoir is divided into two parts: the hot reservoir representing the hot atmosphere (sh and Mhot), and the cold reservoir representing filamentary streams (flt and Mcold). The cold gas feeds the galaxy (ff) directly. The hot gas follows radiative cooling process and then also feeds the galaxy (cool). As described in the text, a galaxy can eject material. The ejecta (wind) are continuously added to the hot-gas reservoir (Mhot) or definitively lost in the IGM according to their velocity distribution. Indeed, the gas stored in the IGM reservoir will never be considered.

On the first hand, we distinguish a cold mode, for which the gas is accreted through filamentary streams. This mode is very efficient. Even if the gas is previously shock-heated (Nelson et al. 2013), it cools with a very high efficiency (tcooltdyn) and falls in the centre of the halo in an approximately free-fall time (tdyn, see Eq. (7)). The first galaxy discs are formed quickly with this accretion mode. This mode dominates the growth of galaxies at high redshifts (z> 3) and the growth of galaxies in the lowest mass haloes (Mvir< 1011 M) at any times.

On the other hand, in more massive haloes (Mvir> 1012 M) and at late times, the accretion is dominated by a hot mode. For these haloes, the cosmological accreted gas falls into the dark-matter structure and is shock-heated to temperatures close to the virial temperature (Tvir). This gas participates in the development of a hot stable atmosphere (T> 106 K) around the central host galaxy.

Consequently, and as shown in Fig. 3, the baryonic cosmological accretion is divided into two phases (Mcold and Mhot). At a given redshift and for a given halo, the hot mass fraction (fsh(Mvir,z), shock heated) is computed following Lu et al. (2011) prescription (their Eqs. (24) and (25)): fsh(Mvir,z)=12[0.1×exp[(z4)2]+0.9]×[1+erf(logMvir11.40.4)]·\begin{eqnarray} \label{shock_heated_fraction} f_{\rm sh}(M_{\rm vir},z) = \dfrac{1}{2}\left[0.1\times \exp\left[-\left(\frac{z}{4}\right)^2\right]+ 0.9\right]~~\nonumber\\ ~~~\times \left[1+{\rm erf}\left(\dfrac{\log M_{\rm vir}-11.4}{0.4}\right)\right]\cdot \end{eqnarray}(14)Figure 4 gives the evolution of fsh with redshift and mass. As explained previously, the cold mode is dominated in the smallest structures (Mvir< 1011 M), and thus fsh(Mvir,z) ≪ 1. Unlike for these small objects, the cosmological accretion onto the largest haloes is dominated by the hot mode, and fsh(Mvir,z) ≃ 1. It is interesting to note that, at high redshift (z> 4), even for the largest haloes (Mvir> 1012 M), a small fraction (10%) of the accreted gas is always attributed to the cold mode. This residual cold fraction decreases with time and tends to zero in the local Universe.

thumbnail Fig. 4

Fraction of hot isotropic gas distributed in the baryonic accretion following (Lu et al. 2011). This fraction depends on the dark-matter halo mass and redshift (colour code).

thumbnail Fig. 5

The total smooth baryonic accretion onto the galaxy. The colour shading shows the total baryonic accretion rate (free-fall coming from filaments and cooling flows) from our reference model, as a function of dark-matter halo mass, and for a set of redshifts. The colour scale indicates the normalised logarithmic density of objects. We add for information the transition mass (1011.4 M, grey vertical arrow) used in the accretion mode repartition (Eq. (14)). The dotted and solid black lines show the average baryonic accretion rate deduced from hydrodynamic simulations by Fakhouri et al. (2010) and Ceverino et al. (2010), respectively. Black circles are for a model with photoionisation prescription by Gnedin (2000) instead of our Okamoto et al. (2008) standard prescription. Grey squares show the accretion for a model without any photoionisation process.

3.2.1. Cold gas accretion, the filamentary part

Before feeding the galaxy, the cold gas falls in the dark-matter potential well and forms some filamentary structures. This cold gas follows a very collimated distribution. The formation rate of these filamentary structures (flt) is deduced from the dark-matter accretion rate (acc, Eq. (1)) and from the fraction of shock-heated gas (fsh, Eq. (14)). We define flt=(1fsh)bar.\begin{eqnarray} \dot{M}_{\rm flt} = (1 - f_{\rm sh})\dot{M}_{\rm bar} . \label{flt_formation_rate} \end{eqnarray}(15)The cold gas is represented by the cold reservoir Mcold (Fig. 3). Once collimated in the filamentary structure, the cold gas falls into the centre of the dark-matter halo, and feeds a galaxy disc with a rate ff (see Fig. 3) close to the free-fall rate (ff). We define ff=εffMcold2tdyn,\begin{eqnarray} \dot{M}_{\rm ff} = \varepsilon_{\rm ff}\dfrac{M_{\rm cold}}{2t_{\rm dyn}} , \label{free-fall-rate} \end{eqnarray}(16)where εff = 1.0 is an efficiency parameter.

The evolution terms of the cold gas reservoir are thus completely known. They follow dMcolddt=cold=fltff.\begin{eqnarray} \dfrac{{\rm d}M_{\rm cold}}{{\rm d}t} = \dot{M}_{\rm cold} = \dot{M}_{\rm flt} - \dot{M}_{\rm ff}. \label{cold_mass_evolution} \end{eqnarray}(17)

3.2.2. Hot gas accretion, the shock-heated part

The other part of the baryonic cosmological accretion on a given dark-matter halo is an isotropically distributed gas that is shock-heated when its falls in the potential well. The shock is localised around the virial radius (rvir) of the dark-matter structure. The accretion rate (sh) of this hot isotropic gas phase is given by sh=fshbar.\begin{eqnarray} \dot{M}_{\rm sh} = f_{\rm sh}\dot{M}_{\rm bar} . \label{hot_cosmo_acc} \end{eqnarray}(18)At any given time, the total hot gas mass (Mhot in Fig. 3) located around the galaxy is the sum of this cosmological shock-heated gas and the ejected gas produced by the galaxy wind (wind in Fig. 3 and Sect. 4.1).

In summary, the time evolution of the hot atmosphere follows the differential equation dMhotdt=hot=sh+windcoolIGM,\begin{eqnarray} \dfrac{{\rm d}M_{\rm hot}}{{\rm d}t} = \dot{M}_{\rm hot} = \dot{M}_{\rm sh} + \dot{M}_{\rm wind} - \dot{M}_{\rm cool} - \dot{M}_{\rm IGM}, \label{hot_mass_evolution} \end{eqnarray}(19)where sh is given by Eq. (18), wind is the ejecta rate (wind) produced by the galaxy (see Sect. 4.1), and cool is the cooling rate computed by taking the radiative cooling of the hot gas into account (see Sect. 4.6). Finally, IGM is an output rate corresponding to the mass fraction (wind + pre-existing hot gas phase), which leaves the hot atmosphere due to its high velocity (see Sect. 4.3 for more details). Each term of this equation is shown in Fig. 3.

3.2.3. Total baryonic accretion on the galaxy

Figure 5 shows the total baryonic accretion on the galaxy (cool + ff) for a set of redshifts (z ∈ [ 0:6 ]). The grey vertical arrow indicates the transition mass (1011.4 M) used in Lu et al. (2011) (Eq. (14)) for the accretion mode repartition. For haloes with mass lower (higher) than this threshold the accretion is dominated by cold mode (cooling process). The colour histogram is built with our reference model. This reference model uses Okamoto et al. (2008) photoionisation prescription (see Sect. 3.1). The decrease in accretion efficiency for massive haloes is clearly visible. For comparison, we plot for haloes with mass lower than Mvir = 1011.4 M, the average accretion on the galaxy computed in the case of Gnedin (2000) photoionisation prescription. As explained in Cousin et al. (2015) the Gnedin (2000) photoionisation prescription is more efficient, and in this case, the average accretion is lower than with our reference model. We also show in the figure the mean accretion on the galaxy in a model without any photoionisation. At high redshift (z> 2), the difference between our reference model and the model without any photoionisation is very small. For z< 2, the differences are larger, especially for the smallest structures (by a factor 10 at z = 1). Again, this confirms that the Okamoto et al. (2008) photoionisation prescription only has a real impact on the smallest structures (Mvir< 1010 M) and at low redshift. Finally, we plot the mean baryonic accretion computed in hydrodynamic simulations by Ceverino et al. (2010, their Eq. (7), ) and Fakhouri et al. (2010, their Eq. (2), ΩbΩmmedian\hbox{$\frac{\Omega_{\rm b}}{\Omega_{\rm m}}\langle\dot{M}\rangle_{\rm median}$}) for comparison. We see that our model is in good agreement with these two simulations.

3.2.4. Some remarks concerning the bimodal accretion prescription

Our accretion model is based on an explicit bimodal baryonic accretion. This prescription allows us to follow the mass of cold and/or hot gas evolving around the galaxy. Indeed, even if the original calculation proposed by White & Frenk (1991) included the cooling and the free-fall timescales (Benson & Bower 2011), it was difficult to follow the mass of gas associated with each of these modes. In the context of our explicit two-phase accretion, it is easier to identify the different galaxy populations associated with a given accretion mode. The total fresh gas accreted on the galaxy (ff+Mcool˙\hbox{$\dot{M}_{\rm ff} + \dot{M_{\rm cool}}$}) in our model has been compared with the prediction from the original calculation proposed by White & Frenk (1991). It differs by less than 10%.

Some recent works based on a new computational hydrodynamic technique (e.g. Nelson et al. 2013), are questioning the existence of cold streams. It seems that the different computational methods (SPH, AMR, or moving mesh) do not give the same results. In parallel, the impact of the numerical resolution is also very important. At the time of writing, the debate is still open.

4. The hot halo phase

We have presented the hot gas accretion mode that is the first input term of the hot gas phase. In a classical galaxy evolution model, SNe and the SMBH activity heat a fraction of the gas. This gas can be ejected from the galaxy (wind) and thus has to be added to the hot gas phase. It is considered as the second input term of the hot gas phase. This section describes the two distinct feedback processes, SN and SMBH. All symbols used in the equations of this section and their definition are summarised in Table 3.

Table 3

Ejecta and heating processes symbols and their definition.

We present in the following (Sects. 4.1 to 4.6) a self-consistent model of the hot atmosphere. We begin (Sect. 4.1) by a complete description of the ejecta rate due to SNe and/or AGN. We continue (Sect. 4.2) by a description of the kinetic and thermal energy carried by the ejected mass. We add to this energy calculation a model of energy transfer from the ejecta to the pre-existing hot gas phase. These energy considerations allow us to compute (Sect. 4.3) an evolving hot halo temperature and the fraction of mass that can stay in equilibrium in the dark-matter halo potential well. The other part of the mass, which leaves the hot atmosphere, is removed from the hot gas reservoir Mhot and is transferred to (see Fig. 3):

  • A passive IGM reservoir if the structure is identified as a main halo,

  • the hot gas reservoir of the host M-H if the structure is a S-H (satellite).

We give the formalism for this hot atmosphere model for two different equations of state. In our reference model, we use the ideal gas case (Sect. 4.5.1). We also give the equations for the polytropic gas case in Appendix B.

SNe and SMBH (see Figs. 8 and 11) generate outflows (wind in Fig. 3). The total galaxy ejecta rate is given by the following relations: wind=wind,sn+wind,\begin{eqnarray} \dot{M}_{\rm wind} = \dot{M}_{\rm wind,sn}+ \dot{M}_{\rm wind,\BH} \label{gal_ejecta} \end{eqnarray}(20)We develop these terms in the following section.

4.1. Ejecta processes

4.1.1. SNe ejecta

The instantaneous ejecta rate produced by SNe (wind,sn) is deduced with the principle of energy conservation. This has already been used (see Dekel & Silk 1986; Kauffmann et al. 1993) and Efstathiou (2000) for a review. We link the kinetic energy per time unit to the power produced by SNe following 12wind,sn,Vwind2=εsnfKin,snEsnηsn,\begin{eqnarray} \dfrac{1}{2}\dot{M}_{\rm wind,sn,}V_{\rm wind}^2 = \varepsilon_{\rm sn}f_{\rm Kin,sn}E_{\rm sn}\eta_{\rm sn}\dot{M}_{\star} , \label{sn_energy_conservation} \end{eqnarray}(21)where

  • (1 − εsn)is the fraction of SNe kinetic energy used to feed the disc velocity dispersion (see Sect. 5.3.1);

  • fKin,snEsnis the SNe energy converted in kinetic energy. We use Esn = 1044 Joules and fKin,sn = 0.3 (Kahn 1975; Aguirre et al. 2001);

  • ηsnis the number of SNe produced per unit of stellar mass ([ηsn]=M-1\hbox{$[\eta_{\rm sn}] = \Msun^{-1}$}). This parameter is linked to the initial mass function IMF2,

  • is the star formation rate.

This formulation of the kinetic energy per unit of time leads to a degeneracy between the ejected mass rate (wind,sn) and the wind velocity (Vwind). To break the degeneracy, we compute the wind velocity using a model extracted from Bertone et al. (2005) and based on Efstathiou (2000) and Shu et al. (2005). The wind velocity can be linked to the star formation rate (Martin 1999) and seems to be independent on the galaxy morphology (Heckman et al. 2000; Frye et al. 2002). We therefore use the following relation (Bertone et al. 2005, their Eq. (9)) Vwind=623(100Myr-1)0.145[kms-1].\begin{eqnarray} V_{\rm wind} = 623\left(\dfrac{\dot{M}_{\star}}{100~\Msun\,{\rm yr}^{-1}}\right)^{0.145}~~\left[{\rm km\,s}^{-1}\right] . \label{Vwind} \end{eqnarray}(22)

4.1.2. Super-massive-black hole activity

An accreting SMBH produces wind. The outflow rate (•,out) is linked to the effective accretion rate on the SMBH (•,acc) by the relation (Ostriker et al. 2010) η=,out,acc·\begin{eqnarray} \eta_{\BH}=\dfrac{\dot{M}_{\rm \BH,out}}{\dot{M}_{\rm \BH,acc}} \cdot \label{agn_out_acc} \end{eqnarray}(23)This parameter η is currently distributed in the range [ 0.1:1 ] (Ostriker et al. 2010). In our reference model we use η = 0.6. If •,inf is the total infall rate (see Sect. 6 and Eq. (88) for explicit prescriptions), the conservation of the mass gives obviously ,inf=,acc+,out,\begin{eqnarray} \dot{M}_{\rm \BH,inf} = \dot{M}_{\rm \BH,acc} + \dot{M}_{\rm \BH,out} , \label{agn_mass_conservation} \end{eqnarray}(24)and therefore ,acc=,inf1+ηand,out=,infη1+η·\begin{eqnarray} \dot{M}_{\rm \BH,acc} = \dfrac{\dot{M}_{\rm \BH,inf}}{1 + \eta_{\BH}}~~\mbox{and}~\dot{M}_{\rm \BH,out} = \dfrac{\dot{M}_{\rm \BH,inf}\eta_{\BH}}{1 + \eta_{\BH}}\cdot \label{agn_in_out_rate} \end{eqnarray}(25)Using the same energy conservation criterion as that used for SNe ejecta (see Sect. 4.1.1 and Eq. (21)), we obtain 12,outVjet2=fKin,,accc2,\begin{eqnarray} \dfrac{1}{2}\dot{M}_{\rm \BH,out} V_{\rm jet}^2 = f_{\rm Kin,\BH}\dot{M}_{\rm \BH,acc}c^2 , \label{agn_energy_conservation} \end{eqnarray}(26)where

  • fKin,•is the fraction of energy produced by the accretion on the SMBH that is converted into kinetic energy. This parameter is not well known but observations and numerical simulations give a range of plausible values, fKin,• ∈ [ 10-4:10-3 ] (Ostriker et al. 2010; Proga et al. 2000; Stoll et al. 2009). In our reference model we use fKin,• = 10-3. This value is compatible with observed outflow velocities, as given for example by Emonts et al. (2005) and Morganti et al. 2005a,b.

  • Vjetis the outflow velocity.

Observations of ionised and neutral gas outflows in radio galaxies show that the SMBH activities and, more precisely, the jet have an impact on the gas content (e.g. Morganti et al. 2005a,b; Emonts et al. 2005; Nesvadba et al. 2006, 2008; Lehnert et al. 2011; Guillard et al. 2012). However, it is still unclear how the gas is affected. To compute the galaxy ejecta rate linked to the SMBH activity (wind,•), we use the momentum conservation between the SMBH jet and the gas in the galaxy. We assume that a fraction of the gas can reach a velocity equal to the galaxy escape velocity (\hbox{$V_{\rm esc,\Gal}$}) thanks to momentum transfer. The total ejecta coming from the SMBH activity is therefore wind,=,infη1+η[1+(cεVesc,𝒢)2fKin,η]\begin{eqnarray} \begin{tiny} \dot{M}_{\rm wind,\BH} = \dfrac{\dot{M}_{\rm \BH,inf}\eta_{\BH}}{1 + \eta_{\BH}}\left[1 + \left(\dfrac{c}{\varepsilon_{\BH}V_{\rm esc,\Gal}}\right)\sqrt{\dfrac{2 f_{\rm Kin,\BH}}{\eta_{\BH}}}\right] \label{agn_tot_ejecta_rate} \end{tiny} \end{eqnarray}(27)where

  • ε = 0.6is the coupling factor between the jet momentum and the gas reservoir, This factor is adjusted to produce an outflow rate of a few solar masses per year during the secular evolution of galaxies and a few tens of solar masses per year during the merger induced activity (e.g. Morganti et al. 2005a,b; Emonts et al. 2005).

  • \hbox{$V_{\rm esc,\Gal}$}is the escape velocity of the galaxy (Eq. (56)).

4.2. Energy of the ejecta, transfers to the hot gas phase

As developed in subsequent sections, we introduce in our model an explicit description of the energetic content of the hot gas surrounding the galaxy. Indeed, accretion and ejecta processes generate energy. A fraction of this energy (kinetic, thermal, and luminous) leaves the galaxy and is distributed in the hot halo gas phase. This energy contributes to the heating of the hot gas.

We follow the internal energy of

  • the hot gas evolving around the galaxy: ETh,atm (28);

  • the hot accreted gas: ETh,acc (29);

  • the wind generated by the galaxy: ETh,wind (32).

With these three terms, it is possible to follow the evolution of the mean temperature of the hot gas phase T\hbox{$\overline{T}$} and estimate the mean temperature Twind\hbox{$\overline{T}_{\rm wind}$} of the wind produced by the galaxy.

In this section, we separately discuss the three terms that are considered as three different internal energy sources. The properties of the hot gas are then computed using the following steps:

  • 1.

    The internal energy of the accreted gas ETh,acc is added to the internalenergy of the pre-existing hot gas ETh,atm.

  • 2.

    We then deduce a mean temperature Tatm\hbox{$\overline{T}_{\rm atm}$} for the pre-existing hot gas phase.

  • 3.

    We deduce a mean temperature Twind\hbox{$\overline{T}_{\rm wind}$} for the wind phase from its internal energy (ETh,wind).

  • 4.

    We deduce a mean temperature Twind\hbox{$\overline{T}_{\rm wind}$} for the wind phase, by taking these two mean temperatures (Twind\hbox{$\overline{T}_{\rm wind}$} and Tatm\hbox{$\overline{T}_{\rm atm}$}) into account in a distribution function, f. We finally compute the escape mass fraction fesc and the effective mean temperature T\hbox{$\overline{T}$}.

4.2.1. Thermal energy of the pre-existing hot phase

At a given time t, the gas mass in the hot pre-existing atmosphere is Mhot. This gas is considered in thermal equilibrium in the dark-matter potential well. If the mean temperature of the gas3 is T\hbox{$\overline{T}$}, and removing the mass that has condensed during the time step (coolΔt, Eq. (54)), the total internal energy of the pre-existing hot gas phase at the end of the time step (t + Δt) is given by ETh,atm=3kBT2μmp[Mhot(t)coolΔt􏽼0􏽻􏽺0􏽽condensedmass].\begin{eqnarray} E_{\rm Th,atm} = \dfrac{3k_{\rm B}\overline{T}}{2\mu m_{\rm p}}\left[M_{\rm hot}(t)-\underbrace{\dot{M}_{\rm cool}\Delta t}_{\rm condensed~mass}\right]. \label{Ethatm} \end{eqnarray}(28)

4.2.2. Thermal energy in the gas accretion

As explained in Sect. 3.2.2, the hot isotropic atmosphere is also fed by cosmological hot (shock-heated) accretion. The thermal energy of the gas accreted during the previous time step Δt is ETh,acc=3kBTvir2μmp[shΔt􏽼0􏽻􏽺0􏽽accretedmass],\begin{eqnarray} E_{\rm Th,acc} = \dfrac{3k_{\rm B} T_{\rm vir}}{2 \mu m_{\rm p}}\left[\underbrace{\dot{M}_{\rm sh}\Delta t}_{\rm accreted~mass}\right] , \label{Ethacc} \end{eqnarray}(29)where Tvir is the virial temperature of the dark-matter halo (Eq. (8)) and shΔt (Eq. (18)) is the cosmological hot gas mass accreted during the time step Δt.

4.2.3. Thermal energy in the ejecta

During a SNe explosion, a fraction fKin,sn of the total energy (Esn) is converted into kinetic energy. The remaining fraction is distributed in luminous energy (photons) and in thermal energy in the ejected gas. We can write the instantaneous non-kinetic power produced by SNe: sn=(1fKin,sn)Esnηsn.\begin{eqnarray} \dot{Q}_{\rm sn} = (1-f_{\rm Kin,sn})E_{\rm sn}\eta_{\rm sn}\dot{M}_{\star} . \label{Qsn} \end{eqnarray}(30)The instantaneous non-kinetic power produced by SMBH activity follows =,accc2(1fKin,)=,infη1+ηc2(1fKin,).\begin{eqnarray} \label{Qagn} \dot{Q}_{\BH} & =& \dot{M}_{\rm \BH,acc}c^2\left(1-f_{\rm Kin,\BH}\right)\nonumber\\ & =& \dot{M}_{\rm \BH,inf}\dfrac{\eta_{\BH}}{1 + \eta_{\BH}}c^2\left(1-f_{\rm Kin,\BH}\right) . \end{eqnarray}(31)Equations (30) and (31) are used to compute the instantaneous power. We then compute the (time-)average on the evolution time step Δt of

  • the non-kinetic power, \hbox{$\langle\dot{Q}_{\rm sn}\rangle$};

  • the SMBH non-kinetic power, \hbox{$\left<\dot{Q}_{\BH}\right>$}.

With these two quantities, we can derive the thermal energy carried by the galactic winds during Δt following ETh,wind=Δt[fTh,snsn+fTh,agn],\begin{eqnarray} E_{\rm Th,wind} = \Delta t\left[f_{\rm Th,sn}\left<\dot{Q}_{\rm sn}\right> + f_{\rm Th,agn}\left<\dot{Q}_{\BH}\right>\right] , \label{Ethwind} \end{eqnarray}(32)where we define fTh,sn and fTh,• as two free parameters describing the fraction of the non-kinetic energy distributed in thermal energy by SNe and SMBH feedback, respectively. The reference values are set to fTh,sn = 0.1 and fTh,• = 10-3. These values have been chosen so as to produce a temperature distribution centred on:

  • 106K for SN wind. Indeed, typical SN remnant emission shows temperatures close to kBT ∈ [ 0.1:0.7 ] keV, which correspond to T ∈ [ 1:8 ] 106 K (Koo et al. 2002; Sasaki et al. 2014);

  • 107K for SMBH wind.

4.2.4. Temperatures

During a time step Δt, the hot gas phase is fed by accretion, cooling, and ejecta processes. In this context, the mean temperature of the hot atmosphere after a time step Δt is given by Tatm=2μmp3kB[ETh,atm+ETh,accMhot(t)+shΔt􏽼0􏽻􏽺0􏽽accretedmasscoolΔt􏽼0􏽻􏽺0􏽽condensedmass]=T(Mhot(t)coolΔt)+TvirshΔtMhot(t)+shΔtcoolΔt·\begin{eqnarray} \overline{T}_{\rm atm} & = &\dfrac{2\mu m_{\rm p}}{3k_{\rm B}}\left[\dfrac{E_{\rm Th,atm} + E_{\rm Th,acc}}{M_{\rm hot}(t)+\underbrace{\dot{M}_{\rm sh}\Delta t}_{\rm accreted~mass}-\underbrace{\dot{M}_{\rm cool}\Delta t}_{\rm condensed~mass}}\right]\nonumber\\ \label{Tatm} & =& \dfrac{\overline{T}\left(M_{\rm hot}(t)-\dot{M}_{\rm cool}\Delta t\right)+ T_{\rm vir}\dot{M}_{\rm sh}\Delta t}{M_{\rm hot}(t)+\dot{M}_{\rm sh}\Delta t-\dot{M}_{\rm cool}\Delta t}\cdot \end{eqnarray}(33)In parallel we consider that the temperature of the wind is constant during the time step Δt. It is given by Twind=2μmp3kBETh,windMwind􏽼0􏽻􏽺0􏽽ejectedmass·\begin{eqnarray} \overline{T}_{\rm wind}= \dfrac{2\mu m_{\rm p}}{3k_{\rm B}}\dfrac{E_{\rm Th,wind}}{\underbrace{M_{\rm wind}}_{\rm ejected~mass}}\cdot \label{Twind} \end{eqnarray}(34)

4.3. Escape fraction and mean temperature

4.3.1. Average velocity and mass of the galaxy wind

On the galaxy scale, the instantaneous total wind velocity (SN and SMBH) is computed using velocity winds of each component weighted by their ejecta rate. We define Vwind=wind,snVwind+wind,Vesc,𝒢wind,sn+wind,·\begin{eqnarray} \small{ V_{\rm wind} = \dfrac{\dot{M}_{\rm wind,sn}V_{\rm wind}+\dot{M}_{\rm wind,\BH}V_{\rm esc,\Gal}}{\dot{M}_{\rm wind,sn}+\dot{M}_{\rm wind,\BH}}\cdot \label{Vwind_gal}} \end{eqnarray}(35)This instantaneous wind velocity is only valid on the galaxy scale. As described in Eq. (A.4) the evolution of the galaxy components (disc and bulge) is performed with an adaptative time step δt shorter than the hot atmosphere time step Δt = ∑ iδti. To take the impact of the hot gas wind on the hot gas atmosphere into account, such as the permanently lost fraction (see Sect. 4.3 and Eq. (41)), it is necessary to estimate the average wind velocity for the time scale used for the hot phase. The wind velocity is therefore computed following Vwind=1Δti=1NVwind(iδt)δt.\begin{eqnarray} \small{ \ds \left<V_{\rm wind}\right> = \dfrac{1}{\Delta t}\sum_{i=1}^N V_{\rm wind}(i\delta t)\delta t . \label{Vwind_hot}} \end{eqnarray}(36)The total mass contained in the wind is therefore given by Mwind=i=1Nwind(iδt)δt.\begin{eqnarray} \small{ \ds M_{\rm wind} = \sum_{i=1}^N \dot{M}_{\rm wind}(i\delta t)\delta t . \label{Mwind_hot}} \end{eqnarray}(37)

4.3.2. Velocity distributions

We assume that the hot atmosphere in which the wind extends is in hydrostatic equilibrium in the dark-matter halo potential well. The hot gas phase has a mean temperature Tatm\hbox{$\overline{T}_{\rm atm}$}. In these conditions, the velocity probability distribution of a particle in this hot gas phase is given by the well-known Maxwell-Boltzman distribution, fMB(v,Tatm)dv\hbox{$f_{\rm MB}(v,\overline{T}_{\rm atm}){\rm d}v$} (see Eq. (C.1)).

Concerning the wind, in the gas frame of reference, the probability velocity distribution is also given by a Maxwell-Boltzman distribution, fMB(v,Twind)dv\hbox{$f_{\rm MB}(v,\overline{T}_{\rm wind}){\rm d}v$}. In the fixed halo (or galaxy) referential, if we take the mean velocity of the wind into account, the velocity probability distribution is also given by the conditional following relation f(v,Twind)dv={\begin{eqnarray} f(v,\overline{T}_{\rm wind}){\rm d}v = \left\{ \begin{array}{ll} f_{\rm MB}\left(v ',\overline{T}_{\rm wind}\right){\rm d}v {:} & {\rm if}\quad v' > 0\\ 0{:} & \mbox{otherwise,} \end{array}\right. \end{eqnarray}(38)with v′ = v − ⟨Vwind.

At any given time, we can consider that the hot gas phase is composed of a pre-existing hot atmosphere in which the hot gas produced by feedback processes is injected. Therefore, the global velocity probability distribution is given by, f(v,Tatm,Twind)dv=(1fwind)fMB(v,Tatm)dv\begin{eqnarray} f^{\dagger}(v,\overline{T}_{\rm atm},\overline{T}_{\rm wind}){\rm d}v &= & \,(1-f_{\rm wind})f_{\rm MB}(v,\overline{T}_{\rm atm}){\rm d}v\nonumber\\ &&+f_{\rm wind}f(v,\overline{T}_{\rm wind}){\rm d}v \label{F_dag} \end{eqnarray}(39)with fwind=MwindMhot+Mwind+shΔt􏽼0􏽻􏽺0􏽽accretedmasscoolΔt􏽼0􏽻􏽺0􏽽condensedmass·\begin{eqnarray} f_{\rm wind} = \dfrac{M_{\rm wind}}{M_{\rm hot}+M_{\rm wind}+\underbrace{\dot{M}_{\rm sh}\Delta t}_{\rm accreted~ mass}-\underbrace{\dot{M}_{\rm cool}\Delta t}_{\rm condensed~mass}}\cdot \label{fwind} \end{eqnarray}(40)

4.3.3. The effective escape fraction

At a given time some particles in the pre-existing hot gas phase and even more in the wind phase can have higher velocities than the escape velocity of the dark-matter structure Vesc,dm (Eq. (9)). We assume that only the dark-matter mass and therefore the dark-matter gravitational well have a direct impact on the hot gas phase confinement. Suto et al. (1998) demonstrate that the self-gravity of the hot gas phase has an impact on the hot gas phase density profile but only at large radius (r>rdm). The density profile of the hot gas is an important quantity in the gas cooling modelling, but if we consider as Suto et al. (1998) that only the core part of the hot atmosphere can cool in a short time, the impact of not taking the self-gravity on cooling processes into account is small.

The mass fraction fesc that can leave the dark-matter potential well is given by fesc=10Vesc,dmf(v,Tatm,Twind)dv.\begin{eqnarray} f_{\rm esc} = 1 - \int_0^{V_{\rm esc,dm}}f^{\dagger}\left(v,\overline{T}_{\rm atm},\overline{T}_{\rm wind}\right){\rm d}v . \label{f_esc} \end{eqnarray}(41)Figure 6 shows the distribution of f(v,Tatm,Twind)dv\hbox{$f^{\dagger}(v,\overline{T}_{\rm atm},\overline{T}_{\rm wind}){\rm d}v$} (Eq. (39)) for two different dark-matter halo masses. In the two cases the wind velocity is close to Vwind = 250 km s-1 but for Mvir = 1010 M, the temperature of the wind is close to T=106\hbox{$\overline{T} =10^6$} K, while for Mvir = 1013 M, it is close to T=107\hbox{$\overline{T} =10^7$} K. As expected, the fraction of the wind mass that can definitively leave the structure is larger for the lowest dark-matter halo mass. Indeed, for such small structures, all the mass contained in the wind leaves the potential well.

thumbnail Fig. 6

Distribution functions of the velocity in the pre-existing hot gas phase and in the wind. The upper and lower panels show the distributions for two different dark-matter haloes with Mvir = 1010 M and Mvir = 1013 M, respectively. In the two cases, the wind speed is close to Vwind = 250 km s-1 but the temperature is close to T=106\hbox{$\overline{T} =10^6$} K in the upper panel and T=107\hbox{$\overline{T}=10^7$} K in the lower panel. The total gas distribution, in black, is the sum of the pre-existing hot gas distribution (red dotted line) and the gas wind distribution (grey dotted line). The blue curve gives the mass fraction. The vertical dotted grey line indicates the escape velocity of the dark-matter structure. The intersection of the mass fraction and this velocity threshold allows us to determine the escape fraction. The red curve shows the new hot gas phase distribution recomputed with the new temperature T\hbox{$\overline{T}$} (Eq. (42)) derived after accretion, cooling, and ejection processes have been taken into account during the previous time step Δt. In the low-mass case, all the mass contained in the wind leaves the structure. In the high-mass case, a fraction of the wind mass can be retained.

After ejection and evaporation, the distribution and the average temperature are recomputed. In Fig. 6 the new f(v,Tatm,Twind)dv\hbox{$f^{\dagger}(v,\overline{T}_{\rm atm},\overline{T}_{\rm wind}){\rm d}v$} distribution is shown in red. For the halo with Mvir = 1013 M, the temperature has increased because the hot wind has participated in heating the hot isotropic atmosphere.

4.3.4. Mass transfer and mass loss

The mass MIGM = fescMhot leaves the hot halo owing to its high velocity. If the galaxy that generates the outflows is hosted by an S-H (satellite), we recompute the ejection processes, following the previous description, in the host M-H referential. Indeed, the halo finder algorithm creates links between S-H and their host M-H. We compute the escape velocity of M-H, and deduce the fraction of ejected mass that can leave its dark-matter potential well. The difference between the mass lost by S-H and M-H is added to the hot reservoir of M-H. The fraction of the mass that can leave M-H is definitively removed from the hot reservoir and added to a passive IGM reservoir (see Fig. 3).

After mass ejection, the average temperature of the gas is modified. This new mean temperature is then given by T=2μmp3kB(1fesc)0Vesc,dmv2f(v,Tatm,Twind)dv.\begin{eqnarray} \overline{T} = \dfrac{2\mu m_{\rm p}}{3k_{\rm B}\left(1-f_{\rm esc}\right)}\int_0^{V_{\rm esc,dm}}v^2f^{\dagger}\left(v,\overline{T}_{\rm atm},\overline{T}_{\rm wind}\right){\rm d}v. \label{T_hot} \end{eqnarray}(42)

4.4. Some remarks concerning the SMBH feedback efficiency

SMBH feedback has been implemented in different models (e.g. Croton et al. 2006; Cattaneo et al. 2006; Bower et al. 2006; Malbon et al. 2007; Somerville et al. 2008; Ostriker et al. 2010; Bower et al. 2012). For example, in Croton et al. (2006; their Eqs. (10) and (11)) and/or Somerville et al. (2008; their Eq. (21)) SMBH activity directly affects the cooling of the hot gas phase. Indeed in these two models, a fraction of the power generated by the accretion is considered as a heating gas rate heat. This heating rate is then used to compute an effective cooling rate, cool,eff = MAX(0,cool,effheat). This approach leads to a very efficient SMBH feedback in which a large fraction of the power generated by the accretion is directly used to reduce the cooling process.

In our approach, SMBH feedback contributes to gas ejection (Eq. (27)), but the impact on the cooling rate is only given by a possible increase in the hot gas phase temperature due to wind (thermal energy in the wind, Eqs. (31) and (32)). In this context, even if the SMBH activity produces high-temperature wind (107 K), our SMBH feedback implementation leads to less efficiency and therefore has a weaker impact on the cooling process.

4.5. The hot gas density profile

To describe the hot isotropic atmosphere located around the galaxy and to take the impact of the cooling process into account, we need to define a density profile. This section is dedicated to describing of i) the density profile and ii) the cooling process. All the parameters and their definitions are listed in Table 4.

Table 4

Symbols and their definition for the hot atmosphere and cooling process.

We first assume that the hot stable mass Mhot (in hydrostatic equilibrium) is enclosed in the virial radius rvir of the dark-matter structure. To model the mass distribution of the hot gas we assume a spherical geometry and thus define Mhot=4π0rvirρg(r)r2dr.\begin{eqnarray} M_{\rm hot} = 4\pi \int_0^{r_{\rm vir}}\rho_{\rm g}(r)r^2{\rm d}r . \label{Mhot_r} \end{eqnarray}(43)The density distribution ρg(r) of the hot atmosphere is deduced from the hydrostatic equilibrium condition (HEC) computed in a potential dominated by the dark-matter, 1ρg(r)dP(r)dr=GMh(<r)r2,\begin{eqnarray} \dfrac{1}{\rho_{\rm g}(r)}\dfrac{{\rm d}P(r)}{{\rm d}r} = -\dfrac{G M_{\rm h}({<}r)}{r^2} , \label{HEC_id} \end{eqnarray}(44)where G is the gravitational constant, P(r) the pressure radial profile and Mh( <r) the dark-matter mass enclosed in the radius r. In our model it is possible to use two different density profiles deduced from an ideal or a polytropic gas. In the reference model, we use the isothermal gas case. For completeness we describe the polytropic gas case in Appendix B.

4.5.1. The isothermal ideal gas case

For an ideal (id) gas, the equation of state links the pressure, the density, and the temperature by the well-known following relation Pid(r)=kBT0μmpρgid(r),\begin{eqnarray} P^{\rm id}(r)=\frac{k_{\rm B} T_0}{\mu m_{\rm p}}\rho_{\rm g}^{\rm id}(r) , \label{T0} \end{eqnarray}(45)where kB is the Boltzman constant, μ the mean molecular weight of the gas, mp the proton mass, Pid(r) the radial pressure profile, ρgid(r)\hbox{$\rho_{\rm g}^{\rm id}(r)$} the density profile, and T0 the central temperature of the gas. In this simple case of an ideal gas we assume, that the radial profile of the temperature T(r) is set to a constant value and that this constant value is equal to the central value T0. Here, T0 is equal to the average value Tatm\hbox{$\overline{T}_{\rm atm}$}. The HEC in the dark-matter halo potential (Suto et al. 1998; Makino et al. 1998; Komatsu & Seljak 2001; Capelo et al. 2012) gives dln(ρgid)dr=GμgmpkBT0Mh(<r)r2=4πGμmpρdmrdmkBT0φ(x)x2,\begin{eqnarray} \dfrac{{\rm dln}(\rho_{\rm g}^{\rm id})}{{\rm d}r} = -\dfrac{G\mu_{\rm g} m_{\rm p}}{k_{\rm B}T_0}\frac{M_{\rm h}({<}r)}{r^2} = -\dfrac{4\pi G\mu m_{\rm p}\rho_{\rm dm}r_{\rm dm}}{k_{\rm B}T_0}\frac{\phi(x)}{x^2} , \end{eqnarray}(46)where ρdm and rdm are the core density and the characteristic radius of the dark-matter halo, respectively and φ(x) is the geometrical function given in Eq. (5).

For an NFW profile (Navarro et al. 1995, 1996, 1997), the hot gas profile that follows the HEC is Suto et al. (1998); Capelo et al. (2012)ρgid(r)=ρg,0id(x)withx=rrdm,\begin{eqnarray} \rho_{\rm g}^{\rm id}(r) = \rho_{\rm g,0}\mathcal{F}^{\rm id}(x)~~\mbox{with}~~x=\frac{r}{r_{\rm dm}}, \label{rhog0} \end{eqnarray}(47)with id(x), a geometrical function defined as id(x)=exp[4πGμmpρdmrdm2kBT0Φ(x)],\begin{eqnarray} \mathcal{F}^{\rm id}(x) = \exp\left[-\dfrac{4\pi G\mu m_{\rm p}\rho_{\rm dm}r_{\rm dm}^2}{k_{\rm B}T_0}\Phi(x)\right] , \label{Fidx} \end{eqnarray}(48)where Φ(x) is another geometrical function based on the previous one φ(x) (Eq. (5)) and defined by Φ(x)=0xφ(u)u2du=1ln(1+x)x·\begin{eqnarray} \Phi(x)=\ds\int_{0}^{x}\dfrac{\phi(u)}{u^2}{\rm d}u = 1 - \dfrac{\ln(1+x)}{x} \cdot \label{Phix} \end{eqnarray}(49)In the case of an ideal gas, the geometrical distribution of the hot gas is completely set by the dark-matter properties. The last parameter that remains to be determined is the gas core density ρg,0. It is extracted from the total mass contained in the hot halo (Eq. (43)) as ρg,0=Mhot[4πrdm30cdmid(x)x2dx]-1.\begin{eqnarray} \rho_{\rm g,0} = M_{\rm hot}\left[4\pi r_{\rm dm}^3\int_0^{c_{\rm dm}} \mathcal{F}^{\rm id}(x)x^2{\rm d}x \right]^{-1}. \label{rho_gas_core_id} \end{eqnarray}(50)

4.6. Cooling process

The cooling and condensation of the hot atmosphere is computed using the classical model initially proposed by White & Frenk (1991). We recall the cooling time function here, which gives the cooling time of the gas mass enclosed in a radius r and which follows a given density ρg(r) and temperature profile T(r), t(r)=1.04μmpT(r)ρg(r)Λ[T(r),Zg]·\begin{eqnarray} t(r) = 1.04\dfrac{\mu m_{\rm p} T(r)}{\rho_{\rm g}(r)\Lambda\left[T(r),Z_{\rm g}\right]}\cdot \label{cooling_time_function} \end{eqnarray}(51)The cooling time function also depends on the mean molecular weight of the gas μ, the proton mass mp, and on the cooling function Λ that depends on the temperature T(r) and metallicity Zg of the hot gas. We use the cooling function Λ from Sutherland & Dopita (1993), for temperature between 104 K and 108.5 K and metallicity between 10-30Z and 100.5Z.

4.6.1. Cooling time, radius, and rate

To compute the cooling rate, and therefore the condensed mass of gas that can fall into the galaxy, we need to define the cooling radius, rcool which the radial extension of the hot atmosphere that can condense during a given time tcool. The cooling radius rcool is the solution of the following equation t(rcool)=tcool,\begin{eqnarray} t(r_{\rm cool}) = t_{\rm cool} \label{cooling_radius_equation} , \end{eqnarray}(52)where t(r) is the cooling time function (Eq. (51)).

Different definitions for tcool exist in the literature (e.g. Cole et al. 2000; Hatton et al. 2003), and thus for rcool (see De Lucia et al. 2010, for a discussion). In our model, we link tcool to a cooling clock:

  • The cooling clock starts when the hot reservoir Mhot receives some mass. The clock runs as long as the hot phase contains gas. Therefore, tcool increases with the evolving time of the hot atmosphere.

  • The cooling clock is set to 0 when a major merger event occurs. We consider a major merger event when the mass ratio of the two hot gas reservoirs (Mhot1\hbox{$M_{\rm hot}^1$} and Mhot2\hbox{$M_{\rm hot}^2$}) MIN ( M hot 1 , M hot 2 ) / MAX ( M hot 1 , M hot 2 ) > 1 / 3. \begin{eqnarray} {\rm MIN}(M_{\rm hot}^1,M_{\rm hot}^2)/{\rm MAX}(M_{\rm hot}^1,M_{\rm hot}^2) >1/3 . \end{eqnarray}(53)

  • After a minor merger, the cooling clock is set to the cooling time tcool of the most massive progenitor4.

The cooling rate of the hot atmosphere is then computed following cool=εcool2tdyn0MIN(rcool,rvir)ρg(r)r2dr.\begin{eqnarray} \dot{M}_{\rm cool} = \dfrac{\varepsilon_{\rm cool}}{2t_{\rm dyn}}\int_0^{{\rm MIN}(r_{\rm cool},r_{\rm vir})}\rho_{\rm g}(r)r^2{\rm d}r. \label{cooling_rate} \end{eqnarray}(54)The mass enclosed in the cooling radius rcool falls into the galaxy with a rate close to the free-fall rate. We assume that the hot atmosphere extends up to the virial radius, so that the cooling radius cannot be larger than this virial radius. The free efficiency parameter is set to εcool = 1.0.

5. Galaxy components

We present the relations and definitions linked to galaxy formation and evolution in the next sections. An eGalICS galaxy (Fig. 7) is composed of two distinct but interacting components: a disc \hbox{$\Disc$} (Fig. 8) and a bulge or bulge-like component (Fig. 11).

thumbnail Fig. 7

Exchanges on the galaxy scale. The galaxy scale is considered as an intermediate scale. A galaxy can be made of two different parts, a disc \hbox{$\Disc$} and a bulge . The fresh gas (ff and Mcool˙\hbox{$\dot{M_{\rm cool}}$}), coming from the cold and the hot reservoirs, only feeds the disc component. The bulge component is created during a merger or by the clump migration process. SNe and the activity of the SMBH produce ejecta (Mwind). These ejecta are considered at an upper scale (see Fig. 3). The gas exchanges associated with the disc and the bulge are described in Figs. 8 and 11.

thumbnail Fig. 8

Exchanges on the disc scale. The disc is the main component. It is fed by the fresh gas coming from the large scales (ff and cool). SNe associated with the stellar population of the disc generate ejecta wind,sn. In a disc, dynamics and gravity are in balance. Under some conditions, giant clumps can be formed (\hbox{$\dot{M}_{\rm \Clumps,in}$}) and/or disrupted (\hbox{$\dot{M}_{\rm \Clumps,out}$}). At any given time, a given fraction \hbox{$M_{\Clumps}$} of the total disc mass is in those clumps. The clumps migrate (\hbox{$\dot{M}_{\rm \Clumps,migr}$}) and feed a pseudo-bulge component ().

5.1. Accretion on the galaxy

At the first steps of galaxy formation, the cold gas (from filamentary collimated structures and/or the isotropic cooling flow) falls into the center of the dark-matter halo. We assume that this cold gas initially forms a thin exponential disc. The gas acquires angular momentum during the mass transfer (Peebles 1969). After its formation, the disc is supported by its angular momentum. This paradigm is based on the prescription given by Blumenthal et al. (1986) and Mo et al. (1998) and has been frequently used in SAMs (e.g. Cole 1991; Cole et al. 2000; Hatton et al. 2003; Somerville et al. 2008). During the secular evolution of the galaxy, some mass may be transferred to a pseudo-bulge component by disc instabilities (see Sects. 5.4.3 and 5.4.4 for more information). After a major merger event (see Sect. 8.3), the remnant galaxy has a spheroidal morphology (see Sect. 5.5). The gas contained in the remnant galaxy and the freshly accreted gas, which will fall into the centre of the halo in the next time step, will then form a new disc.

5.2. The circular and escape velocity of the galaxy

The radial profile of the circular velocity of the galaxy is deduced from the dynamic of all its components, the dark-matter halo \hbox{$\DM$}, disc \hbox{$\Disc$}, and bulge . We obtain Vcirc,𝒢2=Vcirc,𝒟ℳ2+Vcirc,𝒟2+Vcirc,2.\begin{eqnarray} V_{\rm circ,\Gal}^2 = V_{\rm circ,\DM}^2 + V_{\rm circ,\Disc}^2 + V_{\rm circ,\Bulge}^2 . \label{gal_circular_velocity} \end{eqnarray}(55)All the terms will be defined in their corresponding description section (Eqs. (10), (60), and (83)).

The escape velocity of the galaxy is computed following Vesc,𝒢2=2M𝒢(r<r𝒢,50)+Mvir(r<r𝒢,50)r𝒢,50,\begin{eqnarray} V_{\rm esc,\Gal}^2 = \dfrac{2M_{\Gal}\left(r<r_{\Gal,50}\right)+M_{\rm vir}\left(r<r_{\Gal,50}\right)}{r_{\Gal,50}} , \label{gal_escape_velocity} \end{eqnarray}(56)where \hbox{$M_{\Gal}(r<r_{\Gal,50})$} and \hbox{$M_{\rm vir}(r<r_{\Gal,50})$} are the galaxy and dark-matter mass enclosed in the galaxy half-mass radius \hbox{$r_{\Gal,50}$}, respectively.

5.3. The disc \hbox{$\Disc$}

Table 6 summarised the parameters of the disc. In the classical scenario adopted here, the disc component is the first structure formed in a new galaxy. We use the standard approach as proposed by Cole et al. (2000); Hatton et al. (2003), among others.

Table 5

Galaxy symbols and their definition.

Table 6

Disc symbols and their definition.

We assume that the disc component is infinitely thin. Its mass radial distribution is given by M𝒟(r)=M𝒟[1exp(x𝒟)(1+x𝒟)],\begin{eqnarray} M_{\Disc}(r) = M_{\Disc}\left[1-\exp(-x_{\Disc})\left(1+x_{\Disc}\right)\right] , \label{disc_profile} \end{eqnarray}(57)where \hbox{$x_{\Disc}$} is the dimensionless radius \hbox{$x=r/r_{\Disc}$}.

As explained previously, we assume that the gas accretion is only supported by the disc. At a given time n, the disc is defined by its total mass M𝒟n\hbox{$M_{\Disc}^n$} and its exponential radius r𝒟n\hbox{$r_{\Disc}^n$}. During a time step Δt = tn + 1tn, we assume that the new accreted mass acquires angular momentum and therefore forms a very thin disc structure with a mass (ff + coolt and an exponential radius deduced from the dark-matter halo virial radius rvir using the well-know spin parameter λ, racc=λ2rvir\hbox{$r_{\rm acc} = \dfrac{\lambda}{\sqrt{2}}r_{\rm vir}$} (Blumenthal et al. 1986; Sellwood & McGaugh 2005; Macciò et al. 2008; Cole et al. 2008; Antonuccio-Delogu et al. 2010; Muñoz-Cuartas et al. 2011).

The assembly of the pre-existing disc with the accreted mass leads to a secular evolution of the disc exponential radius that follows r𝒟2=M𝒟n+1M𝒟nwindΔt􏽺0􏽽􏽼0􏽻ejectedmassr𝒟n􏽼0􏽻􏽺0􏽽preexistingdisc+(ff+cool)Δtλ2rvir􏽼0􏽻􏽺0􏽽accretion·\begin{eqnarray} r_{\Disc}^2 = \dfrac{M_{\Disc}^{n+1}}{\underbrace{\dfrac{M_{\Disc}^{n} - \overbrace{\dot{M}_{\rm wind}\Delta t}^{\rm ejected~mass}}{r_{\Disc}^n}}_{\rm pre-existing~disc} + \underbrace{\dfrac{(\dot{M}_{\rm ff} + \dot{M}_{\rm cool})\Delta t}{\dfrac{\lambda}{\sqrt{2}}r_{\rm vir}}}_{\rm accretion}} \cdot \label{rd_evolution} \end{eqnarray}(58)In the model, the disc is defined as an object that contains two different parts: 𝒟={\begin{eqnarray} \Disc = \left\{ \begin{array}{ll} M_{\Disc,\sfg}{:}& \mbox{A star-forming (cold) gas component} \\ & \\ M_{\Disc,\star}{:}& \mbox{A stellar population.} \end{array}\right. \label{disc_object} \end{eqnarray}(59)The radial profile of the disc circular velocity is given by (Freeman 1970) Vcirc,𝒟2(r)=GM𝒟2r𝒟3r2[I0K0I1K1]r2r𝒟,\begin{eqnarray} V_{\rm circ,\Disc}^2(r) = \dfrac{G M_{\Disc}}{2r_{\Disc}^3}r^2\left[I_0K_0-I_1K_1\right]_{\frac{r}{2r_{\Disc}}} , \label{Vcirc_d} \end{eqnarray}(60)where Ir and Kr are modified Bessel functions of rank r.

5.3.1. Disc velocity dispersion and dynamical time

In this section, we use a new method based on energy conservation that is used to monitor the average velocity dispersion of the gas in the disc component. At a given time n, we assume that the energy dispersion of the gas \hbox{$\sigma_{v,\Disc}$} is given by Eσvn=12Mg,𝒟σv,𝒟2,\begin{eqnarray} E_{\sigma_v}^n = \frac{1}{2}M_{\sfg,\Disc}\sigma_{v,\Disc}^2, \label{E_vdisp} \end{eqnarray}(61)where \hbox{$M_{\sfg,\Disc}$} is the total mass of the gas in the disc and \hbox{$\sigma_{v,\Disc}$} the average velocity dispersion in the disc component. At time n the kinetic energy contained in the disc rotation is EVn=M𝒟n2r𝒟20rmax=11r𝒟exp(rr𝒟)Vcirc,𝒢2(r)rdr,\begin{eqnarray} E_{V}^n =\dfrac{M_{\Disc}^n}{2r_{\Disc}^2}\int_0^{r_max = 11r_{\Disc}}\exp\left(-\frac{r}{r_{\Disc}}\right)V_{\rm circ,\Gal}^2(r)r{\rm d}r , \label{E_Vrot} \end{eqnarray}(62)where \hbox{$V_{\rm circ,\Gal}(r)$} is the circular velocity of the galaxy (Eq. (55)). Here the integral is computed from 0 to \hbox{$11~r_{\Disc}$}, which encloses 99.99% of the disc mass.

As proposed by Khochfar & Silk (2009) or Ocvirk et al. (2008), we assume that the velocity dispersion is generated by gas infall and SNe energy injection. Since the infall mass is structured through a very thin disc in our model, assuming that the dynamics of the gas is given by the total potential well gives the following definition for the kinetic energy transported by gas fuelling Einf=(cool+ff)Δt011r𝒟exp(rr𝒟)Vcirc,𝒢2(r)rdr,\begin{eqnarray} E_{\rm inf} = \left(\dot{M}_{\rm cool}+\dot{M}_{\rm ff}\right)\Delta t \int_0^{11r_{\Disc}}exp\left(-\frac{r}{r_{\Disc}}\right)V_{\rm circ,\Gal}^2(r)r{\rm d}r , \label{infall_energy} \end{eqnarray}(63)where cool and ff are the cooling rate and the free-fall rate fuelling the galaxy, respectively.

thumbnail Fig. 9

Distribution of the \hbox{$\sigma_{r}/V_{\rm circ,\Gal}(r=r_{\Gal,50})$} ratio for 0 <z< 6 (coloured lines). The mean ratio decreases from 0.6 for high-z structures to 0.1 for local structures. The unstable (high sigma disc) population formed at high z disappears progressively and stable disc structures appear.

In addition we take the SNe kinetic energy injection in the gas into account, Esn,σv=(1ε𝒟)ηsn(1fKin,sn)Esn,𝒟Δt,\begin{eqnarray} E_{{\rm sn},\sigma_v} = (1 - \varepsilon_{\Disc})\eta_{\rm sn}(1-f_{\rm Kin,sn})E_{\rm sn}\dot{M}_{\star,\Disc}\Delta t , \label{Esn_sv} \end{eqnarray}(64)where

  • \hbox{$(1 - \varepsilon_{\Disc})$}is the fraction of SNe kinetic energy used to feed disc velocity dispersion;

  • ηsnis the number of SNe produced per unit of stellar mass formed ([ηsn]=M-1\hbox{$[\eta_{\rm sn}] = \Msun^{-1}$}). This parameter is linked to the initial mass function IMF;

  • fKin,snEsnis the SNe energy converted in kinetic energy. We use Esn = 1044 Joules and fKin,sn = 0.3 (Kahn 1975; Aguirre et al. 2001);

  • \hbox{$\dot{M}_{\star,\Disc}$}is the star formation rate in the disc.

To compute the evolution of the average velocity dispersion, we assume that the variation of the rotational energy between times n and n + 1, (ΔEV=EVn+1EVn\hbox{$\Delta E_{V} = E_{V}^{n+1} - E_{V}^n$}) is supported by gas infall energy (Einf ≥ ΔEV) and that the variation of the energy dispersion (ΔEσv,=Eσvn+1Eσvn\hbox{$\Delta E_{\sigma_v,} = E_{\sigma_v}^{n+1} - E_{\sigma_v}^{n}$}) is then given by ΔEσv,𝒟=(1fdisp)[Esn,σv,𝒟+(EinfΔEV)],\begin{eqnarray} \Delta E_{\sigma_v,\Disc} = \left(1-f_{\rm disp}\right)\left[E_{{\rm sn},\sigma_v,\Disc} + \left(E_{\rm inf}-\Delta E_{V}\right)\right] , \label{dEsv} \end{eqnarray}(65)where fdisp = 0.05 is a dissipation factor. The mean velocity dispersion in the disc is therefore given by σv,𝒟=2(Eσvn+ΔEσv)Mg,𝒟·\begin{eqnarray} \sigma_{v,\Disc} = \sqrt{\dfrac{2\left(E_{\sigma_v}^n +\Delta E_{\sigma_v}\right)}{M_{\sfg,\Disc}}} \cdot \label{sigma_v} \end{eqnarray}(66)Figure 9 shows the evolution with redshift of the normalized distribution of \hbox{$\sigma_{r}/V_{\rm circ,\Gal}(r=r_{\Gal,50})$}. We have assumed an isotropic three-dimensional velocity dispersion, σr=σv,𝒟/3\hbox{$\sigma_r = \sigma_{v,\Disc} /\!\sqrt{3}$}. While high-redshift discs present an average ratio close to 0.6, local discs have ratios close to 0.1. The disturbed disc population disappears progressively with time and a more stable population is generated.

The disc dynamical time is set to the minimal time between the disc’s full rotational time and local velocity dispersion. We use tdyn,𝒟=MIN(2r𝒢,50σr-1;2πVcirc,𝒢-1(r=r𝒢,50)).\begin{eqnarray} t_{\rm dyn,\Disc} = {\rm MIN}\left(2r_{\Gal,50}\sigma_r^{-1};~2\pi V_{\rm circ,\Gal}^{-1}\left(r=r_{\Gal,50}\right)\right) . \label{t_disc} \end{eqnarray}(67)

Table 7

Clumps symbols and their definition.

5.4. Disc instability, bulge-like component formation

For more than one decade, observations (e.g. Cowie et al. 1995; van den Bergh 1996; Elmegreen & Elmegreen 2005; Genzel et al. 2008; Bournaud et al. 2008) and hydrodynamic simulations (e.g. Bournaud et al. 2007; Ceverino et al. 2010, 2012) have shown that there are of gas-rich, turbulent discs. These kinds of discs are unstable, and undergo gravitational fragmentation that forms giant clumps. Such clumps interact and migrate to the centre of the galaxy where they form a bulge-like component (e.g. Elmegreen 2009; Dekel et al. 2009b).

Based on these works, we assume in our model that mass overdensities and low velocity dispersion (Eq. (66)) lead to the formation and the migration of giant clumps. We present an analytic self-consistent model of disc instabilities here. All parameters used in this section are described in Table 7.

5.4.1. Origin of instabilities

Following Toomre (1963, 1964), a disc becomes unstable if its gas surface density (\hbox{$\Sigma_{\Disc,\sfg}$}) is high and its circular velocity and/or mean velocity dispersion is low. The disc stability is controlled by the following parameter Q=σrκπGΣ𝒟,g=r𝒟,902σrκ0.9GM𝒟,g,\begin{eqnarray} Q = \dfrac{\sigma_r\kappa}{\pi G \Sigma_{\Disc,\sfg}}=\dfrac{r_{\Disc,90}^2 \sigma_r\kappa }{0.9 G M_{\Disc,\sfg}} , \label{Toomre} \end{eqnarray}(68)where κ=1M𝒟,g011r𝒟[4Ω2(1+r2Ωdr)]1/2Σ𝒟,g(r)rdr,\begin{eqnarray} \kappa = \dfrac{1}{M_{\Disc,\sfg}}\int_0^{11r_{\Disc}}\left[4\Omega^2\left(1+\dfrac{r}{2\Omega}\dfrac{{\rm d}\Omega}{{\rm d}r}\right)\right]^{1/2}\Sigma_{\Disc,\sfg}(r)r{\rm d}r, \label{kappa} \end{eqnarray}(69)is the mass-weighted epicyclic frequency with Ω the angular velocity and \hbox{$\Sigma_{\Disc,\sfg}(r)$} the gas mass surface density.

In Eq. (68), σr is the mean radial velocity dispersion. As for the disc dynamical time computation (Eq. (67)), we assume an isotropic three-dimensional velocity dispersion, σr=σv,𝒟/3\hbox{$\sigma_r = \sigma_{v,\Disc} /\!\sqrt{3}$}. The gaseous component of the disc \hbox{$M_{\Disc,\sfg}$} is considered as stable or marginally stable if Q>Qcrit = 1.0. We use this stability criterium to evaluate the mass that will form clumps by gravitational instabilities.

5.4.2. The clump component, \hbox{$\Clumps$}

We assume that disc instabilities are linked to the formation of giant clumps. These giant clumps are formed into the disc and migrate to the centre of the disc to form a pseudo-bulge component (see Fig. (8)). As previously, all parameters used in the description of disc instabilities linked to clumps are listed in the dedicated Table 7.

In our model, the clumpy phase (\hbox{$\Clumps$}) is not considered as a new separated component but is considered as an inherent part of the disc. The amount of gas in clumps is therefore given by a fraction of the disc mass. As for the disc component, clumps contain a cold star-forming gas \hbox{$M_{\Clumps,\sfg}$} and a stellar population \hbox{$M_{\Clumps,\star}$}. We monitor the stellar and the gaseous mass in clumps by taking their formation (Sect. 5.4.3), transfer (Sect. 5.4.4), and disruption (Sect. 5.4.5) into account. Star formation in the disc is divided into a clumpy and a homogeneous component.

thumbnail Fig. 10

Evolution of the disc-stability criterion (Toomre parameter Q) and clump properties. In the upper panel is shown the evolution of the Toomre parameter with z (Eq. (68), grey line). On this timeline, vertical arrows indicate a merger event (with short and long arrows for minor and major mergers, respectively). For 11 <z< 1.5, the Toomre parameter is larger than the critical value Q>Qcrit = 1, and therefore the disc is stable, but for 1.5 <z< 0.25 (highlighted by the two vertical black lines), the disc becomes unstable. During the unstable transition, bounded by two merger events, many clumps are formed and migrate. The lower panel shows a zoom of this unstable period. The blue and green lines show the number of clumps (left y-axis) and the average individual mass of these clumps (right y-axis), respectively. After a merger event, disc instabilities lead to the formation of a large number of clumps (42) with small individual masses (\hbox{$M_{\rm \Clumps,ind} < 10^5~\Msun$}). Then, during a series of merger events, clumps are destroyed and formed again, but with higher individual masses. Between two merger events, the number of clumps increases while the individual mass \hbox{$M_{\rm \Clumps,iind}$} remains stable. A decrease in one or two clumps is linked to the migration process. It is interesting to note that the migration process creates a temporary stable disc (Q>Qcrit = 1). In the second part of the unstable period, only one clump is formed but with a much higher mass. Finally, this single clump is disrupted and the disc becomes stable again.

5.4.3. Clump formation

We consider that disc instabilities (and therefore clumps formation) are only due to an excess of mass. When the homogeneous (non-clumpy) gas fraction of the disc (\hbox{$M_{\Disc,\sfg} - M_{\Clumps,\sfg})$} is unstable (Q<Qcrit), we divide the gaseous mass of the disc in two parts, a stable M𝒟,gs\hbox{$M_{\Disc,\sfg}^{\rm s}$} and an unstable part M𝒟,gu\hbox{$M_{\Disc,\sfg}^{\rm u}$}. We have 1Q=0.9G(M𝒟,gs+M𝒟,gu)r𝒟,902σrκ=1Qcrit+0.9GM𝒟,gur𝒟,902σrκ,\begin{eqnarray} \label{unstable_mass} \dfrac{1}{Q'} & = & \dfrac{0.9 G \left(M_{\Disc,\sfg}^{\rm s} + M_{\Disc,\sfg}^{\rm u}\right)}{r_{\Disc,90}^2\sigma_r\kappa}\nonumber\\ & = &\dfrac{1}{Q_{\rm crit}} + \dfrac{0.9 G M_{\Disc,\sfg}^{\rm u}}{r_{\Disc,90}^2\sigma_r\kappa}, \end{eqnarray}(70)with Q is a reduced-Toomre parameter that takes into account only the fraction of the disc homogeneous gas mass (the gaseous mass of the clumps being subtracted in the gas surface density computation). Within this framework, at a given time, a new amount of unstable gas feeds the clumpy part, \hbox{$\Clumps$} (see Fig. 8). The unstable gas mass is given by M𝒟,gu=r𝒟,902σrκ0.9G(1Q1Qcrit)·\begin{eqnarray} M_{\Disc,\sfg}^{\rm u} = \dfrac{r_{\Disc,90}^2\sigma_r\kappa}{0.9G}\left(\dfrac{1}{Q'} - \dfrac{1}{Q_{\rm crit}}\right)\cdot \label{unstable_mass_2} \end{eqnarray}(71)Knowing the new unstable mass in the disc, we can compute the instantaneous formation rate of the clump component (\hbox{$\dot{M}_{\rm \Clumps,in}$}, see Fig. 8) following 𝒞,in=M𝒟,gutdyn,𝒟·\begin{eqnarray} \dot{M}_{\rm \Clumps,in} = \dfrac{M_{\Disc,\sfg}^{\rm u}}{t_{\rm dyn,\Disc}}\cdot \label{clumps_in} \end{eqnarray}(72)Clumps can form stars. The dynamical time for star formation in clumps tdyn,𝒞=r𝒞3GfgM𝒞,indwithfg=M𝒞,gM𝒞,g+M𝒞,,\begin{eqnarray} t_{\rm dyn,\Clumps} = \sqrt{\dfrac{r_{\Clumps}^3}{Gf_{\sfg}M_{\rm \Clumps,ind}}}~~\mbox{with}~~f_{\sfg}=\dfrac{M_{\Clumps,\sfg}}{M_{\Clumps,\sfg}+M_{\Clumps,\star}} , \label{t_clumps} \end{eqnarray}(73)where \hbox{$r_{\Clumps}$} and \hbox{$M_{\rm \Clumps,ind}$} are the characteristic radius and mass of one clump, which we extract from Dekel et al. 2009b (their Eqs. (8) and (9)), r𝒞=π6δr𝒟,90,\begin{eqnarray} r_{\Clumps} = \dfrac{\pi}{6}\delta r_{\Disc,90}, \label{r_clumps} \end{eqnarray}(74)and M𝒞,ind=π236δM𝒟·\begin{eqnarray} M_{\rm \Clumps,ind} = \dfrac{\pi^2}{36}\delta M_{\Disc}\cdot \label{m_clumps} \end{eqnarray}(75)In these relations, δ is the disc mass fraction defined as δ=M𝒟(<r𝒟,90)Mvir(<r𝒟,90)+M𝒟(<r𝒟,90)+M(<r𝒟,90),\begin{eqnarray} \delta = \dfrac{M_{\Disc}({<}r_{\Disc,90})}{M_{\rm vir}({<}r_{\Disc,90}) + M_{\Disc}({<}r_{\Disc,90}) + M_{\Bulge}({<}r_{\Disc,90})} , \label{delta} \end{eqnarray}(76)where \hbox{$M_{\rm vir}({<}r_{\Disc,90})$}, \hbox{$M_{\Disc}({<}r_{\Disc,90})$}, and \hbox{$M_{\Bulge}({<}r_{\Disc,90})$} are the mass, enclosed in \hbox{$r_{\Disc,90}$} (radius for 90% of the disc mass), of the dark-matter halo, disc, and bulge component, respectively. As explained in the following section, the individual mass \hbox{$M_{\rm \Clumps,ind}$} is used as a threshold mass that allows us to compute the clump transfer to the pseudo-bulge component.

5.4.4. Clumps transfer

To take the migrations into account, we define \hbox{$M_{\Clumps '}$} as the mass of clumps that is transferred to the pseudo-bulge. The transfer rate, \hbox{$\dot{M}_{\rm \Clumps,migr}$}, writes as 𝒞,migr=M𝒞M𝒞tmigr,𝒞,\begin{eqnarray} \dot{M}_{\rm \Clumps,migr} = \dfrac{M_{\Clumps} -M_{\Clumps'}}{t_{\rm migr,\Clumps}} , \label{clumps_migr} \end{eqnarray}(77)where \hbox{$M_{\Clumps}$} is the total mass in clumps, and \hbox{$M_{\Clumps '}$} the mass of the clumps that is already in a transfer process. We add \hbox{$t_{\rm migr,\Clumps}$} as the typical clump transfer time scale. Like other clump properties, we compute this characteristic time following Dekel et al. (2009b; their Eq. (19)) tmigr,𝒞=2.1(Qδ)2tdyn,𝒟.\begin{eqnarray} t_{\rm migr,\Clumps} = 2.1\left(\dfrac{Q}{\delta}\right)^2t_{\rm dyn,\Disc}. \label{t_clumps_trans} \end{eqnarray}(78)When the transferred mass reservoir \hbox{$M_{\Clumps '}$} becomes larger than the reference mass \hbox{$M_{\rm \Clumps,ind}$}, the stellar and gas mass contained in this reservoir (\hbox{$\Clumps '$}) is subtracted from the disc and added to the pseudo-bulge component. The transfer is considered as a micro-merger event. When the micro-merger is finished, we reset the transferred mass, and recompute the characteristic individual mass of clumps \hbox{$M_{\rm \Clumps,ind}$}. This new characteristic mass will be the new threshold mass for the next clump migration process.

5.4.5. Clumps disruption

At a given time, a disc \hbox{$\Disc$} may have giant clumps (\hbox{$\Clumps$}). If, at the same time Q>Qcrit then we compute the disruption rate 𝒞,out=M𝒞tdisr,𝒞,\begin{eqnarray} \dot{M}_{\rm \Clumps,out} = \dfrac{M_{\Clumps}}{t_{\rm disr,\Clumps}} , \label{clumps_out} \end{eqnarray}(79)where t𝒞disr\hbox{$t_{\Clumps}^{\rm disr}$} is the disruption time computed as proposed by Dekel et al. (2009b; their Eq. (14)): tdisr,𝒞=1.4tdyn,𝒟Q·\begin{eqnarray} t_{\rm disr,\Clumps} = 1.4\dfrac{t_{\rm dyn,\Disc}}{Q}\cdot \label{t_clumps_disr} \end{eqnarray}(80)In Fig. 10 we show the evolution of the Toomre stability parameter for one disc chosen arbitrarily for z< 11. We also show the evolution of the number and individual mass of clumps for the unstable period. We see that for 0.8 <z< 1.5 the disc instability generates numerous clumps that are destroyed after each merger event and systematically formed again but with higher individual mass. During this period of instabilities, the decrease in the number of clumps has to be linked to the migration process. After z = 0.8 only one giant clumps can be formed, but with a higher mass. For z< 0.5 the disc is stable (Q>Qcrit = 1), the clumps disruption process is therefore turned on and the clumps mass decrease. After the last merger event, the disc does cannot form any clumps owing to the lack of gas.

5.4.6. Summary of mass transfers

The evolution of the mass transfers of the different components of the disc is governed by the following equations: 𝒟={\begin{eqnarray} \Disc = \left\{ \begin{array}{ll} \dot{M}_{\Disc,\sfg} = & + \dot{M}_{\rm ff} + \dot{M}_{\rm cool} + \dot{M}_{\rm wind,\star} \\ & - \dot{M}_{\star,\Disc} - \dot{M}_{\rm wind,sn,\Disc} - \left[M_{\Clumps ',\sfg}\right]_{M_{\Clumps'} > M_{\rm \Clumps,ind}} \\ \dot{M}_{\Disc,\star} = & + \dot{M}_{\star,\Disc} - \dot{M}_{\rm wind,\star} - \left[M_{\Clumps ',\star}\right]_{M_{\Clumps'} > M_{\rm \Clumps,ind}}. \end{array}\right. \end{eqnarray}(81)

thumbnail Fig. 11

Exchanges on the bulge scale. The bulge is considered as a passive component not fed by fresh accreted gas. A bulge () is formed during a merger or by clumps migrations from the disc (\hbox{$\Disc$}). An SMBH can evolve in the bulge. It is fed by a gas reservoir (the torus, \hbox{$\Torus$}) with an accretion rate acc. The torus is fed during a merger or by the clump migration process. The activity of the SMBH generates ejecta (•,out).

Table 8

Bulge symbols and their definition.

5.5. Pseudo-bulge or spheroidal galaxy

In our model, we assume the same properties for a post-major merger galaxy (with spheroidal morphology) and for a pseudo-bulge formed by giant-clumps migration. The mass transfers of these components are illustrated in Fig. 11, and the parameters are listed in Table 8.

As for the disc, we treat the (pseudo-)bulge as a numerical object composed of a cold gas reservoir (Mℬ,g) and a stellar population (Mℬ,). In a pseudo bulge the cold gas reservoir is fed by clump migration. This gas can be used to form new stars. In a spheroidal galaxy, formed during a merger event, the gas reservoir is only fed by the stellar cycle (see Fig. 13). In this context, the amount of gas is very low and therefore no star formation is possible. At each merger, the gas contained in the bulge component is added to the disc component of the remnant galaxy.

The mass distribution of the spheroidal galaxy or the pseudo-bulge is described by a Hernquist (1990) model, M(r)=Mr2(r+r)2,\begin{eqnarray} M_{\Bulge}(r) = M_{\Bulge} \dfrac{r^2}{(r+r_{\Bulge})^2} , \label{bulge_profile} \end{eqnarray}(82)where r is the characteristic radius of the bulge.

5.5.1. Bulge circular and escape velocity

The circular velocity of the bulge or pseudo-bulge is given by (Hernquist 1990) Vcirc,2(r)=GM,(r)r(r+r)2,\begin{eqnarray} V_{\rm circ,\Bulge}^2(r) = \dfrac{G M_{\Bulge,\BH}(r) r}{(r_{\Bulge}+r)^2} , \label{Vcirc_b} \end{eqnarray}(83)with M,(r)=Mr2(r+r)2+M,\begin{eqnarray} M_{\Bulge,\BH}(r) = M_{\Bulge}\dfrac{r^2}{(r+r_{\Bulge})^2} + M_{\BH}, \label{Vesc_b1} \end{eqnarray}(84)where M is the total mass (gas + stars) of the bulge, M the mass of the SMBH, and r the characteristic radius. We give also here the escape velocity of the bulge, Vesc,2(r)=2GM,(r+r)·\begin{eqnarray} V_{\rm esc,\Bulge}^2(r) = \dfrac{2G M_{\Bulge,\BH}}{(r_{\Bulge}+r)}\cdot \label{Vesc_b} \end{eqnarray}(85)

5.5.2. Bulge dynamical time and velocity dispersion

As commonly used in SAMs, we define the dynamical time of the bulge using its velocity dispersion (σv,) instead of its circular velocity, tdyn,=2r,50σv,,\begin{eqnarray} t_{\rm dyn,\Bulge} = \dfrac{2r_{\Bulge,50}}{\sigma_{v,\Bulge}} , \label{t_bulge} \end{eqnarray}(86)where rℬ,50 is the radius that enclosed 50% of the bulge mass, and where we used the definition of the velocity dispersion proposed by Hernquist (1990; their Eq. (10)), computed at the half-mass radius.

5.5.3. Summary of the mass transfers

The physical processes acting in the bulges are linked to clump migration and gas ejection (wind,agn due to SMBH activity). As in the disc, the gas in the bulge is converted into stars and these stars feed the gas reservoir (wind,, Sect. 7 and Fig. 13). We summarize the mass transfer rates with the following equations: ={\begin{eqnarray} \Bulge = \left\{ \begin{array}{ll} \dot{M}_{\Bulge,\sfg} = & \left[M_{\Clumps ',\sfg}\right]_{M_{\Clumps '} > M_{\rm \Clumps,ind}} + \dot{M}_{\rm wind,\star}\\ & - \dot{M}_{\star,\Bulge} - \dot{M}_{\rm wind,sn,\Bulge} - \dot{M}_{\rm wind,\BH}\\ \dot{M}_{\Bulge,\star} = & \left[M_{\Clumps ',\star}\right]_{M_{\Clumps '} > M_{\rm \Clumps,ind}} + \dot{M}_{\star,\Bulge} - \dot{M}_{\rm wind,\star}. \end{array}\right. \end{eqnarray}(87)

Table 9

SMBH symbols and their definition.

6. The SMBH component

An SMBH can be formed and evolve in the centre of the galaxy. This section describes the hypothesis behind the formation of an SMBH and the accretion mechanisms that govern its evolution. The parameters used in this section are listed in a dedicated table (Table 9).

We decided to form the SMBH by converting a fraction of the bulge mass into an SMBH when the bulge mass (M) becomes greater than a threshold Mℬ,• = 106 M. We assume that the formation process of an SMBH is only possible during a merger. Indeed, we assume that the SMBH formation needs powerful processes. Therefore, after a merger, if the total mass of the bulge is higher than Mℬ,• an SMBH with a mass M•,min = 103Msun is created. We subtract the equivalent mass in their respective proportions from the gas phase and from the stellar population. As can be expected, the hypothesis of a constant seed mass as a threshold to create an SMBH generates a scatter in the M versus M relation. This scatter is similar to what is obtained by Croton et al. (2006).

The SMBH evolves in the centre of the bulge and grows by accretion processes. The total accretion rate is given by ,inf=MIN(MAX(𝒯,Bondi),Edd),\begin{eqnarray} \dot{M}_{\rm \BH,inf} = {\rm MIN}\left({\rm MAX}\left(\dot{M}_{\Torus},\dot{M}_{\rm Bondi}\right),\dot{M}_{\rm Edd}\right) , \label{agn_inf} \end{eqnarray}(88)where \hbox{$\dot{M}_{\Torus}$} is an accretion rate linked to the evolution of a gaseous torus formed around the black hole and defined as 𝒯=M𝒯2tdyn,·\begin{eqnarray} \dot{M}_{\Torus}=\dfrac{M_{\Torus}}{2t_{\rm dyn,\Bulge}}\cdot \label{agn_inf_torus} \end{eqnarray}(89)In this last equation, \hbox{$M_{\Torus}$} is the total mass of gas in the torus and tdyn,ℬ the dynamical time of the bulge. During the galaxy evolution processes the torus \hbox{$M_{\Torus}$} is fed by

  • the gas coming from clumps migration (\hbox{$M_{\Clumps ',\sfg}$}),

  • a fraction of the bulge gas (\hbox{$\delta_{\Torus} = 10^{-4}$})5, instantaneously added during a merger event.

In Eq. (88), Bondi is the classical Bondi (1952) accretion rate. We use the same calculation as proposed by Somerville et al. (2008), Bondi=3πGμmp4kBTΛ(T,Z)M,\begin{eqnarray} \dot{M}_{\rm Bondi} = \dfrac{3\pi G \mu m_{\rm p}}{4}\dfrac{k_{\rm B} T_{\BH}}{\Lambda(T_{\BH},Z)}M_{\BH}, \label{Bondi} \end{eqnarray}(90)where T ≃ 107 K is the mean temperature in the accretion torus, Λ the cooling function, and M the SMBH mass. At any time, the accretion onto the SMBH is possible only if the torus contains gas (\hbox{$M_{\Torus} > 0$}) and only in the Eddington accretion rate limit, Edd=4πGmp[1(11+η)(1+fKin,)]σTcM.\begin{eqnarray} \dot{M}_{\rm Edd} = \dfrac{4\pi G m_{\rm p}}{\left[1-\left(\frac{1}{1+\eta_{\BH}}\right)\left(1+f_{\rm Kin,\BH}\right)\right]\sigma_T{\rm c}}M_{\BH}. \label{Edd_acc_rate} \end{eqnarray}(91)We recall that the complete description of ejecta processes driven by the SMBH is given in Sects. 4.1.2 and 4.2.3.

thumbnail Fig. 12

Fraction of structures hosting an SMBH in the centre of their galaxy as a function of the dark-matter halo mass for a set of redshift (coloured lines). The fraction increases with the dark-matter halo mass to reach a value of 1 for haloes more massive than ≃ 4 × 1011 M. The insert panel shows the fraction of active galaxy nuclei as a function of the dark-matter halo mass. An SMBH is considered to be active if its bolometric luminosity is higher than 5% of its Eddington luminosity.

The instantaneous bolometric luminosity produced by the SMBH is deduced from Eq. (31) and is set to L,bol=(1fTh,).\begin{eqnarray} L_{\rm \BH,bol} = \left(1-f_{\rm Th,\BH}\right)\dot{Q}_{\BH}. \label{agn_Lum} \end{eqnarray}(92)We show in Fig. 12 the evolution of the fraction of dark-matter haloes that hosts an SMBH, as a function of the halo mass. The fraction increases with Mh. In the model presented here, all structures with dark-matter halo mass larger than 4 × 1011 M host an SMBH, regardless of the redshift.

We also show the fraction of active SMBH in Fig. 12. An SMBH is considered as active if L•,bol> 0.05LEdd. For redshift z> 2 this fraction is large then >40% but it decreases for structures evolving at lower redshift (z1). This decrease is mainly visible at the highest and lowest mass. Haloes with intermediate mass (1011<Mvir< 1012 M) host most of the active SMBH.

7. The stellar component

In the disc, bulge, and clumps, the gas (\hbox{$M_{\sfg,\Disc}$}, Mg,ℬ, \hbox{$M_{\sfg,\Clumps}$}, respectively) is converted into stars (\hbox{$M_{\star,\Disc}$}, M,ℬ, \hbox{$M_{\star,\Clumps}$}, respectively). We adopt the same definition for the star formation rate in all these components, ,i=ε,iMg,itdyn,i·\begin{eqnarray} \dot{M}_{\star,i} = \varepsilon_{\star,i}\dfrac{M_{\sfg,i}}{t_{{\rm dyn},i}}\cdot \label{sfr} \end{eqnarray}(93)In this formulation, ε,i is a free efficiency parameter, depending on the component. We adopt ε,i = 0.02 in discs and bulge and 0.2 in clumps. These values allow to be reproduced the Kennicutt (1998) law on both galactic and clump scales. Here Mg,i is the total mass of star-forming, gas and tdyn,i is the dynamical time (Eqs. (67), (86), and (73) for the disc, bulge and, clumps, respectively).

thumbnail Fig. 13

Cartoon showing the stellar cycle. Stars are formed from cold gas in the galaxy components (). Stars evolve and die. We follow the evolution of the ejected mass produced by the stellar population wind,.

If the disc contains some clumps, the star formation is divided into two parts: i) the homogeneous gas part in the disc (\hbox{$M_{\sfg,\Disc}$}, Sect. 5.4.3 and Eq. (70)) and ii) the clumpy part (\hbox{$M_{\sfg,\Clumps}$}). The overall disc star formation rate is computed by considering the fraction of mass in each of the two components ,𝒟+,𝒞=ε,𝒟Mg,𝒟tdyn,𝒟+ε,𝒞Mg,𝒞tdyn,𝒞·\begin{eqnarray} \dot{M}_{\star,\Disc} + \dot{M}_{\star,\Clumps} = \varepsilon_{\star,\Disc}\dfrac{{M}_{\sfg,\Disc}}{t_{\rm dyn,\Disc}} + \varepsilon_{\star,\Clumps}\dfrac{{M}_{\sfg,\Clumps}}{t_{\rm dyn,\Clumps}}\cdot \label{sfr_disc_all} \end{eqnarray}(94)As shown in Fig. 13 the evolution of stellar population is closely related to the gas reservoirs. At the end of their life, stars re-injects gas into the galaxy’s interstellar medium. This gas is enriched in metals. As in its predecessor (Hatton et al. 2003), the metals enrichment is followed using the Devriendt et al. (1999) prescription.

8. Mergers

8.1. Major or minor merger?

We define the merger type (major or minor) using the following mass ratio, ηmerger=MIN(M1/2,1;M1/2,2)MAX(M1/2,1;M1/2,2),\begin{eqnarray} \eta_{\rm merger} = \dfrac{{\rm MIN}(M_{1/2,1};~M_{1/2,2})}{{\rm MAX}(M_{1/2,1};~M_{1/2,2})}, \label{eta_merger} \end{eqnarray}(95)where \hbox{$M_{1/2} = M_{\Gal}({<}r_{\Gal,50}) + 2M_{\rm vir}({<}r_{\Gal,50}) $} is the total mass (galaxy and dark-matter halo), enclosed in the galaxy half-mass radius \hbox{$r_{\Gal,50}$}. We consider a minor merger event when ηmerger< 0.25. Table 10 gives symbols and their definition associated with mergers.

Table 10

Merger event symbols and their definition.

8.2. Minor mergers

During a minor merger, we keep the disc and the bulge. The stellar populations of each disc (bulge), are added to form the remnant disc (bulge). In contrast, cold gas reservoirs are added in the remnant disc. The size of the remnant disc (bulge) is set to the larger disc (bulge) progenitor size. The velocities (dispersion and circular) are recomputed with the properties of the remnant galaxy and its dark-matter halo.

We consider that the minor merger of two discs has a strong impact on the evolution and formation of clumps. We decided to destroy all clumps in the remnant disc. New clumps have to be formed in the new disc after the merger, by taking the post-merger dark-matter and disc properties into account.

8.3. Major mergers

In the case of a major merger (ηmerger> 0.25), the galaxy structure and dynamic are completely modified. The stellar population of the remnant galaxy presents an elliptical form modelled by the same Hernquist (1990) mass distribution than the bulge component. In parallel the gas components of the two progenitors are added and a new disc is formed. The characteristics of this new disc are computed in the context of the post-merger dark-matter halo (rvir).

At the merger time, the half-mass bulge radius rℬ,50 is computed using (Hatton et al. 2003), r,50=(M1/2,1+M1/2,2)M1/2,12r𝒢,1+M1/2,22r𝒢,22.5M1/2,1M1/2,2(r𝒢,1+r𝒢,2),\begin{eqnarray} r_{\Bulge,50} = \dfrac{\left(M_{1/2,1} + M_{1/2,2}\right)}{\frac{M_{1/2,1}^2}{r_{\Gal,1}} + \frac{M_{1/2,2}^2}{r_{\Gal,2}} - 2.5\frac{M_{1/2,1} M_{1/2,2}}{\left(r_{\Gal,1}+r_{\Gal,2}\right)}}, \label{merger_bulge_half_mass_radius} \end{eqnarray}(96)where M1 / 2,i has been introduced in Eq. (95). The characteristic bulge radius (r) is then deduced following r,50=r(1+2).\begin{eqnarray} r_{\Bulge,50} = r_{\Bulge}\left(1 + \sqrt{2}\right). \label{merger_bulge_radius} \end{eqnarray}(97)In the case of a micro-merger event, when a clump migrates to the centre, we use the same framework as in Eq. (95), M1 / 2,1 being in this case the bulge plus dark-matter halo mass enclosed in the bulge half mass radius, and \hbox{$M_{1/2,2}\simeq\Clumps_{\rm ind}$} the transferred clump mass.

8.4. SMBH: evolution and/or formation

For the SMBH component, if the two progenitors contain an SMBH, their masses are added, and a fraction \hbox{$\delta_{\Torus} = 10^{-4}$} of the gas contained in the remnant galaxy is added to the SMBH accretion torus (\hbox{$M_{\Torus}$}). If only one of the two bulge components contains an SMBH, the remnant galaxy keeps this SMBH component. If the progenitors do not contain any SMBH, but if the remnant bulge mass becomes higher than the mass threshold, an SMBH is created with a mass M•,min = 103 M.

8.5. Merger boost factor

Mergers are violent events. During a short time-scale galaxy evolution processes are strongly modified and secular evolution laws (efficiencies) are no longer valid. To take these modifications into account we define a merger boost factor following εboost(τ)=100ηmergerηgasexp(ττmerger)·\begin{eqnarray} \varepsilon_{\rm boost}(\tau) = 100\eta_{\rm merger}\eta_{\rm gas}\exp\left(-\dfrac{\tau}{\tau_{\rm merger}}\right)\cdot \label{boost_factor} \end{eqnarray}(98)In this definition, ηmerger is the progenitor mass ratio (Eq. (95)), ηgas the gas fraction in the post-merger structure, τ the time elapsed since the last merger event, and τmerger = 0.05 Gyr is the characteristic merger time scale. This factor is used to modify the efficiencies just after a merger event. We use MAX[ 1,εboost(τ) ] as a correction factor to the star formation law (Eq. (93)), for each component of the galaxy (disc and bulge) and for the SMBH torus accretion rate (Eq. (89)). The formulation used here allows simulating a time-dependent gas compression (decreasing with time), and takes the gas content of the two progenitors into account. The more gas they contain, the higher the gas compression.

8.6. The hot gas phase

During a merger, the hot gas masses of the two haloes are added, and the equilibrium temperature is computed with the thermal energy of each hot atmosphere, Ethem,i=32kBTiMhot,iμmp,\begin{eqnarray} E_{{\rm them},i} = \dfrac{3}{2}k_{\rm B}\overline{T}_i\dfrac{M_{{\rm hot},i}}{\mu m_{\rm p}}, \label{Etherm_merger} \end{eqnarray}(99)which is added to the gravitational interaction energy, Egrav=GMhot,1Mhot,2rhot,1+rhot,2,\begin{eqnarray} E_{\rm grav} = G\dfrac{M_{\rm hot,1}M_{\rm hot,2}}{r_{\rm hot,1}+r_{\rm hot,2}}, \label{Egrav_merger} \end{eqnarray}(100)where rhot,i is the half-mass radius. The equilibrium temperature is then given by T=2μmp3kBEthem,1+Ethem,2+EgravMhot,1+Mhot,2·\begin{eqnarray} \overline{T}=\dfrac{2\mu m_{\rm p}}{3k_{\rm B}}\dfrac{E_{\rm them,1} + E_{\rm them,2} + E_{\rm grav}}{M_{\rm hot,1}+M_{\rm hot,2}}\cdot \label{Thot_merger} \end{eqnarray}(101)

9. Conclusions

We have presented in this paper the physical mechanisms acting on the formation and evolution of a galaxy in a dark-matter halo in a semi-analytical framework. The model follows the dark-matter structuration through merger trees. The growth of dark-matter haloes was followed using a set of properties, including the virial mass and radius and density profile parameters. The dark-matter accretion rate is computed using background particles, which were identified as particles that have never been assigned to any other halo at a previous time step. These particles allowed following the baryonic accretion. Using a photoionisation model, the baryonic accretion rate was deduced from the smooth dark-matter accretion rate, built with the background particles. The baryonic accretion is then divided into two parts:

  • A cold mode, that dominates the growth of small structures. It relies on the fact that the cooling time of the gas accreted by low-mass haloes (Mvir< 5 × 1010 M) is shorter than their dynamical time.

  • A hot mode: for massive haloes (Mvir> 5 × 1010 M), the accreted hot gas has to cool before collapsing.

In the model, the evolution of the hot atmosphere, in equilibrium in the dark-matter potential well, is based on a detailed monitoring of its average temperature and density profile. These two quantities were followed by using energy conservation laws based on SN feedback and hot-gas accretion. In addition, the coupling between the pre-existing hot-gas phase and the wind phase generated by the galaxy, allowed us to compute the gas escape fraction. The cooling rate, computed using the standard White & Frenk (1991) prescription, relies on this explicit monitoring of the hot-gas phase properties.

Coming from the cooling process or directly from the cold accretion mode, the cold gas forms a disc at the centre of the dark-matter halo. The gas in the disc can form stars. In the fiducial approach, we used the Kennicutt law based on the total gas mass and the dynamical time of the disc.

A feature of this new model is the treatment of instabilities in the disc. Based on Dekel et al. (2009b), the model we developed accounts for the formation, migration, and disruption of giant clumps. The dynamic of the disc is followed using the rotational profile that is deduced from the entire mass distribution (dark-matter halo, disc, and bulge). To this rotational profile is added a model to follow the average gas velocity dispersion. The evolution through cosmic time of the ratio \hbox{$\sigma_{r}/V_{\rm cric,\Gal}$} indicates that the population of highly disturbed discs (\hbox{$\sigma_{r}/V_{\rm cric,\Gal} \simeq 0.6$}), formed at high z, decreases progressively with time, and gives birth to a more quiet disc population (\hbox{$\sigma_{r}/V_{\rm cric,\Gal} \le 0.1$}). The mean velocity dispersion, added to the gas surface density computation, allows the stability behaviour of the disc to be derived. Using a method based on the Toomre parameter, we compute at each time-step the unstable mass, which is stored in clumps. The clumpy mass then follows a migration process. The mass that is transferred participates to the formation of a pseudo-bulge component.

This pseudo-bulge, or bulge formed during mergers, can host the formation and the evolution of an SMBH. The accretion of mass onto the SMBH generates energy that is at the origin of mass ejection through wind. Together with the mass ejected due to SNe, the SMBH activity participates in the heating of the hot atmosphere around the galaxy. SN and SMBH are the two main process that allow the star formation to be regulated. Indeed gas ejection induces gas depletion in the galaxy and the heating of the hot atmosphere reduces the cooling rate.

In our model all the up-to-date mechanisms used to regulate the star formation process were considered. Some variations of the fiducial recipes are discussed in detail in Cousin et al. (2015). We showed in this paper that even using strong efficiencies for the photoionisation and feedback processes, the star formation activity remains too efficient in low-mass galaxies. We thus introduce an ad hoc modification of the standard paradigm, based on the presence of a no-star-forming gas component, and an artificial concentration of the star-forming gas in galaxy discs. The no-star-forming gas reservoir generates a delay between the accretion of the gas and star formation. This model agrees much better with the observations of the stellar mass function in the low-mass range than the reference model, and agrees quite well with a large set of observations, including the redshift evolution of the specific star formation rate.

SAMs are currently the best models for following a large number of galaxies in their dark-matter halo. Their predictions, deduced from large scale dark-matter simulations, allow interpreting the recent surveys of galaxy population. However, even if the decoupling between dark-matter evolution (merger tree) and the baryonic physics allows exploring a large set of prescriptions in a short computation time, SAMs cannot follow all the complexity of the galaxy evolution. This applies particularly to the computation of the star formation rate. SAMs generally assume that the gas distribution follows an exponential profile, in a infinitely thin disc. All the gas content in the disc is available to form stars. In contrast, real galaxies host a complex gas structuration; high-z structures are generally gas-rich discs and show a finite scale-height. Galaxy discs are a multi-fluid system, in which multi-scale interactions (2D, 3D Shi et al. 2011) generate instabilities that participate in the ISM structuration (e.g. Jog & Solomon 1984; Hoffmann & Romeo 2012). A better description of these ISM structuration mechanisms is essential for better understanding the star formation regulation processes. In a future work, we plan to follow the gas structuration in the disc from the large (λh) to the small scales (λh), assuming a finite scale-height (h) and using the turbulent cascade description. Assuming that star formation occurs only on the smallest scales, we could estimate the star-forming gas fraction and give a better description of the star formation activity.


1

We use here a fit of the relation given in Fig. 15 of Okamoto et al. (2008).

2

The number of SNe increases with the number of massive stars in the initial mass function (IMF). Thus, top-heavy IMFs have a larger ηsn. In this context, the SNe feedback is also larger.

3

In this context, T\hbox{$\overline{T}$} used here is the result at the previous time step of the Eq. (42).

4

The most massive structure is defined as the structure with the higher mass of hot gas Mhot, which is not necessarily the highest Mvir dark-matter structure.

5

This value has been chosen to reproduce the trend of the M vs. M relation.

6

Over the fraction of the time step already spent.

Acknowledgments

We thank Nicole Nesvabda, Pierre Guillard, and Nicolas Bouché for useful discussions of the SNe and SMBH feedback. We also thank Jérémy Blaizot and Andrea Cattaneo for very insightful discussions.

References

  1. Aguirre, A., Hernquist, L., Schaye, J., et al. 2001, ApJ, 560, 599 [NASA ADS] [CrossRef] [Google Scholar]
  2. Antonuccio-Delogu, V., Dobrotka, A., Becciani, U., et al. 2010, MNRAS, 407, 1338 [NASA ADS] [CrossRef] [Google Scholar]
  3. Aubert, D., Pichon, C., & Colombi, S. 2004, MNRAS, 352, 376 [NASA ADS] [CrossRef] [Google Scholar]
  4. Babul, A., & Rees, M. J. 1992, MNRAS, 255, 346 [NASA ADS] [CrossRef] [Google Scholar]
  5. Baugh, C. M. 2006, Rep. Prog. Phys., 69, 3101 [NASA ADS] [CrossRef] [Google Scholar]
  6. Behroozi, P. S., Wechsler, R. H., & Conroy, C. 2013a, ApJ, 762, L31 [NASA ADS] [CrossRef] [Google Scholar]
  7. Behroozi, P. S., Wechsler, R. H., & Conroy, C. 2013b, ApJ, 770, 57 [NASA ADS] [CrossRef] [Google Scholar]
  8. Benson, A. J., & Bower, R. 2011, MNRAS, 410, 2653 [NASA ADS] [CrossRef] [Google Scholar]
  9. Benson, A. J., Lacey, C. G., Baugh, C. M., Cole, S., & Frenk, C. S. 2002, MNRAS, 333, 156 [NASA ADS] [CrossRef] [Google Scholar]
  10. Benson, A. J., Bower, R. G., Frenk, C. S., et al. 2003, ApJ, 599, 38 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bertone, S., Stoehr, F., & White, S. D. M. 2005, MNRAS, 359, 1201 [NASA ADS] [CrossRef] [Google Scholar]
  12. Béthermin, M., Doré, O., & Lagache, G. 2012, A&A, 537, L5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Blumenthal, G. R., Faber, S. M., Flores, R., & Primack, J. R. 1986, ApJ, 301, 27 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  14. Bondi, H. 1952, MNRAS, 112, 195 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  15. Bournaud, F., Elmegreen, B. G., & Elmegreen, D. M. 2007, ApJ, 670, 237 [Google Scholar]
  16. Bournaud, F., Daddi, E., Elmegreen, B. G., et al. 2008, A&A, 486, 741 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Bower, R. G., Benson, A. J., Malbon, R., et al. 2006, MNRAS, 370, 645 [NASA ADS] [CrossRef] [Google Scholar]
  18. Bower, R. G., Benson, A. J., & Crain, R. A. 2012, MNRAS, 422, 2816 [NASA ADS] [CrossRef] [Google Scholar]
  19. Bullock, J. S., Kravtsov, A. V., & Weinberg, D. H. 2000, ApJ, 539, 517 [NASA ADS] [CrossRef] [Google Scholar]
  20. Capelo, P. R., Coppi, P. S., & Natarajan, P. 2012, MNRAS, 422, 686 [NASA ADS] [CrossRef] [Google Scholar]
  21. Caputi, K. I., Cirasuolo, M., Dunlop, J. S., et al. 2011, MNRAS, 413, 162 [NASA ADS] [CrossRef] [Google Scholar]
  22. Cattaneo, A., Dekel, A., Devriendt, J., Guiderdoni, B., & Blaizot, J. 2006, MNRAS, 370, 1651 [NASA ADS] [CrossRef] [Google Scholar]
  23. Ceverino, D., Dekel, A., & Bournaud, F. 2010, MNRAS, 404, 2151 [NASA ADS] [Google Scholar]
  24. Ceverino, D., Dekel, A., Mandelker, N., et al. 2012, MNRAS, 420, 3490 [NASA ADS] [CrossRef] [Google Scholar]
  25. Cole, S. 1991, ApJ, 367, 45 [NASA ADS] [CrossRef] [Google Scholar]
  26. Cole, S., Aragon-Salamanca, A., Frenk, C. S., Navarro, J. F., & Zepf, S. E. 1994, MNRAS, 271, 781 [NASA ADS] [CrossRef] [Google Scholar]
  27. Cole, S., Lacey, C. G., Baugh, C. M., & Frenk, C. S. 2000, MNRAS, 319, 168 [NASA ADS] [CrossRef] [Google Scholar]
  28. Cole, S., Helly, J., Frenk, C. S., & Parkinson, H. 2008, MNRAS, 383, 546 [NASA ADS] [CrossRef] [Google Scholar]
  29. Couchman, H. M. P., & Rees, M. J. 1986, MNRAS, 221, 53 [NASA ADS] [CrossRef] [Google Scholar]
  30. Cousin, M., Lagache, G., Blaizot, I., Bethermin, M., & Guiderdoni, B., et al. 2015, A&A, 575, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Cowie, L. L., Hu, E. M., & Songaila, A. 1995, AJ, 110, 1576 [NASA ADS] [CrossRef] [Google Scholar]
  32. Croton, D. J., Springel, V., White, S. D. M., et al. 2006, MNRAS, 365, 11 [NASA ADS] [CrossRef] [Google Scholar]
  33. De Lucia, G., Boylan-Kolchin, M., Benson, A. J., Fontanot, F., & Monaco, P. 2010, MNRAS, 406, 1533 [NASA ADS] [Google Scholar]
  34. Dekel, A., & Silk, J. 1986, ApJ, 303, 39 [NASA ADS] [CrossRef] [Google Scholar]
  35. Dekel, A., Birnboim, Y., Engel, G., et al. 2009a, Nature, 457, 451 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  36. Dekel, A., Sari, R., & Ceverino, D. 2009b, ApJ, 703, 785 [Google Scholar]
  37. Devriendt, J. E. G., Guiderdoni, B., & Sadat, R. 1999, A&A, 350, 381 [NASA ADS] [Google Scholar]
  38. Doroshkevich, A. G., Zel’dovich, Y. B., & Novikov, I. D. 1967, Sov. Astron., 11, 233 [NASA ADS] [Google Scholar]
  39. Efstathiou, G. 1992, MNRAS, 256, 43P [Google Scholar]
  40. Efstathiou, G. 2000, MNRAS, 317, 697 [NASA ADS] [CrossRef] [Google Scholar]
  41. Eke, V. R., Frenk, C. S., Baugh, C. M., et al. 2004, MNRAS, 355, 769 [NASA ADS] [CrossRef] [Google Scholar]
  42. Eke, V. R., Baugh, C. M., Cole, S., et al. 2005, MNRAS, 362, 1233 [NASA ADS] [CrossRef] [Google Scholar]
  43. Elmegreen, B. G. 2009, in Galaxy Evolution: Emerging Insights and Future Challenges, eds. S. Jogee, I. Marinova, L. Hao, & G. A. Blanc, ASP Conf. Ser., 419, 23 [Google Scholar]
  44. Elmegreen, B. G., & Elmegreen, D. M. 2005, ApJ, 627, 632 [NASA ADS] [CrossRef] [Google Scholar]
  45. Emonts, B. H. C., Morganti, R., Tadhunter, C. N., et al. 2005, MNRAS, 362, 931 [NASA ADS] [CrossRef] [Google Scholar]
  46. Fakhouri, O., & Ma, C.-P. 2010, MNRAS, 401, 2245 [NASA ADS] [CrossRef] [Google Scholar]
  47. Fakhouri, O., Ma, C.-P., & Boylan-Kolchin, M. 2010, MNRAS, 406, 2267 [NASA ADS] [CrossRef] [Google Scholar]
  48. Faucher-Giguère, C.-A., Kereš, D., & Ma, C.-P. 2011, MNRAS, 417, 2982 [NASA ADS] [CrossRef] [Google Scholar]
  49. Freeman, K. C. 1970, ApJ, 160, 811 [NASA ADS] [CrossRef] [Google Scholar]
  50. Frye, B., Broadhurst, T., & Benítez, N. 2002, ApJ, 568, 558 [NASA ADS] [CrossRef] [Google Scholar]
  51. Genel, S., Bouché, N., Naab, T., Sternberg, A., & Genzel, R. 2010, ApJ, 719, 229 [NASA ADS] [CrossRef] [Google Scholar]
  52. Genzel, R., Burkert, A., Bouché, N., et al. 2008, ApJ, 687, 59 [NASA ADS] [CrossRef] [Google Scholar]
  53. Gnedin, N. Y. 2000, ApJ, 542, 535 [NASA ADS] [CrossRef] [Google Scholar]
  54. Guillard, P., Ogle, P. M., Emonts, B. H. C., et al. 2012, ApJ, 747, 95 [NASA ADS] [CrossRef] [Google Scholar]
  55. Guo, Q., White, S., Li, C., & Boylan-Kolchin, M. 2010, MNRAS, 404, 1111 [NASA ADS] [Google Scholar]
  56. Guo, Q., White, S., Boylan-Kolchin, M., et al. 2011, MNRAS, 413, 101 [NASA ADS] [CrossRef] [Google Scholar]
  57. Hatton, S., Devriendt, J. E. G., Ninin, S., et al. 2003, MNRAS, 343, 75 [NASA ADS] [CrossRef] [Google Scholar]
  58. Heckman, T. M., Lehnert, M. D., Strickland, D. K., & Armus, L. 2000, ApJS, 129, 493 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  59. Helly, J. C., Cole, S., Frenk, C. S., et al. 2003, MNRAS, 338, 903 [NASA ADS] [CrossRef] [Google Scholar]
  60. Henriques, B. M. B., White, S. D. M., Thomas, P. A., et al. 2013, MNRAS, 431, 3373 [NASA ADS] [CrossRef] [Google Scholar]
  61. Hernquist, L. 1990, ApJ, 356, 359 [NASA ADS] [CrossRef] [Google Scholar]
  62. Hoeft, M., Yepes, G., Gottlöber, S., & Springel, V. 2006, MNRAS, 371, 401 [NASA ADS] [CrossRef] [Google Scholar]
  63. Hoffmann, V., & Romeo, A. B. 2012, MNRAS, 425, 1511 [NASA ADS] [CrossRef] [Google Scholar]
  64. Ikeuchi, S. 1986, Ap&SS, 118, 509 [NASA ADS] [CrossRef] [Google Scholar]
  65. Ilbert, O., Salvato, M., Le Floc’h, E., et al. 2010, ApJ, 709, 644 [NASA ADS] [CrossRef] [Google Scholar]
  66. Ilbert, O., McCracken, H. J., Le Fevre, O., et al. 2013, A&A, 556, A55 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Jog, C. J., & Solomon, P. M. 1984, ApJ, 276, 114 [NASA ADS] [CrossRef] [Google Scholar]
  68. Kahn, F. D. 1975, in Int. Cosmic Ray Conf., 11, 3566 [Google Scholar]
  69. Kauffmann, G., White, S. D. M., & Guiderdoni, B. 1993, MNRAS, 264, 201 [NASA ADS] [CrossRef] [Google Scholar]
  70. Kauffmann, G., Colberg, J. M., Diaferio, A., & White, S. D. M. 1999, MNRAS, 303, 188 [NASA ADS] [CrossRef] [Google Scholar]
  71. Kennicutt, Jr., R. C. 1998, ApJ, 498, 541 [NASA ADS] [CrossRef] [Google Scholar]
  72. Kereš, D., Katz, N., Weinberg, D. H., & Davé, R. 2005, MNRAS, 363, 2 [NASA ADS] [CrossRef] [Google Scholar]
  73. Khochfar, S., & Silk, J. 2009, ApJ, 700, L21 [NASA ADS] [CrossRef] [Google Scholar]
  74. Komatsu, E., & Seljak, U. 2001, MNRAS, 327, 1353 [NASA ADS] [CrossRef] [Google Scholar]
  75. Koo, B.-C., Lee, J.-J., & Seward, F. D. 2002, AJ, 123, 1629 [Google Scholar]
  76. Kravtsov, A. V., Gnedin, O. Y., & Klypin, A. A. 2004, ApJ, 609, 482 [NASA ADS] [CrossRef] [Google Scholar]
  77. Lacey, C., & Cole, S. 1994, MNRAS, 271, 676 [Google Scholar]
  78. Larson, R. B. 1974, MNRAS, 169, 229 [NASA ADS] [CrossRef] [Google Scholar]
  79. Leauthaud, A., Tinker, J., Bundy, K., et al. 2012, ApJ, 744, 159 [NASA ADS] [CrossRef] [Google Scholar]
  80. Lehnert, M. D., Tasse, C., Nesvadba, N. P. H., Best, P. N., & van Driel, W. 2011, A&A, 532, L3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Lu, Y., Kereš, D., Katz, N., et al. 2011, MNRAS, 416, 660 [NASA ADS] [CrossRef] [Google Scholar]
  82. Macciò, A. V., Dutton, A. A., & van den Bosch, F. C. 2008, MNRAS, 391, 1940 [NASA ADS] [CrossRef] [Google Scholar]
  83. Makino, N., Sasaki, S., & Suto, Y. 1998, ApJ, 497, 555 [NASA ADS] [CrossRef] [Google Scholar]
  84. Malbon, R. K., Baugh, C. M., Frenk, C. S., & Lacey, C. G. 2007, MNRAS, 382, 1394 [NASA ADS] [Google Scholar]
  85. Martin, C. L. 1999, ApJ, 513, 156 [NASA ADS] [CrossRef] [Google Scholar]
  86. Mo, H. J., Mao, S., & White, S. D. M. 1998, MNRAS, 295, 319 [Google Scholar]
  87. Morganti, R., Oosterloo, T. A., Tadhunter, C. N., van Moorsel, G., & Emonts, B. 2005a, A&A, 439, 521 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  88. Morganti, R., Tadhunter, C. N., & Oosterloo, T. A. 2005b, A&A, 444, L9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  89. Muñoz-Cuartas, J. C., Macciò, A. V., Gottlöber, S., & Dutton, A. A. 2011, MNRAS, 411, 584 [NASA ADS] [CrossRef] [Google Scholar]
  90. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1995, MNRAS, 275, 56 [Google Scholar]
  91. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1996, ApJ, 462, 563 [NASA ADS] [CrossRef] [Google Scholar]
  92. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 [NASA ADS] [CrossRef] [Google Scholar]
  93. Nelson, D., Vogelsberger, M., Genel, S., et al. 2013, MNRAS, 429, 3353 [NASA ADS] [CrossRef] [Google Scholar]
  94. Nesvadba, N. P. H., Lehnert, M. D., Eisenhauer, F., et al. 2006, ApJ, 650, 693 [NASA ADS] [CrossRef] [Google Scholar]
  95. Nesvadba, N. P. H., Lehnert, M. D., De Breuck, C., Gilbert, A. M., & van Breugel, W. 2008, A&A, 491, 407 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  96. Ocvirk, P., Pichon, C., & Teyssier, R. 2008, MNRAS, 390, 1326 [NASA ADS] [Google Scholar]
  97. Okamoto, T., Gao, L., & Theuns, T. 2008, MNRAS, 390, 920 [NASA ADS] [CrossRef] [Google Scholar]
  98. Ostriker, J. P., Choi, E., Ciotti, L., Novak, G. S., & Proga, D. 2010, ApJ, 722, 642 [NASA ADS] [CrossRef] [Google Scholar]
  99. Peebles, P. J. E. 1969, ApJ, 155, 393 [NASA ADS] [CrossRef] [Google Scholar]
  100. Press, W. H., & Schechter, P. 1974, ApJ, 187, 425 [NASA ADS] [CrossRef] [Google Scholar]
  101. Proga, D., Stone, J. M., & Kallman, T. R. 2000, ApJ, 543, 686 [NASA ADS] [CrossRef] [Google Scholar]
  102. Quinn, T., Katz, N., & Efstathiou, G. 1996, MNRAS, 278, L49 [NASA ADS] [CrossRef] [Google Scholar]
  103. Rees, M. J. 1986, MNRAS, 218, 25 [NASA ADS] [CrossRef] [Google Scholar]
  104. Sasaki, M., Heinitz, C., Warth, G., & Pühlhofer, G. 2014, A&A, 563, A9 [Google Scholar]
  105. Sellwood, J. A., & McGaugh, S. S. 2005, ApJ, 634, 70 [NASA ADS] [CrossRef] [Google Scholar]
  106. Shapiro, P. R., Giroux, M. L., & Babul, A. 1994, ApJ, 427, 25 [NASA ADS] [CrossRef] [Google Scholar]
  107. Shi, Y., Helou, G., Yan, L., et al. 2011, ApJ, 733, 87 [NASA ADS] [CrossRef] [Google Scholar]
  108. Shu, C.-G., Mo, H.-J., & Shu-DeMao. 2005, Chin. J. Astron. Astrophys., 5, 327 [Google Scholar]
  109. Silk, J. 2003, MNRAS, 343, 249 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  110. Somerville, R. S. 2002, ApJ, 572, L23 [NASA ADS] [CrossRef] [Google Scholar]
  111. Somerville, R. S., & Kolatt, T. S. 1999, MNRAS, 305, 1 [Google Scholar]
  112. Somerville, R. S., Hopkins, P. F., Cox, T. J., Robertson, B. E., & Hernquist, L. 2008, MNRAS, 391, 481 [NASA ADS] [CrossRef] [Google Scholar]
  113. Somerville, R. S., Gilmore, R. C., Primack, J. R., & Domínguez, A. 2012, MNRAS, 423, 1992 [NASA ADS] [CrossRef] [Google Scholar]
  114. Stoll, R., Mathur, S., Krongold, Y., & Nicastro, F. 2009 [arXiv:0903.5310] [Google Scholar]
  115. Sutherland, R. S., & Dopita, M. A. 1993, ApJS, 88, 253 [NASA ADS] [CrossRef] [Google Scholar]
  116. Suto, Y., Sasaki, S., & Makino, N. 1998, ApJ, 509, 544 [NASA ADS] [CrossRef] [Google Scholar]
  117. Thoul, A. A., & Weinberg, D. H. 1996, ApJ, 465, 608 [NASA ADS] [CrossRef] [Google Scholar]
  118. Toomre, A. 1963, ApJ, 138, 385 [NASA ADS] [CrossRef] [Google Scholar]
  119. Toomre, A. 1964, ApJ, 139, 1217 [NASA ADS] [CrossRef] [Google Scholar]
  120. Tormen, G., Bouchet, F. R., & White, S. D. M. 1997, MNRAS, 286, 865 [NASA ADS] [CrossRef] [Google Scholar]
  121. Tweed, D., Devriendt, J., Blaizot, J., Colombi, S., & Slyz, A. 2009, A&A, 506, 647 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  122. van Daalen, M. P., Schaye, J., Booth, C. M., & Dal la Vecchia, C. 2011, MNRAS, 415, 3649 [NASA ADS] [CrossRef] [Google Scholar]
  123. van de Voort, F., Schaye, J., Booth, C. M., Haas, M. R., & Dal la Vecchia, C. 2011, MNRAS, 414, 2458 [NASA ADS] [CrossRef] [Google Scholar]
  124. van den Bergh, S. 1996, AJ, 112, 2634 [NASA ADS] [CrossRef] [Google Scholar]
  125. van den Bosch, F. C. 2002, MNRAS, 331, 98 [Google Scholar]
  126. Wechsler, R. H., Bullock, J. S., Primack, J. R., Kravtsov, A. V., & Dekel, A. 2002, ApJ, 568, 52 [NASA ADS] [CrossRef] [Google Scholar]
  127. White, S. D. M., & Frenk, C. S. 1991, ApJ, 379, 52 [NASA ADS] [CrossRef] [Google Scholar]
  128. White, S. D. M., & Rees, M. J. 1978, MNRAS, 183, 341 [NASA ADS] [CrossRef] [Google Scholar]
  129. Zhao, D. H., Jing, Y. P., Mo, H. J., & Börner, G. 2009, ApJ, 707, 354 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: An adaptive time-step scheme

In a galaxy, different processes act at the same time but on different time scales. It is not efficient, for example, to compute the evolution of the cold filamentary phase that evolves with a typical dark-matter dynamical time (tdyn = 106 yr) and with a time step dedicated to the ejecta rate of the galaxy (5 × 103 yr). If we use the same time step to follow these two components, we generate numerical errors on the component with the shorter dynamical time.

In this context, each part of the gas cycle (in the baryonic halo phase or in the galaxy) must evolve with a time step that is as close as possible to its dynamical time. Accordingly, we have developed an adaptive time-step scheme. Each component evolves passively by receiving and transfers mass until the evolution of an other component significantly affects its own evolution. We consider that a component is modified if the variation ΔM of its mass (gas or stars) is greater than p% = 10%. This value is low enough to accurately follow star formation and ejection processes in gas-rich galaxies up to z ≃ 4. Naturally, for passive galaxy (like spheroidal at low z), a higherr value could be acceptable.

Appendix A.1: Description

At any given time step each component with a mass M must evolve with an input and output rate (in and out). With these two gas flows, the optimal time step, δtMAX, for a component with a mass M is therefore computed as δtMAX=p%M|inout|·\appendix \setcounter{section}{1} \begin{eqnarray} \delta t^{\rm MAX} = p\%\dfrac{M}{\left|\dot{F}_{\rm in}-\dot{F}_{\rm out}\right|} \cdot \end{eqnarray}(A.1)For example, for the cold gas component Mcold (Fig. (3)), in = flt (Eq. (15)) and out = ff (Eq. (16)); therefore, δtcoldMAX=p%Mcold|fltff|·\appendix \setcounter{section}{1} \begin{eqnarray} \delta t^{\rm MAX}_{\rm cold} = p\%\dfrac{M_{\rm cold}}{\left|\dot{M}_{\rm flt}-\dot{M}_{\rm ff}\right|}\cdot \end{eqnarray}(A.2)At the first time step, the mass is not defined and the reference is then fixed to Mlim,bar=p%fbMlim,dm􏽼0􏽻􏽺0􏽽20×mp.\appendix \setcounter{section}{1} \begin{eqnarray} M_{\rm lim,bar} = p\% \left<f_{\rm b}\right>\underbrace{M_{\rm lim,dm}}_{20\times m_{\rm p}} . \end{eqnarray}(A.3)

Appendix A.2: The hot-gas phase

The case of the hot atmosphere is more complicated. Indeed, the hot phase:

  • receives gas from the intergalactic medium: sh (Eq. (18));

  • follows cooling process: cool (Eq. (54));

  • receives gas coming from the galaxy: wind (Eq. (20)).

This last term is not linked to the mass of the hot phase Mhot but depends on the star formation activity of the galaxy. To follow the evolution of the hot atmosphere with a high accuracy, and therefore compute the optimal time step δthotMAX\hbox{$\delta t^{MAX}_{\rm hot}$}, we must take possible variations in the ejecta rate into account.

In this context, it is important to check whether the evolution of the hot phase is dominated by input (sh + wind) or output (cool) gas flows. In the first case in = sh + wind>out = cool if we take instantaneous values of all input and output rates to compute δthotMAX\hbox{$\delta t^{\rm MAX}_{\rm hot}$}, an increase, in the next time step, in the galaxy ejection rate would lead to a violation, by excess, of the “stability criterium”: ΔMhot/Mhotp%. On the contrary, if the hot phase is dominated by the cooling process in = sh + wind<out = cool, a decrease in the galaxy ejection rate would lead to a larger decrease in the hot gas mass than p%.

We thus have δthotMAX={\appendix \setcounter{section}{1} \begin{eqnarray} \delta t^{\rm MAX}_{\rm hot} = \left\{ \begin{array}{ll} p\%\frac{M_{\rm hot}}{\left|\dot{M}_{\rm sh} + (1+p\%)\dot{M}_{\rm wind} - \dot{M}_{\rm cool}\right|}{:} & {\rm if} \quad \dot{M}_{\rm in} >\dot{M}_{\rm out} \\ p\%\frac{M_{\rm hot}}{\left|\dot{M}_{\rm sh} + (1-p\%)\dot{M}_{\rm wind} - \dot{M}_{\rm cool}\right|}{:} & \mbox{otherwise,} \end{array}\right. \label{dt_hot} \end{eqnarray}(A.4)where we authorize, depending on the case, an increase or a decrease of p% of the instantaneous galaxy ejection rate wind. These two cases are imprinted in two limits (a lower limit: MINwind\hbox{$\dot{M}_{\rm wind}^{\rm MIN}$} and a upper limit: MAXwind\hbox{$\dot{M}_{\rm wind}^{\rm MAX}$}) beyond which the criterion of p% variation is no longer valid.

For the hot-gas phase, with this optimal time step, we can deduce the two lower and upper limits for the galaxy ejection rate. The lower limit reads as MINwind={\appendix \setcounter{section}{1} \begin{eqnarray} \dot{M}^{\rm MIN}_{\rm wind} = \left\{ \begin{array}{ll} p\%\frac{M_{\rm hot}}{\delta t^{\rm MAX}_{\rm hot}}- \dot{M}_{\rm sh} + \dot{M}_{\rm cool}{:} & {\rm if} \quad \dot{M}_{\rm in} >\dot{M}_{\rm out} \\ (1-p\%)\dot{M}_{\rm wind}{:} & \mbox{otherwise} \end{array}\right. \label{M_wind_min} \end{eqnarray}(A.5)and the upper limit as MAXwind={\appendix \setcounter{section}{1} \begin{eqnarray} \dot{M}^{\rm MAX}_{\rm wind} = \left\{ \begin{array}{ll} (1+p\%)dot{M}_{\rm wind} {:} & \mbox{\rm if}\quad \dot{M}_{\rm in} <\dot{M}_{\rm out} \\ p\%\frac{M_{\rm hot}}{\delta t^{\rm MAX}_{\rm hot}}- \dot{M}_{\rm sh} + \dot{M}_{\rm cool}{:} & \mbox{otherwise.} \end{array}\right. \label{M_wind_max} \end{eqnarray}(A.6)Between these two limits, the stability criterion is followed. Figure A.1 illustrates this optimal time step computation. In the case chosen for the figure, the evolution of the hot gas phase is dominated by the input flux (sh + wind>cool). The optimal time step for the hot gas phase evolution is therefore computed by assuming an increase in the galaxy ejecta rate.

Appendix A.3: The stop and restart mechanism

The two galaxy ejection rate limits are MAXwind\hbox{$\dot{M}^{\rm MAX}_{\rm wind}$} and MINwind\hbox{$\dot{M}^{\rm MIN}_{\rm wind}$}. Above MAXwind\hbox{$\dot{M}^{\rm MAX}_{\rm wind}$} and below MINwind\hbox{$\dot{M}^{\rm MIN}_{\rm wind}$}, the evolution of the hot gas phase violates the p% criterion. These two limits are transferred to the galaxy evolution procedure. During the evolution of the galaxy components, if the average6 of the ejecta rate is above or below these two limits, the evolution process is stopped. The properties of the hot gas phase and of the galaxy are computed again at this time, and a new time step starts from this point.

Appendix A.4: Cold and hot phase coupling

In the new eGalICS algorithm, it is therefore possible to compute the optimal time step for the cold (δtcoldMAX\hbox{$\delta t^{\rm MAX}_{\rm cold}$}) and the hot phases (δthotMAX\hbox{$\delta t^{\rm MAX}_{\rm hot}$}). The same computation is performed for all galaxy components, such as for the disc and the bulge components, which composed a galaxy. We present in this section the coupling that exists between the hot and cold phases.

thumbnail Fig. A.1

Diagram illustrating the computation of the optimal time step of the hot gas phase (δthotMAX\hbox{$\delta t^{\rm MAX}_{\rm hot}$}). The two horizontal dash lines materialise the maximum variation, p%. In the case shown here, the evolution of the hot gas phase is dominated by input flows: sh + wind>cool. The optimum time step for the hot-gas phase evolution is therefore computed assuming an increase in the galaxy ejecta rate (1 + p%)wind.

thumbnail Fig. A.2

Scheme presenting the coupled evolution of the cold and the hot phases, and the galaxy. During the evolution, the properties of cold and hot phases allow some barriers tcoldSTOP\hbox{$t_{\rm cold}^{\rm STOP}$} and thotSTOP\hbox{$t_{\rm hot}^{\rm STOP}$} to be established. The evolution of the galaxy is computed following these barriers.

We assume that ff and cool are constant during δtcoldMAX\hbox{$\delta t^{\rm MAX}_{\rm cold}$} and δthotMAX\hbox{$\delta t^{\rm MAX}_{\rm hot}$}, respectively. As shown in Fig. A.2 and as explained previously, the computation of δtcoldMAX\hbox{$\delta t^{\rm MAX}_{\rm cold}$} and δthotMAX\hbox{$\delta t^{\rm MAX}_{\rm hot}$} allows some barriers to be established: tcoldSTOP\hbox{$t_{\rm cold}^{\rm STOP}$} and thotSTOP\hbox{$t_{\rm hot}^{\rm STOP}$}. These barriers mark the evolution of the cold and the hot phase. If the galaxy has already evolved until tgal = tcold = thot, the next galaxy time step is computed following δtgal=MIN(tcoldSTOPthotMAX)tgal.\appendix \setcounter{section}{1} \begin{eqnarray} \delta t_{\rm gal} = {\rm MIN}\left(t_{\rm cold}^{\rm STOP},\delta t^{\rm MAX}_{\rm hot}\right)- t_{\rm gal} . \end{eqnarray}(A.7)During the first part of the evolution presented in Fig. A.2 the time step of the galaxy is limited by the hot gas phase evolution. Indeed, the cold time step is longer (tcoldSTOP>thotSTOP\hbox{$t_{\rm cold}^{\rm STOP} > t_{\rm hot}^{\rm STOP}$}). The hot gas phase and the galaxy therefore evolve on the same time grid. Later on, the galaxy evolution is limited by the cold-phase evolution (tcoldSTOP<thotSTOP,2\hbox{$t_{\rm cold}^{\rm STOP} < t_{\rm hot}^{\rm STOP,2}$}). Indeed, ff has to be recomputed after tcoldSTOP\hbox{$t_{\rm cold}^{\rm STOP}$}. Owing to the strong coupling between the galaxy and the hot gas phase (ejection processes wind), the evolution of these two components are synchronized at tcoldSTOP\hbox{$t_{\rm cold}^{\rm STOP}$}. In these conditions, even if the properties of the hot-gas phase should allow an evolution over a longer time step (until thotSTOP,2\hbox{$t_{\rm hot}^{\rm STOP,2}$}), its evolution is stopped. The three components are updated, their time evolution are then similar (tgal = tcold = thot), and the evolution process is started again.

thumbnail Fig. A.3

Histograms of the time step duration for the galaxy components (disc or bulge in grey), the hot or cold baryonic halo phase (blue), and the dark-matter halos (green). For information, the 94 different snapshots of the N-body simulation used to follow the dark-matter component are distributed with an average time step of 0.15 Gyr.

Appendix A.5: Effective time step duration

Figure A.3 shows the distribution of the time step of the dark-matter component, the cold or/and hot phase and of the galaxy component (disc and bulge). The histograms are built following the analysis of a dark-matter halo and its galaxy between z = 11 and z = 0. The time step durations used most for the hot or cold phases are identical to those used for the dark-matter. Indeed, the bulk of the blue histogram coincides with the green histogram. Minimum values coincides with the synchronisation process (tcoldSTOP<thotSTOP,2\hbox{$t_{\rm cold}^{\rm STOP} < t_{\rm hot}^{\rm STOP,2}$} in Fig. A.2). Concerning the distribution of time steps used for the galaxy components (disc or bulge), the maximum value is linked to the stellar evolution process (δt = 10-3 Gyr, Devriendt et al. 1999). The second bump of the grey histogram corresponds to merger time-scale events (τmerger = 5 × 10-4 Gyr). Lower values are again linked to the synchronisation process (clump transfers) between the disc and the bulge.

Appendix B: The polytropic hot gas case

For a polytropic gas (poly), the temperature distribution is non-isothermal. The pressure, density, and temperature profiles are linked with the following relations Tpoly(r)=T0[ρgpoly(r)ρg,0]Γ1,Ppoly(r)=P0[ρgpoly(r)ρg,0]Γ,\appendix \setcounter{section}{2} \begin{eqnarray} T^{\rm poly}(r)=T_0\left[\dfrac{\rho_{\rm g}^{\rm poly}(r)}{\rho_{\rm g,0}}\right]^{\Gamma -1} , P^{\rm poly} (r)=P_0\left[\dfrac{\rho_{\rm g}^{\rm poly}(r)}{\rho_{\rm g,0}}\right]^{\Gamma} , \label{HEC_poly} \end{eqnarray}(B.1)where T0, P0, and ρg,0 are the central temperature, pressure, and density, respectively, and Γ is the polytropic index of the gas. The HEC in the dark-matter halo potential (Suto et al. 1998; Makino et al. 1998; Komatsu & Seljak 2001; Capelo et al. 2012) gives ρgpoly(r)=ρg,0poly(x,T0)withx=rrdm,\appendix \setcounter{section}{2} \begin{eqnarray} \rho_{\rm g}^{\rm poly}(r) = \rho_{\rm g,0}\mathcal{F}^{\rm poly}(x,T_0)~~\mbox{with}~~x=\frac{r}{r_{\rm dm}} , \label{rho_poly} \end{eqnarray}(B.2)with poly(x,T0), a geometrical function defined as poly(x,T0)=[1ΓΓ14πGμmpρdmrdm2kBT0Φ(x)]1Γ1,\appendix \setcounter{section}{2} \begin{eqnarray} \mathcal{F}^{\rm poly}(x,T_0) = \left[1-\dfrac{\Gamma}{\Gamma-1}\dfrac{4\pi G\mu m_{\rm p}\rho_{\rm dm}r_{\rm dm}^2}{k_{\rm B}T_0}\Phi(x)\right]^{\frac{1}{\Gamma-1}} , \label{F_poly} \end{eqnarray}(B.3)where Φ(x) is the geometrical function defined by Eq. (49) (and see Table B.1 for parameter definitions).

Table B.1

Symbols and their definition used for the polytropic gas description.

In the polytropic gas case, T0 is not equal to the mean temperature Tatm\hbox{$\overline{T}_{\rm atm}$} of the gas phase but to the mass-average temperature Tpoly=Mhot-10rvirρgpoly(r)Tpoly(r)4πr2dr=4πρg,0T0rdm3Mhot0cdm[poly(x,T0)]Γx2dx.\appendix \setcounter{section}{2} \begin{eqnarray} \label{eq_avg_temp_poly} \overline{T}^{\rm poly} & = M_{\rm hot}^{-1}\int_0^{r_{\rm vir}}\rho_{\rm g}^{\rm poly}(r)T^{\rm poly}(r)4\pi r^2{\rm d}r\nonumber\\ & = \dfrac{4\pi \rho_{\rm g,0}T_0 r_{\rm dm}^3}{M_{\rm hot}}\int_0^{c_{\rm dm}}\left[\mathcal{F}^{\rm poly}(x,T_0)\right]^{\Gamma}x^2{\rm d}x . \end{eqnarray}(B.4)Following Eq. (50), the core density ρg,0 is given by ρ0,g=Mhot[4πrdm30cdmpoly(x,T0)x2dx]-1.\appendix \setcounter{section}{2} \begin{eqnarray} \rho_{\rm 0,g} = M_{\rm hot}\left[4\pi r_{\rm dm}^3\int_0^{c_{\rm dm}} \mathcal{F}^{\rm poly}(x,T_0)x^2{\rm d}x \right]^{-1} . \label{rho_gas_core_poly} \end{eqnarray}(B.5)Using Eqs. (B.4) and (B.5), we can compute an implicit relation between T\hbox{$\overline{T}$} and T0Tpoly=T00rvirr0[poly(x,T0)]Γx2dx0rvirr0poly(x,T0)x2dx·\appendix \setcounter{section}{2} \begin{eqnarray} \ds{\overline{T}^{\rm poly} = T_0\dfrac{\int_0^{\frac{r_{\rm vir}}{r_0}}\left[\mathcal{F}^{\rm poly}(x,T_0)\right]^{\Gamma}x^2{\rm d}x}{\int_0^{\frac{r_{\rm vir}}{r_0}}\mathcal{F}^{\rm poly}(x,T_0)x^2{\rm d}x}}\cdot \label{T0_Tpoly} \end{eqnarray}(B.6)where poly(x,T0) (Eq. (B.3)) is a function of the central temperature T0. Since Tpoly\hbox{$\overline{T}^{\rm poly}$} is a strictly increasing function of T0 for a set of dark-matter properties (rdm and ρdm), the implicit relation given above can be inverted using a dichotomy-research algorithm. We can thus compute the central temperature T0 of the polytropic gas connected to the mean temperature Tatm\hbox{$\overline{T}_{\rm atm}$} deduced from energy exchanges (see Sect. (4.2)).

The geometrical profile defined by Eq. (B.3) must be strictly positive to ensure density (Eq. (B.2)) and temperature (Eq. (B.1)) consistency. This condition leads to a minimum threshold value of the central temperature T0min\hbox{$T_0^{\rm min}$}. Indeed all other parameters are linked to dark-matter halo properties and are already fixed. This critical central temperature T0min\hbox{$T_0^{\rm min}$} depends on the geometrical function Φ(x) (Eq. (49)), which is a strictly increasing function of x. As limx+Φ(x)=1,\appendix \setcounter{section}{2} \begin{eqnarray} \underset{x\to+\infty}{\rm lim}\Phi(x) = 1 , \end{eqnarray}(B.7)then T0min>ΓΓ14πGμmpρdmrdm2kB·\appendix \setcounter{section}{2} \begin{eqnarray} T_0^{\rm min} > \dfrac{\Gamma}{\Gamma-1}\dfrac{4\pi G\mu m_{\rm p}\rho_{\rm dm}r_{\rm dm}^2}{k_{\rm B}}\cdot \label{T0_min} \end{eqnarray}(B.8)If the central temperature T0 of the polytropic gas derived from the average temperature Tatm\hbox{$\overline{T}_{\rm atm}$} is lower than this minimum threshold value T0min\hbox{$T_0^{\rm min}$}, the hot atmosphere is not considered stable enough to activate the cooling processes. The gas contained in this hot phase again falls on the galaxy with a feeding rate close to the free-fall rate (as defined for the filamentary structure). We set in this case cool=εffMhot2tdyn,\appendix \setcounter{section}{2} \begin{eqnarray} \dot{M}_{\rm cool} = \varepsilon_{\rm ff}\dfrac{M_{\rm hot}}{2t_{\rm dyn}} , \end{eqnarray}(B.9)where Mhot is the total mass of gas contained in the unstable hot phase, t0 the dynamical time of the dark-matter structure and εff = 1.0 a free parameter (the same as used in Eq. (16)).

Appendix C: Maxwell-Boltzman distribution

We use the following formulation of the well-known Maxwell-Boltzman velocity distribution FMB(v,T)dv=4π(μmp2πkBT)3/2v2exp(μmpv22kB)dv.\appendix \setcounter{section}{3} \begin{eqnarray} F_{\rm MB}(v,T){\rm d}v = 4\pi\left(\dfrac{\mu m_{\rm p}}{2\pi k_{\rm B}T}\right)^{3/2}v^2\exp\left(-\dfrac{\mu m_{\rm p}v^2}{2k_{\rm B}}\right){\rm d}v . \label{MB_dist} \end{eqnarray}(C.1)

Appendix D: Data repository

Outputs from the model are available at the CDS under the reference m1. Data are distributed under *.fits format and are therefore compatible with the TOPCAT software (http://www.star.bris.ac.uk/~mbt/topcat/).

All Tables

Table 1

Dark-matter symbols and their definition.

Table 2

Baryonic accretion symbols and their definition.

Table 3

Ejecta and heating processes symbols and their definition.

Table 4

Symbols and their definition for the hot atmosphere and cooling process.

Table 5

Galaxy symbols and their definition.

Table 6

Disc symbols and their definition.

Table 7

Clumps symbols and their definition.

Table 8

Bulge symbols and their definition.

Table 9

SMBH symbols and their definition.

Table 10

Merger event symbols and their definition.

Table B.1

Symbols and their definition used for the polytropic gas description.

All Figures

thumbnail Fig. 1

Merger tree showing the evolution of the dark-matter smooth accretion. We only show the fifty most massive branches. The junctions with less massive branches are indicated with small arrowheads. The main branch (on the left) is built up through a smooth accretion process acc (Eq. (1)) and multiple mergers with other haloes. The evolution of the smooth dark-matter accretion through the merger tree is shown with colour. A halo without any dark-matter accretion is shown with a black symbol. Circles are used for main haloes (M-H), while squares are used for subhaloes (S-H). The size of the symbols is proportional to the dark-matter halo virial mass. It is clearly visible that before a merger event with a larger structure, a halo is systematically identified as a substructure of this more massive halo. During this subhalo phase, the dark-matter accretion rate strongly decreases and may even stop.

In the text
thumbnail Fig. 2

Evolution of the effective baryonic fraction (fb: Eq. (11)) with dark-matter halo mass and redshift. We use the Okamoto et al. (2008) prescription with a reionisation redshift zreion = 9 (see text for more details). The vertical grey line indicates our mass resolution limit (20 dm-particules). At this mass resolution and for a redshift of about 2, the effective baryonic fraction fb applied to the accretion flux is ΩbΩm\hbox{$\frac{\Omega_{\rm b}}{\Omega_{\rm m}}$}.

In the text
thumbnail Fig. 3

Large-scale exchanges. On a large scale, dark-matter and baryons are coupled. During the accretion process and while the dark-matter smooth accretion (acc) feeds the dark-matter halo (Mh), the associated baryons (bar) feed a baryonic reservoir (Mbar). In the context of bimodal accretion, this baryonic reservoir is divided into two parts: the hot reservoir representing the hot atmosphere (sh and Mhot), and the cold reservoir representing filamentary streams (flt and Mcold). The cold gas feeds the galaxy (ff) directly. The hot gas follows radiative cooling process and then also feeds the galaxy (cool). As described in the text, a galaxy can eject material. The ejecta (wind) are continuously added to the hot-gas reservoir (Mhot) or definitively lost in the IGM according to their velocity distribution. Indeed, the gas stored in the IGM reservoir will never be considered.

In the text
thumbnail Fig. 4

Fraction of hot isotropic gas distributed in the baryonic accretion following (Lu et al. 2011). This fraction depends on the dark-matter halo mass and redshift (colour code).

In the text
thumbnail Fig. 5

The total smooth baryonic accretion onto the galaxy. The colour shading shows the total baryonic accretion rate (free-fall coming from filaments and cooling flows) from our reference model, as a function of dark-matter halo mass, and for a set of redshifts. The colour scale indicates the normalised logarithmic density of objects. We add for information the transition mass (1011.4 M, grey vertical arrow) used in the accretion mode repartition (Eq. (14)). The dotted and solid black lines show the average baryonic accretion rate deduced from hydrodynamic simulations by Fakhouri et al. (2010) and Ceverino et al. (2010), respectively. Black circles are for a model with photoionisation prescription by Gnedin (2000) instead of our Okamoto et al. (2008) standard prescription. Grey squares show the accretion for a model without any photoionisation process.

In the text
thumbnail Fig. 6

Distribution functions of the velocity in the pre-existing hot gas phase and in the wind. The upper and lower panels show the distributions for two different dark-matter haloes with Mvir = 1010 M and Mvir = 1013 M, respectively. In the two cases, the wind speed is close to Vwind = 250 km s-1 but the temperature is close to T=106\hbox{$\overline{T} =10^6$} K in the upper panel and T=107\hbox{$\overline{T}=10^7$} K in the lower panel. The total gas distribution, in black, is the sum of the pre-existing hot gas distribution (red dotted line) and the gas wind distribution (grey dotted line). The blue curve gives the mass fraction. The vertical dotted grey line indicates the escape velocity of the dark-matter structure. The intersection of the mass fraction and this velocity threshold allows us to determine the escape fraction. The red curve shows the new hot gas phase distribution recomputed with the new temperature T\hbox{$\overline{T}$} (Eq. (42)) derived after accretion, cooling, and ejection processes have been taken into account during the previous time step Δt. In the low-mass case, all the mass contained in the wind leaves the structure. In the high-mass case, a fraction of the wind mass can be retained.

In the text
thumbnail Fig. 7

Exchanges on the galaxy scale. The galaxy scale is considered as an intermediate scale. A galaxy can be made of two different parts, a disc \hbox{$\Disc$} and a bulge . The fresh gas (ff and Mcool˙\hbox{$\dot{M_{\rm cool}}$}), coming from the cold and the hot reservoirs, only feeds the disc component. The bulge component is created during a merger or by the clump migration process. SNe and the activity of the SMBH produce ejecta (Mwind). These ejecta are considered at an upper scale (see Fig. 3). The gas exchanges associated with the disc and the bulge are described in Figs. 8 and 11.

In the text
thumbnail Fig. 8

Exchanges on the disc scale. The disc is the main component. It is fed by the fresh gas coming from the large scales (ff and cool). SNe associated with the stellar population of the disc generate ejecta wind,sn. In a disc, dynamics and gravity are in balance. Under some conditions, giant clumps can be formed (\hbox{$\dot{M}_{\rm \Clumps,in}$}) and/or disrupted (\hbox{$\dot{M}_{\rm \Clumps,out}$}). At any given time, a given fraction \hbox{$M_{\Clumps}$} of the total disc mass is in those clumps. The clumps migrate (\hbox{$\dot{M}_{\rm \Clumps,migr}$}) and feed a pseudo-bulge component ().

In the text
thumbnail Fig. 9

Distribution of the \hbox{$\sigma_{r}/V_{\rm circ,\Gal}(r=r_{\Gal,50})$} ratio for 0 <z< 6 (coloured lines). The mean ratio decreases from 0.6 for high-z structures to 0.1 for local structures. The unstable (high sigma disc) population formed at high z disappears progressively and stable disc structures appear.

In the text
thumbnail Fig. 10

Evolution of the disc-stability criterion (Toomre parameter Q) and clump properties. In the upper panel is shown the evolution of the Toomre parameter with z (Eq. (68), grey line). On this timeline, vertical arrows indicate a merger event (with short and long arrows for minor and major mergers, respectively). For 11 <z< 1.5, the Toomre parameter is larger than the critical value Q>Qcrit = 1, and therefore the disc is stable, but for 1.5 <z< 0.25 (highlighted by the two vertical black lines), the disc becomes unstable. During the unstable transition, bounded by two merger events, many clumps are formed and migrate. The lower panel shows a zoom of this unstable period. The blue and green lines show the number of clumps (left y-axis) and the average individual mass of these clumps (right y-axis), respectively. After a merger event, disc instabilities lead to the formation of a large number of clumps (42) with small individual masses (\hbox{$M_{\rm \Clumps,ind} < 10^5~\Msun$}). Then, during a series of merger events, clumps are destroyed and formed again, but with higher individual masses. Between two merger events, the number of clumps increases while the individual mass \hbox{$M_{\rm \Clumps,iind}$} remains stable. A decrease in one or two clumps is linked to the migration process. It is interesting to note that the migration process creates a temporary stable disc (Q>Qcrit = 1). In the second part of the unstable period, only one clump is formed but with a much higher mass. Finally, this single clump is disrupted and the disc becomes stable again.

In the text
thumbnail Fig. 11

Exchanges on the bulge scale. The bulge is considered as a passive component not fed by fresh accreted gas. A bulge () is formed during a merger or by clumps migrations from the disc (\hbox{$\Disc$}). An SMBH can evolve in the bulge. It is fed by a gas reservoir (the torus, \hbox{$\Torus$}) with an accretion rate acc. The torus is fed during a merger or by the clump migration process. The activity of the SMBH generates ejecta (•,out).

In the text
thumbnail Fig. 12

Fraction of structures hosting an SMBH in the centre of their galaxy as a function of the dark-matter halo mass for a set of redshift (coloured lines). The fraction increases with the dark-matter halo mass to reach a value of 1 for haloes more massive than ≃ 4 × 1011 M. The insert panel shows the fraction of active galaxy nuclei as a function of the dark-matter halo mass. An SMBH is considered to be active if its bolometric luminosity is higher than 5% of its Eddington luminosity.

In the text
thumbnail Fig. 13

Cartoon showing the stellar cycle. Stars are formed from cold gas in the galaxy components (). Stars evolve and die. We follow the evolution of the ejected mass produced by the stellar population wind,.

In the text
thumbnail Fig. A.1

Diagram illustrating the computation of the optimal time step of the hot gas phase (δthotMAX\hbox{$\delta t^{\rm MAX}_{\rm hot}$}). The two horizontal dash lines materialise the maximum variation, p%. In the case shown here, the evolution of the hot gas phase is dominated by input flows: sh + wind>cool. The optimum time step for the hot-gas phase evolution is therefore computed assuming an increase in the galaxy ejecta rate (1 + p%)wind.

In the text
thumbnail Fig. A.2

Scheme presenting the coupled evolution of the cold and the hot phases, and the galaxy. During the evolution, the properties of cold and hot phases allow some barriers tcoldSTOP\hbox{$t_{\rm cold}^{\rm STOP}$} and thotSTOP\hbox{$t_{\rm hot}^{\rm STOP}$} to be established. The evolution of the galaxy is computed following these barriers.

In the text
thumbnail Fig. A.3

Histograms of the time step duration for the galaxy components (disc or bulge in grey), the hot or cold baryonic halo phase (blue), and the dark-matter halos (green). For information, the 94 different snapshots of the N-body simulation used to follow the dark-matter component are distributed with an average time step of 0.15 Gyr.

In the text

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

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

Initial download of the metrics may take a while.