19.3. Galaxy growth by gas accretion

\label{chapter-merger-gas}

The extended Press-Schechter formalism or dark-matter-only \(N\)-body simulations are highly effective at determining the merger and accretion history of dark-matter halos. But, without embellishments, they do not distinguish between dark matter and baryons—stars and gas—and treat all matter as collisionless, interacting only through gravity. In real galaxies, however, the baryons are always observed to be much more centrally concentrated than the dark matter. So some form of dissipation must occur for the baryons to cool to the centers of halos.

Before virialization, baryons and dark matter trace each other, but they start to behave differently during virialization. While the dark matter remains collisionless as shells of matter start to cross in the formation of a halo, gas flows cannot intersect, pressure becomes important, and kinetic energy is exchanged for internal energy. This happens first at the center of the forming halo, where the kinetic energy is fully transformed into internal energy, heating up the gas, and causing the bulk velocity of the gas to vanish. Under the right conditions, a shock then propagates outward from the center and the entire gaseous halo’s bulk velocity is transformed into internal energy as the gas flows through the shock front (Rees & Ostriker 1977, Birnboim & Dekel 2003). A shock forms because the gas is moving supersonically: the virializing gas moves at approximately the circular velocity \(v_c\) of the forming halo—\(\mathcal{O}(100\,\mathrm{km\,s}^{-1})\) for the halos that host large galaxies—but is relatively cold, with a sound speed \(c_s \approx 11\,\mathrm{km\,s}^{-1}\,\sqrt{T/ 10^4\,\mathrm{K}} \ll v_c\) (Equation 17.39). However, for the shock to propagate outwards all the way to the virial radius, the pressure behind the shock must remain high enough; this is not the case if the gas behind the shock cools too fast. In that case, the shock stalls near the center, well before reaching the halo’s edge. When the shock is able to propagate throughout the halo, the gaseous halo is shock-heated to high temperature out to the virial radius (see below), gaining a large amount of internal energy and accompanying pressure support that must be radiated away for the gas to cool into a star-forming disk. The shock front further heats any accreted gas after the formation of the halo and its location at all times essentially coincides with the virial radius of the growing halo. When the shock is unable to propagate, gas can flow more easily to the central region of the halo and only shock-heats near the center, both in the initial formation of the halo and in subsequent accretion. The first scenario where gas is shock-heated throughout the halo is known as hot-mode accretion, while the latter scenario is referred to as cold-mode accretion. We discuss these in turn in more detail in what follows. Because radiative cooling is important in both scenarios, we first take a brief detour to describe the basics of cooling in gaseous halos. Finally, not all of the gas that is able to cool turns into stars because of various feedback processes, which we discuss in the final part of this section.

19.3.1. Radiative cooling in gaseous halos

\label{chapter-merger-gas-cooling}

The gas in shock-heated halos cools through a number of two-body interactions that cause the emission of photons; the gas in gaseous halos is largely transparent to these photons, such that they efficiently carry away energy from the gas. In addition to cooling, the gas is heated through photo-ionization and Compton scattering by the ultraviolet background and this has important effects on the net cooling at different redshifts. We start by going through the various relevant cooling processes and then discuss the effect of the ultraviolet background to arrive at a realistic cooling picture relevant to galaxy formation.

There are four main cooling processes in gaseous galactic halos. These all involve an electron and a neutral or ionized atom. The electron is characterized as ‘bound’ when it is bound to a nucleus and ‘free’ otherwise and the processes are categorized by whether or not the electron starts out or ends up free or bound. Discussing these processes in detail is outside of the scope of this book and we refer the reader to specialized texts such as Rybicki & Lightman (1979) for more information.

  • Free-free emission, also known as Bremsstrahlung: A free electron passing by an ionized nucleus is significantly decelerated and, thus, radiates. Free-free emission dominates the cooling at temperatures when all relevant elements are largely ionized (\(\approx 10^6\,\mathrm{K}\) for primordial gas and \(\approx 3\times 10^7\,\mathrm{K}\) for enriched gas).

  • Free-bound emission, or recombination: A photon is emitted as an electron recombines with an ionized nucleus. A further photon can be emitted when the electron is captured in an excited state and decays to the ground state.

  • Bound-free emission, or collisional ionization: An electron and an atom interacting causes the latter to be ionized, radiating a photon in the process.

  • Bound-bound emission, or collisional excitation: Similar to collisional ionization, but an electron in the atom is merely excited to a higher level and a photon is emitted when this electron drops back to its lower level.

The cooling rate \(\mathcal{C}(T,n_H,Z)\) is the resulting net energy loss per unit time and unit volume and it depends on the temperature \(T\), the number density \(n_H\), and the composition \(Z\) of the gas. Because the emission processes described above all involve the interaction of two particles, their rates all scale as the number density squared (or, better, as the density of electrons times the density of atoms, but in most relevant situations these two are proportional to each other). In practice we use the number density \(n_H\) of hydrogen atoms; the total number density is proportional to this with a factor that depends on the composition. Thus, we can write \(\mathcal{C}(T,n,Z)\) in terms of a cooling function \(\Lambda(T,Z)\) as (e.g., Sutherland & Dopita 1993) \begin{equation}\label{eq-hierarchgal-coolingfunc-def} \mathcal{C}(T,n_H,Z) = n_H^2\,\Lambda(T,Z)\,, \end{equation} because \(\Lambda(T,Z)\) then does not depend on the number density. This will no longer be the case when we take into account the ultraviolet background radiation.

Computing cooling functions requires data on atomic levels and reaction cross sections for many elements and they require assumptions about how the levels and ionization states are populated. It is generally a good assumption to assume that the gas is in ionization equilibrium where the densities of all ions are equal to their equilibrium values. Then we can compute cooling functions for different temperatures and compositions using codes like Cloudy (Ferland et al. 1998; Ferland et al. 2017). When only the collisional processes discussed above are considered, the equilibrium assumed is collisional ionization equilibrium, but we can also include photo-ionization and other photon scattering processes. We have done this for some parameter settings relevant to galaxy formation and we load and discuss the results here.

In detail, we have computed cooling rates \(\mathcal{C}\) (which can easily be translated to cooling functions), using Cloudy version 17.03 (Ferland et al. 2017). We computed cooling functions assuming collisional ionization equilibrium for primordial or enriched gas using input parameter files like:

database H2
database CHIANTI levels maximum
database stout levels maximum
database lamda levels maximum
abundances "primordial.abn" // or 'metals scale factor 1' for solar abundances
coronal t=3 vary
grid, range from 3 to 10 with 0.025 dex increments sequential
hden 2
stop zone 1
set dr 0
set WeakHeatCool 0.001
cosmic rays background
save grid ".grd"
save cooling ".col"

We have added the cosmic-ray background here to allow the code to converge at low temperatures, but in practice it does not affect the result in the regime at \(T > 10^4\,\mathrm{K}\) in which we are interested. To add photo-ionization due to the ultraviolet background from Haardt & Madau (2012) (and black-body radiation from the CMB), we add:

table HM12 redshift [redshift]
cmb redshift [redshift]

where [redshift] is the redshift at which we evaluate the background radiation and we additionally change the hden line such that it uses the (log-base-ten) mean virial hydrogen density at the relevant redshift (see below).

The following code reads the pre-computed cooling functions that we will explore:

[6]:
import sys
from pathlib import Path
import os
import shutil
import tempfile
import subprocess
import gzip
from astropy.table import QTable
_CACHE_BASEDIR= Path(os.getenv('HOME')) / '.galaxiesbook' / 'cache'
_CLOUDY_DIR= Path('cloudy')
_CACHE_CLOUDY_DIR= _CACHE_BASEDIR / _CLOUDY_DIR
_ERASESTR= "                                                                                "
def read_cooling(metals='primordial',
                 read_fracs=False,scale_factor=1e-4,
                 hr=False,
                 photo=False,redshift=0):
    """Read cooling functions computed using Cloudy v17.03

    metals= ('primordial') either 'primordial' or a Z value

    read_fracs= (False) if True, also return the fractional contributions to the cooling from different processes

    scale_factor= (1e-4) 1/nH^2: We computed the basic cooling curves at nH = 100 / cm^3, so need to scale by 1/nH^2 = 1e-4 by default

    hr= (False) if True, read from the high-resolution Z grid

    photo= (False) if True, include photoionization from the UV background and the CMB; in this case, also specify the redshift

    redshift= (0) redshift when using photo=True
    """
    if metals == 'primordial':
        metals_str= 'primordial'
    else:
        metals_str= f"{metals}_solar"
    if photo:
        photo_str= "photo_"
        redshift_str= f"_redshift_{redshift}"
    else:
        photo_str= ""
        redshift_str= ""
    basename= f'cooling_{photo_str}{metals_str}{redshift_str}.col'
    dirname= 'cooling_funcs'
    if hr:
        dirname+= '_hr'
    elif photo:
        dirname+= '_photo'
    filePath= _CACHE_CLOUDY_DIR / dirname /basename
    if not filePath.exists():
        download_cloudy_file( Path(dirname) / basename)
    # Now read the file
    T= []
    C= []
    tmp_fracs= []
    with open(filePath) as infile:
        for line in infile:
            if line[0] == '#': continue
            entries= line.split()
            T.append(float(entries[1]))
            C.append(float(entries[3])-float(entries[2])) # net cooling - heating
            if not read_fracs: continue
            ii= 4
            this_fracs= {}
            while ii < len(entries):
                # Find end of this frac entry
                jj= ii+1
                key= entries[ii]
                while jj < len(entries):
                    try:
                        float(entries[jj])
                    except ValueError:
                        key+= f"_{entries[jj]}"
                        jj+= 1
                    else:
                        break
                while jj < len(entries):
                    try:
                        float(entries[jj])
                    except ValueError:
                        break
                    else:
                        jj+= 1
                if key in this_fracs:
                    this_fracs[key]+= float(entries[jj-1])
                else:
                    this_fracs[key]= float(entries[jj-1])
                ii= jj
            tmp_fracs.append(this_fracs)
    T= numpy.array(T)
    Lambda= numpy.array(C)*scale_factor
    if not read_fracs:
        return T,Lambda
    # Post-processes fractions
    # magic to flatten the keys list
    frac_keys= set([k for ks in [list(f.keys()) for f in tmp_fracs] for k in ks])
    cool_fracs= {}
    for key in frac_keys:
        this_frac= numpy.zeros_like(C)
        for ii in range(len(C)):
            try:
                this_frac[ii]= tmp_fracs[ii][key]
            except KeyError:
                pass
        cool_fracs[key]= this_frac
    return T,Lambda,cool_fracs
def download_cloudy_file(filename):
    filePath= _CACHE_CLOUDY_DIR / filename
    if filePath.exists(): return None
    try:
        # make all intermediate directories
        os.makedirs(filePath.parent)
    except OSError: pass
    sys.stdout.write('\r'+"Downloading file %s ...\r" \
                         % (filePath.name))
    sys.stdout.flush()
    # Safe way of downloading
    downloading= True
    interrupted= False
    file, tmp_savefilename= tempfile.mkstemp()
    os.close(file) #Easier this way
    ntries= 1
    downloadPath= Path('https://galaxiesbook.org/data') \
                    / _CLOUDY_DIR / filename
    while downloading:
        try:
            subprocess.check_call(['curl','-L',downloadPath,
                                   '-o',tmp_savefilename,
                                   '--fail','--silent','--show-error',
                                   '--connect-timeout','10',
                                   '--retry','3'])
            shutil.move(tmp_savefilename,filePath)
            downloading= False
            if interrupted:
                raise KeyboardInterrupt
        except:
            if not downloading: #Assume KeyboardInterrupt
                raise
            elif ntries > 5:
                raise IOError('File %s does not appear to exist on the server ...' % (filePath.name))
        finally:
            if Path(tmp_savefilename).exists():
                os.remove(tmp_savefilename)
        ntries+= 1
    sys.stdout.write('\r'+_ERASESTR+'\r')
    sys.stdout.flush()
    return None

and the following function plots the cooling function:

[7]:
def plot_cooling(metals='primordial',
                 scale_factor=1e-4,
                 hr=False,
                 photo=False,redshift=0,
                 *plot_args,**plot_kwargs):
    T,Lambda= read_cooling(metals=metals,
                           read_fracs=False,
                           scale_factor=scale_factor,
                           hr=hr,photo=photo,
                           redshift=redshift)
    loglog(T,Lambda,*plot_args,**plot_kwargs)
    xlabel(r'$T\,(\mathrm{K})$')
    return None

First we look at the cooling function for gas with a primordial composition of only hydrogen and helium, specifically a mass fraction of helium of \(Y = 0.2477\) (Planck Collaboration et al. 2014). We show the overall cooling function in Figure 19.12 and also display contributions from the different collisional processes discussed above split by the identity of the involved atom by parsing the output from Cloudy.

[8]:
def plot_cooling_components_primordial():
    T,Lambda,cool_fracs= read_cooling(metals='primordial',
                                      read_fracs=True)
    loglog(T,Lambda,'k',lw=2.)
    for key in cool_fracs:
        if numpy.all(cool_fracs[key] < 1e-2): continue # too small
        if key == 'eeff' or key == 'H2cX': continue # small
        # Determine linestyle
        if 'ion' in key:
            ls= ':'
            lw= 3.
        elif 'lin' in key:
            ls= '--'
            lw= 2.
        elif 'rcol' in key:
            ls= '-.'
            lw= 2.
        else:
            ls= '-'
            lw= 1.5
        # Determine color
        if 'H_He' in key:
            color= get_cmap("tab10")(0)
        elif 'H_H' in key:
            color= get_cmap("tab10")(1)
        elif 'HeHe' in key:
            color= get_cmap("tab10")(2)
        else:
            color= get_cmap("tab10")(3)
        loglog(T,Lambda*cool_fracs[key],
                    color=color,
                    label=rf"{key}",
                    ls=ls,lw=lw)
    xlim(8e3,3e8)
    ylim(1e-25,5e-22)
    xlabel(r'$T\,(\mathrm{K})$')
    ylabel(r'$\Lambda\,(\mathrm{erg\,s}^{-1}\,\mathrm{cm}^{3})$')
    # Legend
    import matplotlib.lines as mlines
    excitation_line= mlines.Line2D([], [], color='k',marker=None,ls='--',lw=2.,
                                   label=r'$\mathrm{excitation}$')
    ionization_line= mlines.Line2D([], [], color='k',marker=None,ls=':',lw=3.,
                                   label=r'$\mathrm{ionization}$')
    recombination_line= mlines.Line2D([], [], color='k',marker=None,ls='-.',lw=2.,
                                      label=r'$\mathrm{recombination}$')
    freefree_line= mlines.Line2D([], [], color=get_cmap("tab10")(3),
                                marker=None,ls='-',lw=1.5,
                                label=r'$\mathrm{free\!-\!free}$')
    gca().add_artist(legend(handles=[excitation_line,
                                     ionization_line,
                                     recombination_line,
                                     freefree_line],
                            loc='upper right',
                            fontsize=18.,
                            frameon=False))
    HH_line= mlines.Line2D([], [],color=get_cmap("tab10")(1),
                           marker=None,ls='-',lw=2.,
                           label=r'$\mathrm{H}$')
    HeHe_line= mlines.Line2D([], [],color=get_cmap("tab10")(2),
                           marker=None,ls='-',lw=2.,
                           label=r'$\mathrm{He}$')
    HHe_line= mlines.Line2D([], [],color=get_cmap("tab10")(0),
                           marker=None,ls='-',lw=2.,
                           label=r'$\mathrm{He}^+$')
    gca().add_artist(legend(handles=[HH_line,
                                     HeHe_line,
                                     HHe_line],
                            loc=(0.65,0.2),
                            fontsize=18.,
                            frameon=False))
    return None
figure(figsize=(6,6))
plot_cooling_components_primordial()
../_images/chapters_IV-03.-Hierarchical-Galaxy-Formation_100_0.svg

Figure 19.12: The cooling function of primordial gas.

We see that the cooling function has four main features: a cut-off below \(T \approx 10^4\,\mathrm{K}\); peaks at \(T \approx 2\times 10^4\,\mathrm{K}\) and \(T \approx 10^5\,\mathrm{K}\); and an increasing cooling function at high temperatures. The cut-off at \(T \approx 10^4\,\mathrm{K}\) is caused by the fact that, at these temperatures, hydrogen and helium are entirely neutral, such that there are no free electrons to participate in the cooling processes. The high-temperature behavior is dominated by free-free emission with its \(\Lambda(T) \propto T^{1/2}\) behavior, because the other cooling processes stop occurring when hydrogen and helium are fully ionized. The decomposition of the cooling function into its constituents demonstrates that the peaks are caused by collisional excitation of hydrogen (for the lower temperature peak) and of ionized helium; these processes stop occurring when hydrogen (helium) becomes fully (doubly) ionized. Collisional ionization and collisional recombination as well as processes involving neutral helium are sub-dominant contributors at all temperatures.

For enriched gas, the cooling function is enhanced by the fact that metals have many more excited levels to participate in collisional excitation. As an example, we compute the cooling function using Cloudy’s default set of solar abundances and split it into contributions from different elements. Parsing Cloudy's output to perform this split and plotting the result is done using the following functions:

[9]:
import copy
def parse_cool_fracs_metals(cool_fracs):
    cool_fracs= copy.deepcopy(cool_fracs)
    out_cool_fracs= {}
    metals= ['Ne','Na','Mg','Al','Si','P',
             'Cl','Ar','K','Ca','Ti','V','Cr','Mn','Fe','Co','Ni',
             'C','N','O']
    # For ordering below
    metals_in_order= ['C','N','O','F','Ne','Na','Mg','Al','Si','P','S','Cl',
                      'Ar','K','Ca','Sc','Ti','V','Cr','Mn','Fe','Co','Ni']
    keys_to_delete= []
    for key in cool_fracs.keys():
        for metal in metals:
            if metal in key:
                if metal in out_cool_fracs.keys():
                    out_cool_fracs[metal]+= cool_fracs[key]
                else:
                    out_cool_fracs[metal]= cool_fracs[key]
                keys_to_delete.append(key)
                break
        # Do 'S' separately
        if key == 'S' or key.endswith('H_S') or key.endswith('HeS'):
            if 'S' in out_cool_fracs.keys():
                out_cool_fracs['S']+= cool_fracs[key]
            else:
                out_cool_fracs['S']= cool_fracs[key]
            keys_to_delete.append(key)
        # Do H/He
        if 'H_H' in key or 'HeHe' in key or 'H_He' in key:
            if 'HHe' in out_cool_fracs.keys():
                out_cool_fracs['HHe']+= cool_fracs[key]
            else:
                out_cool_fracs['HHe']= cool_fracs[key]
            keys_to_delete.append(key)
    for key in keys_to_delete:
        del cool_fracs[key]
    # Copy over those that weren't changed
    for key in cool_fracs.keys():
        out_cool_fracs[key]= cool_fracs[key]
    # Re-order for plotting: elements in order and HHe and FF last
    metals_in_order.extend(['HHe','FF_c'])
    return {k: out_cool_fracs[k] for k in metals_in_order
            if k in out_cool_fracs.keys()}
def plot_cooling_components_metals(metals,include_primordial=False,
                                   scale_factor=1e-4):
    T,Lambda,cool_fracs= read_cooling(metals=metals,read_fracs=True,
                                      scale_factor=scale_factor)
    cool_fracs= parse_cool_fracs_metals(cool_fracs)
    loglog(T,Lambda,'k',lw=2.5)
    for key in cool_fracs:
        if numpy.all(cool_fracs[key] < 1e-2): continue # too small
        if key == 'eeff' or key == 'H2cX': continue # small
        if key == 'HHe':
            color= 'k'
            ls= '--'
            label= r'$\mathrm{H/He}$'
        elif key == 'FF_c':
            color= 'k'
            ls= '-.'
            label= r'$\mathrm{free\!-\!free}$'
        else:
            color= None
            ls= '-'
            label=rf"$\mathrm{{{key}}}$"
        loglog(T,Lambda*cool_fracs[key],
                    label=label,
                    color=color,ls=ls,lw=1.5)
    if include_primordial:
        plot_cooling(metals='primordial',label=r'$Z=0$',
                     ls='-',lw=4.,color='0.8',zorder=-100)
    xlim(8e3,3e8)
    ylim(7e-25,4e-21)
    xlabel(r'$T\,(\mathrm{K})$')
    ylabel(r'$\Lambda\,(\mathrm{erg\,s}^{-1}\,\mathrm{cm}^{3})$')
    legend(frameon=False,fontsize=12.,ncol=4)
    return None

and then we plot the cooling function in the left panel of the following figure:

[10]:
from matplotlib.gridspec import GridSpec
figure(figsize=(11.5,5.5))
gs= GridSpec(1,3,width_ratios=[0.5,0.5,0.02])
subplot(gs[0])
# You can also use metals=0.001, 0.01, 0.1, or 10.
plot_cooling_components_metals(metals=1.,include_primordial=True)
ylim(3e-24,9e-21)
# Also plot the cooling function as a function of metallicity
subplot(gs[1])
def read_cooling_hr_metals():
    basename= "cooling_hr_metals.csv"
    dirname= "cooling_funcs_hr"
    filePath= _CACHE_CLOUDY_DIR / dirname /basename
    if not filePath.exists():
        download_cloudy_file( Path(dirname) / basename)
    return numpy.loadtxt(filePath)
metals_hr= read_cooling_hr_metals()
from matplotlib.cm import viridis as cmap
MHs= numpy.log10(numpy.array(metals_hr,dtype='float'))
# Everything below -3 is pretty much to same, so no point
# in downloading all those curves
indx= (MHs >= -3.)+(MHs % 1 == 0.)
metals_hr= metals_hr[indx]
MHs= MHs[indx]
for Zstr,MH in zip(metals_hr,MHs):
    plot_cooling(metals=Zstr,hr=True,
                 color=cmap((MH-numpy.amin(MHs))/(numpy.amax(MHs)-numpy.amin(MHs))),
                 lw=4.5,zorder=-MH)
# Add 'contours' of MH = -1, 0, 1
for Z in [1./10.,1.,10.]:
    plot_cooling(metals=Z,lw=2.5,color='0.6',zorder=100)
xlim(8e3,3e8)
ylim(3e-24,9e-21)
gca().set_yticklabels([])
# Hack colorbar
subplot(gs[2])
import matplotlib
matplotlib.colorbar.ColorbarBase(gca(), orientation='vertical',
                               cmap='viridis',
                               norm=matplotlib.colors.Normalize(numpy.amin(MHs),numpy.amax(MHs)),
                               label=r'$[\mathrm{M/H}]$')
#img= imshow(numpy.array([[numpy.amin(MHs),numpy.amax(MHs)]]),cmap=cmap,aspect=2.)
#img.set_visible(False)
#colorbar(shrink=1.,label=r'$[\mathrm{M/H}]$')
tight_layout(w_pad=1.);
../_images/chapters_IV-03.-Hierarchical-Galaxy-Formation_107_0.svg

Figure 19.13: The cooling function and its components for solar-metallicity gas (left panel) and at different metallicities (right panel).

We see that cooling is significantly increased by up to almost two orders of magnitude at all temperatures \(10^4\,\mathrm{K} \lesssim T \lesssim 3\times 10^7\,\mathrm{K}\). There is an especially large effect near \(T \approx 10^5\,\mathrm{K}\) because of the many excited levels of oxygen and carbon, respectively the third and fourth most common element by mass in the solar composition. Heavier atoms require higher temperatures to be fully ionized and collisional excitation provides significant cooling up to \(T \approx 3\times 10^7\,\mathrm{K}\) due to heavy elements such as iron (the sixth most common element by mass in the solar composition). At higher temperatures, all elements are again fully ionized and the cooling function simplifies to the free-free component.

The effect of metallicity on the cooling function can be seen in more detail in the right panel of Figure 19.13, which displays the cooling function for metallicities \(Z\) in the range \(10^{-6}\,Z_\odot \leq Z \leq 10\,Z_\odot\); to help guide the eye, we plot the cooling functions for \(Z = 0.1\,Z_\odot\), \(Z_\odot\), and \(10\,Z_\odot\) as curves. We see that for low metallicities, \(Z \lesssim 10^{-2}\,Z_\odot\), the cooling function is essentially unchanged from the primordial one. But at higher metallicities there is a significant increase. As the intergalactic medium in the Universe gets enriched by galactic outflows, the gas that flows into galaxies and provides the fuel for star formation has an easier time cooling.

So far we have ignored the effect of photo-ionization and Compton heating by the ultraviolet background. But because of the ionizing radiation from, first, massive stars and, later, quasars, there is a significant ultraviolet background at all redshifts at which galaxies form (e.g., Haardt & Madau 1996). Photo-ionization causes heating and changes the cooling rate, because it increases the ionization fraction of the gas, increasing the importance of free-free emission and changing the rate at which the other collisional cooling processes happen. Because of the presence of heating, we now define the cooling function using the effective cooling rate as \(\mathcal{C}-\mathcal{H} = n_H^2\,\Lambda(T,Z,n_H)\), where \(\mathcal{H}\) is the heating rate. Because photo-ionization involves only a single atom, its rate depends on the number density \(n\) of gas particles, while the cooling rate from the collisional processes depends on \(n_H^2\). The resulting cooling function is then no longer independent of \(n\), which means that we have to specify the number density in its computation. For detailed calculations of the cooling of gaseous halos that takes into account internal density gradients, one would have to compute the full three-dimensional \(\Lambda(T,Z,n_H)\) behavior. But to gain an overall understanding of the cooling of gaseous halos, we can simply compute \(\Lambda(T,Z,n_H)\) at the mean hydrogen density of galactic halos, which in the standard picture of spherical collapse and virialization only depends on redshift, but not on halo mass or other properties (see Chapter 17.3). Thus, to illustrate the effect of the ultraviolet background, we compute \(\Lambda(T)\) for primordial gas as a function of redshift, including the ultraviolet background from Haardt & Madau (2012) for that redshift and using the mean hydrogen density within a virialized halo obtained using the virial overdensity from Equation (17.104), computed for the Planck Collaboration et al. (2020b) cosmological parameters.

The following function computes this mean hydrogen density:

[11]:
from astropy.cosmology import Planck18
from astropy.constants import m_p
from astropy import units
def mean_hydrogen_density(redshift,log10=True):
    mean_matter= Planck18.critical_density(redshift)*Planck18.Om(redshift)
    def Deltav(redshift):
        x= Planck18.Om(redshift)-1.
        return (18.*numpy.pi**2.+82.*x-39*x**2.)/Planck18.Om(redshift)
    # Assume primordial composition, so H is 75% by mass
    out= (0.75*Planck18.Ob(redshift)/Planck18.Om(redshift)*Deltav(redshift)*mean_matter/m_p).to(1/units.cm**3)
    if log10:
        return numpy.log10(out.to_value(1/units.cm**3))
    else:
        return out

and we then plot the resulting cooling function that includes the effect of the ultraviolet background:

[12]:
from matplotlib.cm import plasma as cmap
from matplotlib import colors, cm
figure(figsize=(8,6))
redshifts= numpy.arange(0,16,1)
for redshift in redshifts:
    plot_cooling(photo=True,redshift=redshift,
                 scale_factor=1./mean_hydrogen_density(redshift,log10=False)**2.,
                 color=cmap(redshift/numpy.amax(redshifts)))
plot_cooling(metals='primordial',color='#1f77b4',lw=7.5,zorder=-100)
ylabel(r'$\Lambda\,(\mathrm{erg\,s}^{-1}\,\mathrm{cm}^{3})$')
xlim(1e4,1e8)
ylim(1e-24,5e-22)
ax= gca()
# Add colorbar
vmin, vmax= redshifts[0], redshifts[-1]
tcmap= colors.ListedColormap([cmap(f) for f in numpy.linspace(0.,1.,len(redshifts))])
dtick= (redshifts[1]-redshifts[0])
sm = cm.ScalarMappable(cmap=tcmap,norm=Normalize(vmin=vmin-0.5*dtick,vmax=vmax+0.5*dtick))
sm._A = []
cbar= colorbar(sm,use_gridspec=True,ax=ax,format=r'$%.0f$',ticks=redshifts)
cbar.set_label(r'$\mathrm{redshift}$')
../_images/chapters_IV-03.-Hierarchical-Galaxy-Formation_114_0.svg

Figure 19.14: The cooling function of primordial gas in the presence of the ultraviolet background.

We have included the primordial cooling curve without the ultraviolet background as the thick line. We see that at high redshift, the ultraviolet background significantly increases the contribution from free-free emission at high temperatures, but leaves the cooling by collisional excitation at low temperatures unchanged. Both of these behaviors result from the high density of halos at high redshift: the density is so high that collisional excitation is efficient, but Compton scattering of the ultraviolet-background photons off the free electrons increases the cooling by free-free emission. At lower redshifts where the mean hydrogen density inside halos is smaller, photo-ionization significantly impedes cooling by collisional excitation, removing the low-temperature, hydrogen peak in the cooling function entirely at \(z \lesssim 10\) and significantly decreasing the second, ionized-helium peak at low redshift; cooling by free-free emission at low redshift is approximately the same as in the absence of an ultraviolet background. Photo-ionization affects the cooling function of enriched gas in a similar manner, reducing the cooling function by up to an order of magnitude at \(10^4\,\mathrm{K} \lesssim T \lesssim 10^6\,\mathrm{K}\) (Wiersma et al. 2009).

Thus, we see that the ultraviolet background has a significant effect on how gaseous halos cool, making it difficult at low redshifts to cool halos with temperatures below \(10^5\,\mathrm{K}\). However, we caution once more that when the ultraviolet background is important, it is doubly important to take the density structure of the gaseous halo into account. Even if the entire halo is unable to cool efficiently, cooling in the higher-density center may still proceed quickly, because the effect of photo-ionization is much smaller at higher densities.

19.3.2. Hot-mode and cold-mode accretion

\label{chapter-merger-gas-hotcoldaccretion}

In both hot- and cold-mode accretion, the gas that forms a halo or is accreted onto it comes from the intergalactic medium (IGM). After recombination, the temperature of the intergalactic gas is determined by the trade-off between adiabatic cooling and Compton heating from the cosmic microwave background (CMB). The resulting temperature evolution can be computed using tools such as RECFAST (Seager et al. 1999). The result is that until \(z \approx 500\), the gas temperature tracks the behavior of the CMB, \(T_\mathrm{gas} = T_\mathrm{CMB} = T_{0,\mathrm{CMB}}/a\), where \(T_{0,\mathrm{CMB}} = 2.7255\,\mathrm{K}\) is the CMB temperature at redshift zero (Fixsen 2009). At lower redshifts, the gas temperature starts dropping faster as \(T_\mathrm{gas} \propto a^{-2}\) until the ionizing radiation from the first massive stars heats the IGM at \(z \approx 20\) to \(40\); at that point, the gas temperature had dropped to a few tens of Kelvin, but stars then quickly heat the gas to \(\approx 1000\,\mathrm{K}\) by \(z=10\) and all the way up to a few times \(10^4\,\mathrm{K}\) at redshifts \(z \approx 2\) when massive galaxies form (Becker et al. 2011; Rudie et al. 2012). From the discussion of radiative cooling above, it is clear that this gas is unable to cool and, thus, this is the temperature of the gas as it starts virializing inside a halo.

To determine whether the virializing gas behaves differently from the collisionless dark matter, we need to consider its hydrodynamical properties. As we discussed in Chapter 17.1.4.1, gas pressure becomes important on length scales smaller than the Jeans length \(\lambda_J\) from Equation (17.55). We can also conveniently express this criterion using an equivalent Jeans mass \(M_J\) below which pressure is important and gravitational collapse is halted. The Jeans mass is the mass within a sphere of radius \(\lambda_J/2\) and is, thus, given by \begin{align} M_J & = {\pi^{5/2}\over 6}\,\left({q+2\over q}\,{k_B \over \mu\,G\,m_p}\right)^{3/2}\,{T_\mathrm{gas}^{3/2}\over\overline{\rho}_m^{1/2}} \approx 10^{10}\,M_\odot\,\left({T_\mathrm{gas} \over 10^4\,\mathrm{K}}\right)^{3/2}\,\left({\overline{\rho}_m \over 10^{-28}\,\mathrm{g\,cm}^{-3}}\right)^{-1/2}\,,\label{eq-hierarchgal-MJ-vals} \end{align} where \(\overline{\rho}_m\) is the mean matter density at the relevant redshift. In the second equation, we have used reference values for the temperature of the gas and the mean matter density of the Universe that are approximately those at \(z=2\). Because the temperature of the intergalactic medium is lower and the mean matter density is higher at higher redshifts, the Jeans mass is smaller at earlier times. So pressure does not affect the evolution of baryonic overdensities on scales relevant for galaxy formation (except for perhaps the lowest mass galaxies). For non-cosmological gas clouds, one obtains the same expression for the Jeans mass that separates regimes in which the cloud either gravitationally collapses or is supported by pressure (Jeans 1902). This applies when the halo collapses. Then, its mean density increases by a factor equal to the virial overdensity \(\Delta_v \approx 200\), so the Jeans mass decreases by an order of magnitude. Thus, the Jeans mass is much less than the virializing halo mass of all but the smallest galaxies and pressure is therefore only immediately important at the center of a collapsing halo.

In the hot-mode accretion picture introduced by Rees & Ostriker (1977), Silk (1977), and Binney (1977), the pressurized center drives a shock through the virializing halo that fully thermalizes the bulk velocity of the gas out to the virial radius (Gull & Northover 1975; Lea 1976). The resulting virial temperature can be estimated by equating the infall kinetic energy, \(E_\mathrm{kin} = M_\mathrm{gas}\,v^2_{\mathrm{gas}} / 2\) to the final thermal energy \(E_\mathrm{therm} = 3\,M_\mathrm{gas}\,k_B\,T/(2\mu\,m_p)\), where we have assumed a mono-atomic gas. We can estimate the infall kinetic energy using the virial theorem, assuming that the gas approximately virializes before being shocked. Using the scalar virial theorem from Equation (14.18) and approximating the potential energy of the gas as that of a homogeneous density sphere with radius \(R = r_{\mathrm{vir}}\) (see Equation 14.31; note that here we need to replace \(M^2 \rightarrow M_\mathrm{gas}\,M_\mathrm{vir}\), because the gas is only a fraction of the total mass), we find that \begin{align} T_\mathrm{vir} & = {1\over 5}\,{\mu\,m_p\over k_B}\,{G\,M_\mathrm{vir}\over r_\mathrm{vir}}\,. \end{align} However, the pre-factor depends on the density structure assumed to compute the potential energy and it is common to instead use \(1/2\) as this corresponds approximately to an isothermal halo. Also substituting in the circular velocity \(v^2_\mathrm{vir} = G\,M_\mathrm{vir}/ r_\mathrm{vir}\) at the virial radius, we then find \begin{align}\label{eq-hierarchgal-tvir} T_\mathrm{vir} & = 3.6\times10^5\,\mathrm{K}\,\left({v_\mathrm{vir} \over 100\,\mathrm{km\,s}^{-1}}\right)^2\,, \end{align} where we have used \(\mu = 0.59\) as the mean molecular weight, because at these temperatures the primordial gas is fully ionized. When shock-heated to this temperature, the Jeans mass from Equation (19.16) increases by two orders of magnitude or more and pressure becomes important for the entire halo.

After being shock heated, the gaseous halo reaches hydrostatic equilibrium in which the gravitational force is balanced by pressure gradients (see Equation 17.28, where \(\vec{u} = \nabla \vec{u} = 0\) in equilibrium). However, because at this point \(T_\mathrm{vir} > 10^4\,\mathrm{K}\) for most forming halos, the radiative cooling processes discussed in the previous subsection occur and this affects the hydrostatic equilibrium. Exactly what happens depends on how fast cooling proceeds. We can define the cooling time \(t_\mathrm{cool}\) as the time it takes to remove all of the internal energy of a gas of temperature \(T\), number density \(n\), and composition \(Z\) using radiative cooling. The internal energy per unit mass is \(u=3\,k_B\,T/[2\mu\,m_p]\) (Equation 17.34) or per unit volume \(u_V = \rho\,u = 3\,n\,k_B\,T/2\). The cooling time is then \(t_\mathrm{cool} \equiv u_V/\mathcal{C}\) or \begin{align}\label{eq-hierarchgal-tcool} t_\mathrm{cool} & = {3\,n\,k_B\,T \over 2\,n_H^2\,\Lambda(T,Z)}= 0.44\,\mathrm{Gyr}\,\left({T \over 3\times 10^5\,\mathrm{K}}\right)\,\left({10^{-3}\,\mathrm{cm}^{-3}\over n_H}\right)\,\left({10^{-23}\,\mathrm{erg\,s}^{-1}\,\mathrm{cm}^3\over \Lambda[T,Z]}\right)\,, \end{align} where in the second equation we have used that \(n_H = 0.445\,n\) for a fully ionized primordial gas (\(n_H/n = x/[2x+3y]\), where \(x= 0.924\) and \(y=0.076\) are the primordial abundances of hydrogen and helium by number). The reference values of \(T\), \(n_H\), and \(\Lambda(T,Z)\) used in the second equation are typical of a \(v_{\mathrm{vir}} = 100\,\mathrm{km\,s}^{-1}\) halo collapsing at \(z \approx 2\) (including the ultraviolet background) and we see that in this case the cooling time is about half a Gyr.

To determine the effect of radiative cooling on the hydrostatic equilibrium of gaseous halos, we then need to compare \(t_\mathrm{cool}\) to two other time scales. The first is the age of the Universe \(t_0 \approx t_H = 1/H(z)\), where \(t_H\) is the Hubble time. If \(t_\mathrm{cool} \gg t_0\), radiative cooling is unimportant and the halo remains locked in hydrostatic equilibrium. However, when \(t_\mathrm{cool} \lesssim t_0\), radiative cooling causes the halo to contract as it loses pressure support. This contraction can either be slow, if the halo has time to re-establish hydrostatic equilibrium upon losing significant amounts of energy, or fast, if the halo fully loses pressure support. Which regime one finds oneself in depends on the ratio of \(t_\mathrm{cool}\) to the time it would take the halo to collapse in the absence of pressure support, known as the free-fall time \(t_\mathrm{ff}\). The free-fall time is (see Equation 17.82; \(t_c = 2\,t_\mathrm{ff}\) in that equation) \begin{equation} t_\mathrm{ff} = \sqrt{\frac{3\pi}{32\,G\,\bar{\rho}}} = {1\over\sqrt{32}}\,t_\mathrm{dyn}\,, \end{equation} where we have also expressed it in terms of the dynamical time from Equation (2.28). For the gaseous halo, the relevant \(\bar{\rho} = \Delta_v(z)\,\bar{\rho}_m(z)\) is the mean halo overdensity at the collapse redshift \(z\) (see Chapter 17.3.1). To compute \(t_\mathrm{cool}/t_\mathrm{ff}\), we can first express the cooling time from Equation (19.19) in terms of \(t_\mathrm{ff}\) as \begin{align}\label{eq-hierarchgal-tcool-in-terms-of-tff} t_\mathrm{cool} & = {16\over \pi}\,\mu\,\left({n\over n_H}\right)^2\,k_B\,G\,m_p\,\left({\Omega_m \over \Omega_b}\right)\,\left({T \over \Lambda(T,Z)}\right)\,t^2_\mathrm{ff}\,, \end{align} and then we find that \begin{align}\label{eq-hierarchgal-tcool-over-tff} {t_\mathrm{cool}\over t_\mathrm{ff}} & = 1\times \left({T \over 3\times 10^5\,\mathrm{K}}\right)\,\left({10^{-23}\,\mathrm{erg\,s}^{-1}\,\mathrm{cm}^3\over \Lambda[T,Z]}\right)\,\left({t_\mathrm{ff} \over 0.71\,\mathrm{Gyr}}\right)\,, \end{align} where we have again used \(\mu=0.59\) and \(n_H/n = 0.445\) for a fully ionized primordial gas and \(\Omega_m/\Omega_b = 6.32\) from Planck Collaboration et al. (2020b). The reference free-fall time in this equation is that at redshift \(z \approx 1.6\).

When \(t_0 > t_\mathrm{cool} > t_\mathrm{ff}\), cooling happens, but the gas is able to re-equilibrate by heating up as it contracts due to the loss of pressure support. The gas, thus, settles into a new hydrostatic equilibrium of a slightly denser, hotter halo. However, when \(t_\mathrm{cool}/ t_\mathrm{ff} \ll 1\), cooling happens so fast that the gas is unable to recover from the loss of pressure support. Instead, the gas cools quickly to \(10^4\,\mathrm{K}\), significantly lowering its Jeans mass and thus losing all pressure support. The gaseous halo then collapses in a free-fall time. Note that halos that originally start out in the \(t_0 > t_\mathrm{cool} > t_\mathrm{ff}\) regime can enter the free-fall-time-collapse regime later if they enter the region where \(t_\mathrm{cool}/ t_\mathrm{ff} < 1\) as part of their slow evolution. In such halos, the gas collapses at a delay after the formation of the halo. In either case, radiative cooling is efficient at removing internal energy, but not at removing angular momentum. As the gas contracts, it therefore spins up, eventually settling into a rotating disk.

Because halos at higher redshift are denser, they cool faster, but they also have shorter free-fall times. The dependence of the ratio \(t_\mathrm{cool}/ t_\mathrm{ff}\) on temperature and free-fall time is determined by the complex behavior of the cooling function. In Figure 19.15, we display \(t_\mathrm{cool}/ t_\mathrm{ff}\) as a function of \(T\) and \(t_\mathrm{ff}\) for gas with a primordial composition; we indicate the \(t_\mathrm{cool}/ t_\mathrm{ff}=1\) boundary using a black curve.

[28]:
# Cooling function interpolation
from scipy import interpolate
T,Lambda= read_cooling(metals='primordial')
Lambda_int= interpolate.InterpolatedUnivariateSpline(numpy.log10(T),numpy.log10(Lambda),k=3)
def tcool_over_tff(T,tff,Lambda_int=Lambda_int):
    return (T/3e5)*(10.**(-23.-Lambda_int(numpy.log10(T))))*tff/0.71
# Evaluate ratio on grid
Ts= numpy.geomspace(6e3,1e8,201)
tffs= numpy.geomspace(1e3,1e-2,201)
figure(figsize=(9,8))
galpy_plot.dens2d(numpy.log10(tcool_over_tff(Ts[:,None],tffs[None,:],
                                             Lambda_int=Lambda_int)).T,
                  origin='lower',
                  xrange=[numpy.log10(Ts[0]),numpy.log10(Ts[-1])],
                  yrange=[numpy.log10(tffs[0]),numpy.log10(tffs[-1])],
                  aspect=-(numpy.log10(Ts[-1])-numpy.log10(Ts[0]))\
                      /(numpy.log10(tffs[-1])-numpy.log10(tffs[0])),
                  xlabel=r'$T\,(\mathrm{K})$',
                  ylabel=r'$t_\mathrm{ff}\,(\mathrm{Gyr})$',
                  colorbar=True,cmap='coolwarm',shrink=1.,
                  zlabel=r'$\log_{10}\left(t_\mathrm{cool}/t_\mathrm{ff}\right)$',
                  vmin=-2.,vmax=2.,gcf=True,
                  contours=True,levels=[0.])
# Add Z=Zsolar curve
T,Lambda= read_cooling(metals=1.)
Lambda_int= interpolate.InterpolatedUnivariateSpline(numpy.log10(T),numpy.log10(Lambda),k=3)
galpy_plot.dens2d(numpy.log10(tcool_over_tff(Ts[:,None],tffs[None,:],
                                             Lambda_int=Lambda_int)).T,
                  origin='lower',
                  xrange=[numpy.log10(Ts[0]),numpy.log10(Ts[-1])],
                  yrange=[numpy.log10(tffs[0]),numpy.log10(tffs[-1])],
                  justcontours=True,overplot=True,
                  contours=True,levels=[0.],
                  cntrcolors='0.3')
from matplotlib.ticker import FuncFormatter
gca().set_xticks([4.,5.,6.,7.,8.])
gca().xaxis.set_major_formatter(FuncFormatter(\
        lambda y,pos: rf"$10^{{{y:.0f}}}$"))
gca().yaxis.set_major_formatter(FuncFormatter(\
        lambda y,pos: rf"$10^{{{y:.0f}}}$" if not y == 0 else '$1$'))
# Alternative scales
ax_orig= gca()
if True:
    # Vvir on alt x
    orig_xlim= ax_orig.get_xlim()
    ax_twiny= ax_orig.twiny()
    old_ticks= ax_orig.get_xticks()
    ax_twiny.set_xticks(old_ticks)
    def velocity_label(x,pos):
        v= numpy.sqrt(10.**x/(3.5*1e5))*100.
        # Round to nearest 10,100...
        lscale= numpy.floor(numpy.log10(v))
        v= numpy.round(v/10.**lscale)*10.**lscale
        return rf"${{{v:.0f}}}$"
    ax_twiny.xaxis.set_major_formatter(FuncFormatter(velocity_label))
    ax_twiny.set_xlim(orig_xlim) # has to be ~last
    ax_twiny.set_xlabel(r"$v_\mathrm{vir}\,(\mathrm{km\,s}^{-1})$")
# z on alt y, by hand...
# stuff to compute z from tff
from astropy.constants import G
def Deltav(redshift):
    x= Planck18.Om(redshift)-1.
    return (18.*numpy.pi**2.+82.*x-39*x**2.)/Planck18.Om(redshift)
zs= numpy.linspace(0.,20.,101)[::-1]
tff_int= interpolate.InterpolatedUnivariateSpline(\
    numpy.log10(((3.*numpy.pi/G/32/(Planck18.critical_density(zs)*Planck18.Om(zs)*Deltav(zs)))**0.5).to_value(u.Gyr)),
    zs,k=3)
sca(ax_orig)
for z in [0,1,2,4,6,8,10,15]:
    if z == 15:
        annotate(rf'$z$',(numpy.log10(Ts[-1])-0.1,
                        numpy.log10(((3.*numpy.pi/G/32/(Planck18.critical_density(z)*Planck18.Om(z)*Deltav(z)))**0.5).to_value(u.Gyr))),
                 size=17.,annotation_clip=False)
    else:
        annotate(rf'${z}$',(numpy.log10(Ts[-1])-0.1,
                        numpy.log10(((3.*numpy.pi/G/32/(Planck18.critical_density(z)*Planck18.Om(z)*Deltav(z)))**0.5).to_value(u.Gyr))),
                 size=17.,annotation_clip=False)
t= text(4.2,-1.6,r'$\mathrm{cooling\ effective}$',
        color='k',fontsize=20.)
t.set_bbox(dict(facecolor='w',alpha=0.75,edgecolor='none'))
t= text(6.25,2.75,r'$\mathrm{cooling\ ineffective}$',
        color='k',fontsize=20.)
t.set_bbox(dict(facecolor='w',alpha=0.75,edgecolor='none'))
text(5.1,-0.3,r'$Z=0$',fontsize=18.)
text(5.5,0.3,r'$Z=Z_\odot$',fontsize=18.);
../_images/chapters_IV-03.-Hierarchical-Galaxy-Formation_118_0.svg

Figure 19.15: The efficiency of cooling in galactic gaseous halos.

To place this figure in the galaxy-formation context, we have also converted \(T\) to \(v_\mathrm{vir}\) using Equation (19.18) (top \(x\) axis) and \(t_\mathrm{ff}\) to the redshift at which this is the free-fall time within the virial radius (on the right \(y\) axis). Thus, if the gaseous halo is shock-heated to \(T_\mathrm{vir}\), this figure shows whether or not the gas can efficiently cool. To demonstrate the effect of increased cooling of enriched gas, we also show the \(t_\mathrm{cool}/ t_\mathrm{ff}=1\) boundary computed using the \(Z=Z_\odot\) cooling function (all \(\Lambda[T,Z]\) here ignore the ultraviolet background; see below). Real halos have density gradients and Figure 19.15 demonstrates that cooling at any temperature is more efficient at the higher densities (shorter \(t_\mathrm{ff}\)) at the centers of halos.

Thus, we see that shock-heated small galaxies with \(v_\mathrm{vir}\) down to \(20\,\mathrm{km\,s}^{-1}\) can efficiently cool at all redshifts regardless of the composition of the gas in the absence of the ultraviolet background; however, the cooling cut-off at \(T \approx 10^4\,\mathrm{K}\) means that galaxies with \(v_\mathrm{vir} \lesssim 20\,\mathrm{km\,s}^{-1}\) cannot efficiently cool, unless they have significant amounts of molecular hydrogen, which allows cooling at \(T < 10^4\,\mathrm{K}\). Large galaxies with \(v_\mathrm{vir} \gtrsim 100\,\mathrm{km\,s}^{-1}\) cannot cool efficiently at the low redshifts at which they typically collapse unless their gas is enriched. Gas in clusters of galaxies with \(v_\mathrm{vir} \approx 1000\,\mathrm{km\,s}^{-1}\) cannot efficiently cool at any redshift, even if it is enriched. The implications of these behaviors for galaxy formation, however, depend on whether the gaseous halo was shock-heated to the virial temperature in the first place. As we will discuss below, gas in low-mass halos is likely not shock-heated, but rather flows efficiently to the center where is then heated and cools at higher density. Gas in clusters, however, is shock heated, so the discussion here implies that this gas cannot cool. And indeed, clusters contain hot, \(T \approx 10^7\,\mathrm{K}\) gas that can be observed through X-ray emission (e.g., Forman et al. 1979).

Combining the hot-mode accretion picture with the hierarchical merging and accretion behavior in the \(\Lambda\)CDM Universe in which we live, leads to the over-cooling problem (White & Rees 1978). In the hierarchical merging picture, the first objects to form are small halos that merge together at later times to form larger halos. Because a significant amount of mass in large halos was at some earlier redshift part of a \(20\,\mathrm{km\,s}^{-1} \lesssim v_{\mathrm{vir}} \lesssim 50\,\mathrm{km\,s}^{-1}\) halo in which gas cooling is highly efficient in the cooling picture discussed in the previous paragraph, much of the mass cools effectively and should form stars. Thus, we expect that much of the gas in the Universe is transformed into stars, contrary to what is observed (see the discussion in the paragraph following Equation 18.19 in Chapter 18.2) . Indeed, determinations of the stellar-mass–halo-mass relation (see Figure 18.16) show that small galaxies are instead highly inefficient at turning gas into stars. One important ingredient in solving this over-cooling problem is the ultraviolet background, which we have so far ignored in the discussion in this section. As we saw in the previous sub-section, photo-ionization caused by the ultraviolet background strongly reduces the efficiency of the collisional-excitation cooling that dominates at \(10^4\,\mathrm{K} < T < 10^6\,\mathrm{K}\), in particular almost entirely halting cooling at \(10^4\,\mathrm{K} < T < 10^5\,\mathrm{K}\) at \(z\lesssim 10\). Thus, in the presence of the ultraviolet background, small galaxies with \(v_\mathrm{vir} \lesssim 50\,\mathrm{km\,s}^{-1}\) are unable to cool their gas and they are only able to do so and form stars at high redshift when the effect of the ultraviolet background is smaller. Observed stellar populations in Milky Way dwarf galaxies are old, consistent with this picture (e.g., Weisz et al. 2014).

Until the early 2000s, the hot-mode accretion picture was the accepted formation scenario for all galaxies, even though early cosmological, hydrodynamical simulations that included radiative cooling demonstrated that the gas in halos does not necessarily shock heat to the virial temperature (e.g., Fardal et al. 2001). Birnboim & Dekel (2003) investigated the conditions under which a shock is able to propagate in the presence of radiative cooling using analytic methods and simple simulations and concluded that only halos with mass \(M_\mathrm{halo} \gtrsim 10^{11}\,M_\odot\) are able to sustain enough pressure behind the shock to drive it all the way to the virial radius. In lower mass halos, radiative cooling is too efficient in the shock-heated gas: the gas behind the shock front quickly loses pressure support and is then unable to drive the shock to the virial radius. With a more realistic setup, Dekel & Birnboim (2006) derive a critical halo mass separating the hot- and cold-mode regimes of \(M_\mathrm{halo} \approx 6\times 10^{11}\,M_\odot\). High-resolution numerical simulations confirm this critical mass (e.g., Keres et al. 2005; Dekel & Birnboim 2006). Numerical simulations also show that in the few galaxies with \(M_\mathrm{halo} \gtrsim 6\times 10^{11}\,M_\odot\) that form at \(z \gtrsim 2\), some cold gas flows are able to penetrate deep inside the halo even in the presence of a virial shock that heats most of the gaseous halo.

Gas accreted in the cold mode reaches the inner halo quickly without ever being significantly heated, reaching the same high density characterizing galactic-disk gas that is the end-point of hot-mode accretion. Thus, star formation in cold-mode halos (\(M_\mathrm{halo} \lesssim 6\times 10^{11}\,M_\odot\)) commences more efficiently and is more likely to be bursty, affecting the subsequent evolution of the galaxy through, for example, the amount of feedback from supernovae that bursty star formation generates. Gas accreted later in the cold mode similarly flows to the central galaxy quickly. However, the effect of photo-ionization by the ultraviolet background remains, heating the gas in low-mass halos, thus still preventing it from cooling and forming stars. The critical halo mass separating cold- and hot-mode accretion is close to the mass at which the stellar-mass–halo-mass relation (SMHMR) peaks (see Figure 18.16) , indicating that the transition to inefficient hot-mode accretion is part of the reason for the high-mass turn-over in the SMHMR.

In this section, we have focused on the question of how gas from the intergalactic medium forms into or accretes onto a disk. Once the gas settles into a disk, additional cooling and fragmentation into molecular clouds and clumps within those clouds are necessary for stars to form. We have discussed star formation in disks of galaxies in some detail in Chapter 18.5 and we refer the reader to specialized monographs for more details on the star-formation process (e.g, Krumholz 2017). We have also concentrated on what happens to gas as a halo first collapses or as it is accreted into a halo and flows to the center. However, when gas accretes onto an existing galaxy, other heating processes occur, such as thermal and mechanical feedback from supernova explosions or active galactic nuclei, and these affect the extent to which accreted gas is able to cool into a star-forming disk. We discuss these in the next subsection.

19.3.3. Feedback and galaxy formation

\label{chapter-merger-gas-feedback}

In the previous subsection, we saw that gas in galaxies with \(v_\mathrm{vir} \lesssim 50\,\mathrm{km\,s}^{-1}\) is unable to cool because of the ultraviolet background and that the shock-heated gas in groups and clusters of galaxies with \(v_\mathrm{vir} \gtrsim 500\,\mathrm{km\,s}^{-1}\) remains hot because of the inefficiency of radiative cooling at high virial temperatures. But galaxies in between these two extremes, including all large galaxies and a large fraction of the observed galaxy population as a whole, should be able to cool mildly-enriched gas efficiently in both the cold and hot modes of accretion (e.g., Figure 19.15). Thus, galaxy formation should be highly efficient for most galaxies, turning all of the accreted gas into stars in a relatively short time. However, galaxy formation is observed to be highly inefficient. For example, determinations of the stellar-mass–halo-mass relation (SMHMR; see Chapter 18.2.3 and in particular Figure 18.16) show that the fraction of a galaxy’s mass contained in stars is much less than the universal baryon fraction, \(M_*/M_{\mathrm{halo}}\ll \Omega_b/\Omega_m\), at the \(M_{\mathrm{halo}} \approx 10^{12}\,M_\odot\) peak of the SMHMR and that \(M_*/M_{\mathrm{halo}}\) becomes another order of magnitude less at lower and higher masses. The Kennicutt-Schmidt law disussed in Chapter 18.5 demonstrates that only a fraction of the available gas in galaxy disks is converted to stars (e.g., \(\approx 10\%\) of the gas is turned into stars in a dynamical time; see Figure 18.27). At the molecular-cloud level, star formation appears to be highly inefficient, turning only about 1% of the molecular gas into stars in a free-fall time (Krumholz et al. 2012). These observations demonstrate that galaxy formation is inefficient at all scales: getting the gas from the halo into the disk and keeping it there, getting the disk gas to form star-forming clouds, and forming stars from dense molecular gas in the short period before the clouds are disrupted, are all difficult to achieve. The main reason behind all of these inefficiencies is currently believed to be various types of energetic feedback processes (although magnetic-field pressure may provide additional support against collapse; e.g., Parker 1969).

There are two main sources of feedback in galaxies: massive stars and active galactic nuclei (AGN; see Chapter 18.3). Massive stars primarily affect their surrounding interstellar medium (ISM) when they explode as supernovae, but intense radiation from massive stars before they explode can also have significant effects. For AGN, radiative feedback from the central region as well as feedback from the highly-energetic jets launched by the AGN that reach far into the host galaxy’s halo are both important. For both massive stars and AGN, feedback can be either preventive (keeping gas in the hot halo rather than cooling into the cold ISM, similar to how the ultraviolet background does this) or ejective (removing gas from the cold, star-forming ISM); while both types likely occur to some degree, ejective feedback dominates in typical star-forming galaxies (Somerville & Davé 2015).

When massive stars explode as Type II supernovae, they inject energy and momentum into the surrounding ISM in a process known as supernova feedback. The effect of the large amount of energy produced in supernova explosions (\(\approx 10^{51}\,\mathrm{erg}\)) was early on considered to be potentially important for the structure of galaxies (Larson 1974; Dekel & Silk 1986; Babul & Rees 1992), but the small scales on which supernova explosions occur, and where they affect their immediate surroundings, have been (and remain!) hard to resolve in numerical simulations. Therefore, investigations necessarily need to make assumptions about the structure of the ISM near the supernovae and how the explosion interacts with this structure (Efstathiou 2000). The specific assumptions made significantly affect the outcome (Navarro & White 1993); this remains an issue in current simulations. One problem is that stars form in dense environments where cooling times are short, so energy added to the local ISM through explosive heating is efficiently radiated away, preventing any significant feedback from occurring (Katz 1992). Direct momentum transfer—pushing rather than heating—between the explosion and the ISM is therefore likely the main mechanism through which supernova feedback happens.

Supernovae are, however, not the only source of stellar feedback. Radiation pressure, jets (Murray et al. 2010), and radiative heating (Murray et al. 2011) from massive stars can disrupt molecular clouds before their first supernova explosion and drive enriched winds deep into the halo. Simulations incorporating such early-stellar feedback demonstrate that it provides an important contribution to the overall regulation of star formation, in particular preventing more star formation at \(z \gtrsim 1\) than supernova feedback alone would and leading to better agreement with the evolution of the SMHMR at high redshift (Hopkins et al. 2011; Stinson et al. 2013; Hopkins et al. 2014). Cosmic rays injected in supernova and stellar-wind shocks can also significantly suppress star formation in higher-mass galaxies at late times (\(M_{\mathrm{halo}} \gtrsim 10^{11}\,M_\odot\) and \(z \lesssim 2\); Hopkins et al. 2020).

The combination of supernova and early-stellar feedback causes the SMHMR to decline at halo masses below its \(M_{\mathrm{halo}} \approx 10^{12}\,M_\odot\) peak: the smaller the galaxy, the more easily feedback processes can remove star-forming gas, leading to a declining \(M_*/M_{\mathrm{halo}}\) towards small halo masses. As discussed in Chapter 6.8, supernova feedback also affects the inner dark-matter structure of small galaxies, turning the inner cusp into a core through the non-adiabatic potential changes induced by the outflows.

As discussed in Chapters 16.1 and 18.3, various lines of evidence indicate that all massive galaxies contain a supermassive black hole (SMBH) that likely grows during brief epochs of AGN activity. The energetic radiation and jets emanating from the AGN pushes and heats its host galaxy’s gas similar to supernova explosions (Babul & White 1991; Silk & Rees 1998). As for stellar feedback, one main question is whether this AGN feedback is primarily momentum- or energy (heating)-based, with direct momentum transfer again likely playing an outsize role (e.g., King 2003; Murray et al. 2005). Because AGN activity manifests itself in different ways, its feedback effects also vary and two main modes are usually distinguished: quasar mode feedback for feedback arising from the radiation produced within the AGN (which can extend all the way to kpc scales) and radio mode feedback for feedback arising from the radio jets launched by the AGN. Because radiation in the quasar mode originates in the high-density center where cooling times are short, direct momentum transfer must again dominate feedback. Jets—whose kinetic energy can significantly exceed the total luminosity of the AGN—can penetrate far into the lower-density hot halo and significantly heat the gas and heating therefore plays a significant role in radio mode feedback (McNamara & Nulsen 2007).

The suppression of star formation by AGN feedback is largest for galaxies with masses above the SMHMR \(M_{\mathrm{halo}} \approx 10^{12}\,M_\odot\) peak mass, where there is not enough energy in supernova explosions to sufficiently suppress star formation (Benson et al. 2003); AGN feedback is therefore the likely culprit for the high-mass cut-off in the galaxy stellar mass function (see, e.g., Equation 18.28). Because AGN feedback pushes away and heats the gas that could fuel the growth of the central SMBH, AGN feedback is also important in regulating the growth of the SMBH and establishing the M–sigma relation from Equation (16.6) (see, e.g., Di Matteo et al. 2005; Sijacki et al. 2007; Somerville et al. 2008).

While feedback produces \(\mathcal{O}(1)\) effects on galaxies, directly observing feedback in action is difficult. One environment where feedback can be somewhat directly constrained is the circumgalactic medium, the general name for the reservoir of gas inside of a galaxy that is outside of the disk (and what we have more informally referred to as the gaseous halo here). The circumgalactic medium contains diffuse halo gas, cold filaments of accreting gas for galaxies in the cold-mode accretion regime, and gas expelled from the disk through feedback processes (Tumlinson et al. 2017). Galactic winds or outflows can be directly observed (Veilleux et al. 2005), either by directly staring down the barrel of the outflow and observing blueshifted metal absorption in a galaxy’s own spectrum (e.g., Martin 2005; Weiner et al. 2009; Steidel et al. 2010) or by looking for absorption by gas in a galaxy’s halo in the spectrum of a background quasar (e.g., Chen et al. 1998; Tripp et al. 2000; Stocke et al. 2006; Lehner et al. 2015). Biconical outflows coming from the center are also directly seen in the Milky Way (Fox et al. 2015). Such observations demonstrate that outflows are ubiquitous in star-forming galaxies and that they extend out to 100 kpc. Metal absorption in the outflows allows one to determine that the outflows are enriched in metals (e.g., Prochaska et al. 2017) and that they are, thus, important in removing metals from the ISM. This removal is required by the observed mass–metallicity relation (see the discussion of leaky-box models in Chapter 11; Tremonti et al. 2004). A direct accounting of the amount of metals contained in stars and the ISM also falls far short of the total amount of metals that have been produced as a by-product of cosmic star formation (Peeples et al. 2014), and the CGM contains a significant fraction of these “missing metals”.

As an essential ingredient of galaxy formation, stellar and AGN feedback is included in all modern galaxy formation simulations. The small scales involved in star formation and in the interaction of stellar feedback and the surrounding ISM, and the poorly-understood physical mechanisms at play in AGN feedback, however, means that feedback has to be approximated using physically-informed prescriptions on the scales that cannot be resolved and/or fully understood. Such assumptions are referred to as sub-grid physics. The specific sub-grid prescriptions used generally depend on the resolution of the simulation, with lower resolution simulations often having to make quite unrealistic and ad hoc assumptions (such as artificially turning off cooling after a supernova explosion to avoid that the explosion’s heating is simply radiated away; Gerritsen & Icke 1997) and high-resolution simulations typically needing to make fewer assumptions that are informed by more controlled, non-cosmological, high-resolution simulations (e.g., Hopkins et al. 2018). The detailed physics and the various effects of feedback on galaxies remains one of the biggest open problems in galaxy formation.