IQ Collaboratory. II. The Quiescent Fraction of Isolated, Low-mass Galaxies across Simulations and Observations

, , , , , , , , , , , , and

Published 2021 July 5 © 2021. The American Astronomical Society. All rights reserved.
, , Citation Claire M. Dickey et al 2021 ApJ 915 53 DOI 10.3847/1538-4357/abc014

Download Article PDF
DownloadArticle ePub

You need an eReader or compatible software to experience the benefits of the ePub3 file format.

0004-637X/915/1/53

Abstract

We compare three major large-scale hydrodynamical galaxy simulations (EAGLE, Illustris-TNG, and SIMBA) by forward modeling simulated galaxies into observational space and computing the fraction of isolated and quiescent low-mass galaxies as a function of stellar mass. Using SDSS as our observational template, we create mock surveys and synthetic spectroscopic and photometric observations of each simulation, adding realistic noise and observational limits. All three simulations show a decrease in the number of quiescent, isolated galaxies in the mass range M* = 109−10 M, in broad agreement with observations. However, even after accounting for observational and selection biases, none of the simulations reproduce the observed absence of quiescent field galaxies below M* = 109 M. We find that the low-mass quiescent populations selected via synthetic observations have consistent quenching timescales, despite an apparent variation in the late-time star formation histories. The effect of increased numerical resolution is not uniform across simulations and cannot fully mitigate the differences between the simulations and the observations. The framework presented here demonstrates a path toward more robust and accurate comparisons between theoretical simulations and galaxy survey observations, while the quenching threshold serves as a sensitive probe of feedback implementations.

Export citation and abstract BibTeX RIS

1. Introduction

The transition of galaxies from star-forming into quiescence provides key insights into the physical processes that drive galaxy evolution. Major extragalactic surveys like the Sloan Digital Sky Survey (SDSS; York et al. 2000) have been instrumental in quantifying the differing properties of galaxies in these two categories, including morphology, environment, and color (e.g., Blanton et al. 2003; Kauffmann et al. 2003; Brinchmann et al. 2004; Blanton & Moustakas 2009; Moustakas et al. 2013). The processes that drive a galaxy's transformation can be broadly divided into two categories: external processes, which occur in high-density environments, and processes internal to the galaxy, which are often thought to be correlated to galaxy mass (e.g., Peng et al. 2010).

Understanding how and when feedback mechanisms operate within galaxies and in which regimes they are most effective remains a major focus of both observational and theoretical studies of galaxy evolution. Observational galaxy surveys have been used to characterize the populations of quiescent and star-forming galaxies as a function of stellar mass (Kauffmann et al. 2004; Geha et al. 2012), environment (Wetzel et al. 2012), and redshift (Brammer et al. 2009), as well as key correlations between galaxy properties (e.g., the star-forming sequence; Daddi et al. 2007; Elbaz et al. 2007; Noeske et al. 2007; Salim et al. 2007) that inform our understanding of the underlying forces that drive galaxy evolution.

Large-scale cosmological hydrodynamic galaxy formation simulations reproduce the observed universe in a qualitative manner (Genel et al. 2014; Vogelsberger et al. 2014; Schaye et al. 2015; Davé et al. 2017, 2019; Pillepich et al. 2018a). Each simulation uses a distinct set of approximations for the complex physics of underlying galaxy evolution, including star formation, heating and cooling of gas, black hole formation and growth, feedback from active galactic nuclei (AGNs), and stellar feedback (for an overview, see Somerville & Davé 2015; Vogelsberger et al. 2020). While it has been shown that many different subgrid models can produce relatively consistent pictures of galaxy evolution (Naab & Ostriker 2017), the variations between simulations present a unique opportunity to explore how feedback mechanisms shape galaxy evolution.

Many studies have focused on comparing the simulated distributions of quantities such as galaxy masses, colors, or star formation rates (SFRs) to observations (e.g., Genel et al. 2014; Torrey et al. 2014; Vogelsberger et al. 2014; Schaye et al. 2015; Somerville & Davé 2015; Sparre et al. 2015; Trayford et al. 2017; Davé et al. 2017; Nelson et al. 2018; Pillepich et al. 2018a; Donnari et al. 2019, 2020, 2021). These works primarily focus on comparing a singular simulation to a particular set of observations. Studies using multiple simulations require the careful construction of a consistent framework for intersimulation comparison.

In the observed universe, the quiescent fraction for isolated galaxies is zero at M* < 109 M (Geha et al. 2012) in the SDSS volume, suggesting that feedback in this mass regime is highly sensitive to stellar mass. Most studies of feedback in low-mass galaxies have focused on environmental quenching in the context of either high-density environments such as clusters or interactions with a Milky Way–like host galaxy (e.g., Wetzel et al. 2013, 2015; Sales et al. 2015; Fillingham et al. 2018).

Internal feedback may influence galaxy evolution down to M* ∼ 109 M via either AGN (Penny et al. 2018; Dickey et al. 2019; Koudmani et al. 2019) or stellar (El-Badry et al. 2016) feedback. However, the efficiency and interplay of different feedback mechanisms below this mass limit remain uncertain. Large-scale hydrodynamic simulations present a unique opportunity to create an "observed" magnitude-limited galaxy survey and explore how different subgrid feedback implementations shape the distribution.

The Isolated & Quiescent (IQ) Collaboratory 20 aims to bridge the gap between simulations and observations of star-forming and quiescent galaxies to better characterize internal quenching processes. In Hahn et al. (2019, hereafter Paper I), we began the process of comparing the star-forming sequence across a set of simulations.

Our initial analysis in Paper I highlighted the large variation in the apparent quiescent populations of each simulation, but we did not make a direct comparison between simulations and observations. In that work, the SFRs from the simulations are obtained "directly" from the simulation output, as opposed to being derived from "observables" like Hα- or UV+IR-derived SFRs or indices such as the Hα equivalent width (EW) and Dn 4000, an index that measures the strength of the 4000 Å break and loosely traces stellar age. Moreover, even within the observational results, selecting different tracers for, e.g., SFRs can give rise to significant variation (e.g., Salim et al. 2007; Kennicutt & Evans 2012; Speagle et al. 2014; Flores Velázquez et al. 2021; Katsianis et al. 2020).

In addition to deriving galaxy observables from a simulation, fully transforming a simulated volume into a true observational analog requires the careful creation of complete, volume-corrected samples. This is particularly challenging when we are trying to reproduce observations comparing star-forming and quiescent populations, which often suffer from differing levels of incompleteness within the same survey.

Our goal in this work is to create an "apples-to-apples" comparison of the quiescent fraction of isolated low-mass (M* < 1010 M) galaxies as a function of stellar mass using observations of the local universe and a set of simulations, all of which produce realistic cosmological volumes using distinct implementations of subgrid physics. We focus our study on a selection of three major cosmological hydrodynamic simulations (EAGLE, Illustris-TNG, and SIMBA) and compare them to the observed population of low-mass quiescent galaxies in the SDSS.

In Section 2, we present the simulations used in this work, and in Section 3, we discuss the observations that serve as our point of comparison. In Section 4, we create mock observations from each simulation, including synthetic spectra, realistic noise, and SDSS-like incompleteness. Section 5 compares the "observed" simulations to the distribution of low-mass quiescent galaxies in SDSS. In Section 7, we explore the star formation histories of the quiescent galaxies in the simulations, and we review our findings in Section 8.

2. Galaxy Formation Simulations

In this work, we focus on three large-scale hydrodynamic cosmological simulations (EAGLE, Illustris-TNG, and SIMBA), which we will compare to observations. We provide brief descriptions of each simulation below.

2.1. EAGLE

The Virgo Consortium's Evolution and Assembly of GaLaxies and their Environment (EAGLE) project 21 (Crain et al. 2015; Schaye et al. 2015; McAlpine et al. 2016) has a volume of (100 Mpc)3 (comoving), with dark matter and baryonic particle masses of 9.6 × 106 and 1.8 × 106 M. The simulation uses ANARCHY, a modified version of the Gadget3 N-body SPH code (Springel et al. 2001b; Springel 2005; Schaller et al. 2015). The subgrid model for feedback from massive stars and AGNs is based on thermal energy injection in the interstellar medium (Dalla Vecchia & Schaye 2012). The subgrid parameters were calibrated to reproduce the z = 0 stellar mass function and galaxy sizes.

Previous works that have reproduced observables and examined the quenched population of the EAGLE simulation include Schaye et al. (2015), Furlong et al. (2015), Trayford et al. (2015), and Trayford et al. (2017), and the galaxy–black hole relations are discussed in McAlpine et al. (2017, 2018), Bower et al. (2017), and Habouzit et al. (2021). Trčka et al. (2020) showed that mock spectral energy distributions (SEDs; see Camps et al. 2018 for the public release of the SEDs) based on EAGLE galaxies and radiative transfer calculations using skirt (Baes et al. 2011) are in overall agreement with galaxy spectra from the observed sample and highlight the need for consistent comparisons between simulations and observations. Trayford et al. (2017) used the same results of the skirt radiative transfer code to determine quenched fractions based on mock UVJ colors, along with the distribution of Hα flux and Dn 4000. They found that the passive fraction varies significantly depending on the definition of quenched. We will compare our results to previous conclusions in Section 6.1.1.

2.2. Illustris-TNG

The Next Generation Illustris project (Illustris-TNG or TNG 22 ; Marinacci et al. 2018; Naiman et al. 2018; Nelson et al. 2018; Pillepich et al. 2018a; Springel et al. 2018) is the successor to the original Illustris project (Genel et al. 2014; Vogelsberger et al. 2014), with significant updates in the subgrid models and physics included in the simulation. We use the TNG100 simulation, which has a volume of (110.7 Mpc)3, and dark matter and baryonic mass resolutions of 7.6 × 106 and 1.4 × 106 M. The TNG is run using the AREPO moving-mesh code (Springel 2010), which is based on the Gadget code (Springel et al. 2001b; Springel 2005). Adjustments in the TNG model include the addition of magnetohydrodynamics, updated stellar feedback prescriptions, and the transition from thermal "bubbles" in the intergalactic medium to a black hole–driven kinetic wind for the low accretion rate black hole feedback mode (Weinberger et al. 2017; Pillepich et al. 2018b). The color bimodality of low-redshift galaxies is shown to compare well to SDSS (Nelson et al. 2018). Other papers describe the size evolution of quiescent galaxies (Genel et al. 2018) and correlations between galaxy properties and supermassive black holes and AGN feedback (Weinberger et al. 2018; Habouzit et al. 2019, 2021; Li et al. 2020; Davies et al. 2020; Terrazas et al. 2020). Donnari et al. (2019) showed quenched fractions based on UVJ selection and the star-forming sequence and compared these to UVJ-selected observed quenched fractions from COSMOS/UltraVISTA (Muzzin et al. 2013). Donnari et al. (2020, 2021) explored the quenched fraction derived using SFRs and distances from the galaxy star-forming sequence. Donnari et al. (2021) in particular explored how systematic uncertainties affect a single set of observations, which we will discuss further in Section 6.1.2.

2.3. SIMBA

SIMBA (Davé et al. 2019) is a suite of cosmological simulations built on the GIZMO meshless finite-mass hydrodynamics code (Hopkins 2015, 2017), also based on the Gadget code (Springel et al. 2001b; Springel 2005), and forms the next generation of the MUFASA (Davé et al. 2016) simulations with novel black hole growth and feedback subgrid models. The fiducial run has a volume of (143 Mpc)3 and dark matter and baryonic mass resolutions of 9.6 × 107 and 1.8 × 107 M, respectively.

SIMBA includes a model for on-the-fly dust production and destruction (broadly following McKinnon et al. 2017), and star formation is regulated with two-phase kinetic outflows, which were tuned to predictions from the Feedback in Realistic Environments (FIRE) simulations (Hopkins et al. 2014; Muratov et al. 2015; Anglés-Alcázar et al. 2017b; Hopkins et al. 2018). Black hole growth in SIMBA is based on the torque-limited accretion model (Hopkins & Quataert 2011; Anglés-Alcázar et al. 2013, 2015, 2017a) linking the black hole accretion to the properties of the galaxy's inner gas disk. The AGN feedback consists of kinetic bipolar outflows, modeled after observed outflows of AGNs, and X-ray feedback input in the surrounding gas similar to Choi et al. (2012).

Previous work using SIMBA has already discussed the radial density profiles of quenched galaxies (Appleby et al. 2020), the weak correlation between galaxy mergers and quenching (Rodríguez Montero et al. 2019), and the connection between quenching and black hole growth (Davé et al. 2019; Thomas et al. 2019; Habouzit et al. 2021). We compare our results to relevant studies in Section 6.1.3.

3. Observations

Following Geha et al. (2012), we build our sample of observed galaxies from SDSS using the NASA/Sloan Atlas (NSA; Blanton et al. 2011), 23 which includes all galaxies within z < 0.055 in the SDSS footprint. The NSA is a rereduction of SDSS DR8 (Aihara et al. 2011) optimized for nearby low-luminosity objects (Yan 2011; Yan & Blanton 2012), including an improved background subtraction technique (Blanton et al. 2011) and the addition of near- and far-UV photometry from GALEX. The catalog includes emission line fluxes and EWs for all galaxies. We calculate stellar masses as in Mao et al. (2021), using the Bell et al. (2003) relation to determine the mass-to-light ratio from the g − r color. See Section 5.3 for a more in-depth discussion on the effects of uncertainty in the stellar mass measurements on the resultant quiescent fraction.

As in Geha et al. (2012), we consider a low-mass galaxy to be quiescent if it has little to no star formation and is dominated by older stellar populations. For the first condition, we use the Hα EW, which traces recent (<10 Myr) star formation, and require Hα EW < 2 Å. To probe galaxy age, we rely on Dn 4000 (Balogh et al. 1999), an index that quantifies the strength of the 4000 Å break in the spectrum and traces the light-weighted age of the stellar population. The Dn 4000 index is an indirect measure of intermediate ( ∼1 Gyr) star formation (Hamilton 1985; Brinchmann et al. 2004; Moustakas et al. 2006). We use the empirical relation of Geha et al. (2012), Dn 4000 $\gt 0.6+0.1{\mathrm{log}}_{10}({M}_{* }/{M}_{\odot })$, to select galaxies with older stellar populations.

To quantify the environment of each low-mass galaxy, we use dhost, which is defined as the projected distance to the nearest massive neighbor (M* > 2.5 × 1010 M${M}_{{K}_{s}}\lt -23;$ hereafter the host galaxy). We search for potential host galaxies within 7 Mpc and 1000 km s−1 of each low-mass galaxy using the Two Micron All Sky Survey (2MASS) Extended Source Catalog to ensure we do not miss any potential host galaxies that lie outside the SDSS imaging footprint. In the few cases where we do not identify a host within 7 Mpc, we use dhost = 7 Mpc.

We define galaxies as isolated when dhost   >  1.5 Mpc. This is an empirical choice based on the behavior of the quiescent fraction as a function of dhost (see Figure 4 of Geha et al. 2012 or Figure 5 in this work). The quiescent fraction only shows a dependence on environment for dhost < 1.5 Mpc, so we consider low-mass galaxies beyond this threshold to be isolated. Using a less strict definition (e.g., galaxies are isolated beyond dhost  =  1 Mpc) shifts the quenching threshold to a slightly lower mass, while a more conservative definition does not change the observed threshold. For galaxies with stellar masses above M* = 1010 M, we use the central satellite designation from the group catalog of Tinker et al. (2011).

In the left column of Figure 1, we show the distribution of isolated galaxies observed in the local universe as a function of Dn 4000 and SDSS r magnitude in the mass range M* = 108–10 M. Galaxies in isolation are only observed to be quiescent above M* = 109 M (Figure 2).

Figure 1.

Figure 1. Top row: distribution of absolute r magnitude (Mr ) as a function of Dn 4000 from SDSS for observations and derived from the noiseless synthetic spectra from each simulation (from left to right: TNG, EAGLE, and SIMBA) for isolated galaxies (based on dhost) in the stellar mass range M* = 108–10 M. Galaxies are color-coded as star-forming (blue) or quiescent (orange) based on the noise-free Hα EW and Dn 4000. Bottom row: distribution of apparent r magnitude (mr ) as a function of Dn 4000, derived from the noise-added synthetic spectra along a single random sight line through each simulation box. Galaxies are color-coded as star-forming (blue) or quiescent (orange) based on the noise-added Hα EW and Dn 4000. Galaxies below the gray dashed line fall below the SDSS spectroscopic limit and would not be selected for spectroscopic follow-up in SDSS. As such, they are not included in the calculation of fq along a given sight line.

Standard image High-resolution image
Figure 2.

Figure 2. The Dn 4000 vs. Hα EW from SDSS for observations and derived from the noise-added synthetic spectra from each simulation (from left to right: TNG, EAGLE, and SIMBA) for isolated galaxies (based on dhost) in three stellar mass bins (from top to bottom: M* = 108.5−9.0, 109.0–9.5, and 109.5–10.0 M). Galaxies are observed as quiescent if they fall into the orange region of Dn 4000–Hα EW parameter space. The Dn 4000 limit depends on stellar mass and is equal to the solid (dashed) vertical line at the low-mass (high-mass) end of the bin for each row. For the simulations, galaxies are color-coded by their instantaneous SFRs, and galaxies with SFR = 0 M yr−1 are shown in white.

Standard image High-resolution image

4. Mock Observations

To create a true apples-to-apples comparison between the observations and simulations, we must account for the effects of observational incompleteness, finite signal-to-noise ratio (S/N), and a lack of precise distance information for individual galaxies in a wide-field survey.

To that end, we create mock surveys and synthetic observations of each simulation, adding realistic noise and observational limits. We use SDSS as our observational template. We select the mock observational limits and injected noise to match SDSS. We apply the same methodology to all simulations, as described below for a single simulation box. However, this method can also be generalized and applied to semianalytic models and zoom-in simulations, as well as adapted to match other surveys and observations.

4.1. Synthetic Spectra

We begin by generating synthetic spectra for all galaxies in a given simulation using Flexible Stellar Population Synthesis (FSPS; Conroy et al. 2009; Conroy & Gunn 2010) and the Python interface python-fsps (Foreman-Mackey et al. 2014). Our method is described in full in T. K. Starkenburg et. al (2020, in preparation), but in brief, it is as follows.

For each galaxy, we bin the total stellar mass formed by formation age (t) and stellar metallicity (Z). We use a uniform t, Z grid across all simulations to standardize the effects of otherwise variable time resolution and particle size. Age bin size increases linearly as a function of look-back time based on the minimum age steps of the underlying SSP models from t = 0 to 13.75 Gyr. The metallicity grid spans ${\mathrm{log}}_{10}(Z/{Z}_{\odot })=-2.2$ to 0.5 with ΔZ = 0.3 at low metallicity, with bin resolution increasing to ΔZ = 0.1 near Z = 0. We use a Chabrier (2003) initial mass function throughout and include the asymptotic giant branch dust emission model of Villaume et al. (2015). To cover the most recent star formation and avoid the resolution effects of stellar particles, we set the SFRs younger than 15 Myr equal to the instantaneous SFR from the galaxy's gas particles with metallicity equal to the current mass-weighted metallicity of the star-forming gas.

We generate the spectrum of a simple stellar population at each point in the t, Z grid. The sum of the single stellar population spectra, weighted by the stellar mass formed, produces the galaxy spectrum. The nebular emission lines are calculated independently for each stellar metallicity bin using the FSPS nebular emission line prescription based on a CLOUDY lookup table (Byler et al. 2017), with the gas metallicity fixed to that of the metallicity bin.

Dust forms a crucial component to consider when building mock observables and analyzing observational data, and different dust models can strongly affect the results and conclusions. However, for lower-mass galaxies with relatively low gas masses and metallicities, dust is expected to affect the results less strongly. We therefore apply the well-known and often-used two-component dust model from Charlot & Fall (2000), which consists of a dust screen with a power-law dust attenuation curve with index Γ = −0.7 and a normalization of τ(5500 Å) = 0.33. Stars younger than 30 Myr are additionally attenuated following an identical power-law attenuation curve with τ(5500 Å) = 1.0. We discuss the effects of using different dust models in building mock spectra in T. K. Starkenburg et al. (2020, in preparation).

Hahn et al. (2021) used approximate Bayesian computation to infer an empirical dust model using the same set of mock galaxy spectra by fitting SDSS r-band magnitudes, g − r colors, and near-/far-UV colors. We have tested how the results in the present work change when using the best-fit dust empirical model. While the Hα EW measurements can be significantly affected by changing the dust model, Dn 4000 is relatively unchanged, as the index is not strongly sensitive to dust attenuation.

Using the synthetic spectra, we generate SDSS g- and r-band magnitudes, as well as Dn 4000 and Hα EW for all galaxies (Figure 1, top row). These quantities are noise-free and derived from spectra that encompass the total stellar light of each galaxy.

4.2. Mock Surveys

In generating synthetic photometric and spectroscopic quantities for each simulated galaxy, we are much closer to observational analogs for each simulation. However, the distribution of absolute magnitude and Dn 4000 shown in the top row of Figure 1 for each simulation are not directly comparable to the corresponding distribution in SDSS (top left), as the simulations are unaffected by observational noise and each sample contains every galaxy in the simulation volume. To accurately match the observations, we need to create mock surveys of each simulation box. 24

Our method for creating mock surveys is as follows.

  • 1.  
    Place an "observer" at a random location 10 Mpc outside the simulation volume. This distance is selected to match the lower distance limit on galaxies in the NSA catalog from SDSS.
  • 2.  
    Calculate apparent magnitudes, radial velocities, and projected distances for all galaxies as the observer would see them "on sky." For each galaxy, the radial velocity is the sum of the peculiar velocities along the observer's sight line and the recessional velocities given the distance between the galaxy and observer.
  • 3.  
    Convolve the synthetic spectra with the SDSS instrumental line-spread profile (modeled as a Gaussian with σ = 70 km s−1) and resample the spectra to match the SDSS wavelength resolution and coverage.
  • 4.  
    Add realistic spectral noise. The average S/N of SDSS spectra is dependent on galaxy color and apparent magnitude and also a function of wavelength. To reproduce the noise characteristics of SDSS, we bin galaxies from SDSS as a function of g − r color and apparent r magnitude. In each bin, we randomly draw 50 spectra and calculate the average S/N at each wavelength. In each mock observer's frame, noise is added to the synthetic spectra based on galaxy color, apparent magnitude, and wavelength by drawing from a normal distribution σ(λg − rmr ) such that the S/N matches the SDSS model. For simulated galaxies that fall outside the populated regions of the SDSS color–magnitude diagram, we select the closest bin in color–magnitude space.
  • 5.  
    Calculate stellar masses from the g − r colors following the prescription of Mao et al. (2021). Notably, this can result in stellar mass estimates that appear lower than the resolution limits of the simulation.
  • 6.  
    Remeasure Dn 4000 and Hα EW from the noise-added and instrumentally broadened spectra. Figure 2 shows the distribution of noise-added Dn 4000 and Hα EW for the simulations in three stellar mass bins, as compared to SDSS. We follow the definition of Balogh et al. (1999) for Dn 4000 measurements and Yan et al. (2006) for Hα EW measurements.
  • 7.  
    Remove "unobservable" galaxies from the sample. Spectroscopy was only acquired for galaxies in SDSS above the limiting magnitude mr  = 17.7. To accurately match SDSS, simulated galaxies that have apparent magnitudes fainter than 17.7 along a particular sight line must be removed from the sample.

In Figure 1, we show the effects of adding realistic noise, velocity resolution, and completeness cuts to isolated, M* < 1010 M galaxies in each of the three simulations. The top row shows absolute magnitude (Mr ) and Dn 4000 calculated from the noise-free synthetic spectra, while the bottom row shows Dn 4000 measured from the noise-added spectra and the apparent magnitudes (mr ) as seen by the observer along a random sight line.

The gray region in each panel in the bottom row of Figure 1 represents the regime in which galaxies are too faint to be selected for spectroscopic follow-up with SDSS. Of the three simulations, SIMBA produces the fewest "unobservable" galaxies. This is likely driven by the fact the SIMBA does not resolve galaxies all the way down to M* = 108 M, leading to fewer faint galaxies that are then preferentially removed from the observed sample (see Section 5.2 for a more detailed discussion on the effects of resolution in each simulation). The TNG has the faintest population of quiescent galaxies in absolute magnitude (at least in part because it extends to the lowest stellar masses), and correspondingly fewer quiescent galaxies from this simulation are observable with an SDSS-like survey.

In Figure 2, we show the distribution of Dn 4000 and Hα EW in three stellar mass bins for the observations and for the simulations as measured from the noise-added spectra. We color-code the simulations by instantaneous SFR, highlighting the need for synthetic observations to select the analogous quiescent population from simulations. As shown in Geha et al. (2012), no quiescent galaxies are observed in the SDSS volume below M* = 109 M.

4.3. Isolation Criteria

As in the observations, we select galaxies as isolated based on dhost, the projected distance between each simulated galaxy and its most nearby massive neighbor. As in the observations, for each galaxy, we identify potential hosts (galaxies with M* > 2.5 × 1010 M within 7 Mpc in projected distance and 1000 km s−1 in radial velocity). For galaxies with M* > 2.5 × 1010 M, we also require potential hosts to be more massive than the galaxy. Galaxies that have no potential hosts within 1.5 Mpc (dhost > 1.5 Mpc) are considered to be isolated.

In Figure 3, we compare our isolation criteria to the central/satellite classification from the simulations themselves. For all three simulations, that classification is based on a halo finder algorithm that uses the underlying dark matter structure (not available in observations or our mock observations). The specific algorithm varies between EAGLE/TNG and SIMBA. To define halos, both EAGLE and TNG use SUBFIND (Springel et al. 2001a) to find overdense, gravitationally bound (sub)structures within a larger connected structure that is found through a friends-of-friends (FOF; Davis et al. 1985) group-size halo finder. In SIMBA, galaxies are identified as FOF groups of stars and star-forming gas with a spatial linking length of 0.0056 times the mean interparticle spacing, while halos are identified as FOF groups with a linking length of 0.2 times the mean interparticle spacing. Galaxies and halos are cross-matched in postprocessing using the yt-based package CAESAR.

Figure 3.

Figure 3.  Purity (the fraction of galaxies selected as isolated by the dhost criterion that are centrals; solid lines) and completeness (the fraction of the population of centrals that appear as isolated when using dhost; dashed lines) of the dhost isolation criteria compared to the central/satellite classification from the simulations themselves. Shaded regions represent the 1σ variation across 10 randomly oriented sight lines. Using the dhost criterion selects a sample of galaxies that is relatively pure (>85% of isolated galaxies are centrals) but somewhat incomplete (only 55%–70% of centrals are selected as isolated). We use the true simulation-based stellar masses in this figure to provide maximum clarity on the purity and completeness of the dhost criterion.

Standard image High-resolution image

For all three simulations, the most massive subhalo (TNG and EAGLE) or galaxy (SIMBA) in a larger group halo is then classified as the central galaxy, while all other subhalos are classified as satellites.

For each simulation, we show the purity (the fraction of galaxies selected as isolated by the dhost criterion that are centrals; solid lines) and completeness (the fraction of the population of centrals that appear as isolated when using dhost; dashed lines) as a function of mass. For all simulations, purity declines with decreasing stellar mass, though all samples remain relatively pure even at low stellar mass. At M* = 109 M, a sample of isolated galaxies selected with dhost will be ∼85% centrals (∼15% of the "observed" isolated sample of galaxies are actually satellites, as defined by the halo finder).

Completeness is not strongly dependent on stellar mass below M* = 1010 M but decreases above this threshold. Completeness varies more between simulations, possibly due in part to the use of different halo finders. At M* = 109 M, 55%–70% of the centrals selected by the halo finder will be observed as isolated, reflecting the more restrictive nature of the dhost criterion. Recreating the observational selection criteria for isolation is a critical step in comparing between the population of isolated, quiescent galaxies selected from observations and those generated in simulations (see also Figure 6 and the discussion at the end of Section 5.1).

4.4. Quenching Criteria

As in the observations, simulated galaxies are considered to be quiescent if Hα EW  < 2 Å and Dn 4000 $\gt 0.6+0.1{\mathrm{log}}_{10}({M}_{* }/{M}_{\odot })$ (orange shaded region in Figure 2). Because Dn 4000 probes the luminosity-weighted stellar age of a galaxy, it is difficult to make a corresponding selection in SFR space, as highlighted by the color-coding of galaxies in Figure 2.

The star-forming and quiescent contours shown in Figure 1 can overlap because the Dn 4000 criterion is mass-dependent for each galaxy. Measuring Dn 4000 from the realistically noisy synthetic spectra is crucial to accurately match the sample selection from observations.

4.5. Volume Correction

To calculate the quiescent fraction fq , we weight each galaxy by the inverse of the total survey volume over which it could be observed given the SDSS spectroscopic magnitude limit ($1/{V}_{\max }$).

The quiescent fraction of isolated galaxies is then

Equation (1)

where Nq and NSF are the number of quiescent and star-forming galaxies in isolation in each mass bin, respectively.

5. The Quiescent Fraction of Isolated Galaxies

In Section 5.1, we compare the quiescent fraction of isolated galaxies from each simulation to observations from SDSS. In Sections 5.25.4, we discuss the effects of simulation resolution, observationally biased stellar mass, and finite aperture on the quiescent fraction.

5.1. Comparing between Simulations and Observations

In Figure 4, we show the median isolated quiescent fractions for each simulation from 25 randomly placed sight lines (blue, green, and purple points and dashed lines) as compared to SDSS (black circles and solid lines).

Figure 4.

Figure 4. Median quiescent fractions of isolated galaxies as a function of stellar mass for SDSS (black circles), Illustris-TNG (green squares), EAGLE (blue diamonds), and SIMBA (purple triangles). The quiescent fractions for the simulations are the medians of 25 randomly placed sight lines around each simulation box. Shaded regions represent the combination of binomial uncertainty on the quiescent fraction and the variance across sight lines for each simulation. Error bars on the SDSS quiescent fraction are adapted from Geha et al. (2012). Inset: same data shown in log-scale. Black arrows represent 1σ upper limits for the SDSS data in bins where the number of isolated quiescent galaxies is zero, and the number of isolated quiescent galaxies (Nq) in the lowest-mass bin is indicated for each simulation.

Standard image High-resolution image

At intermediate stellar masses (M* = 109.5−10.5 M), the simulations all show the quiescent fraction decreasing with decreasing stellar mass, qualitatively reproducing the isolated quenching threshold seen in SDSS (Geha et al. 2012). In observations, this is thought to be driven by the waning efficiency of internal feedback mechanisms. The threshold for TNG galaxies appears to be at a higher stellar mass than is seen in observations, while EAGLE's quiescent fraction appears to depend less strongly on stellar mass in this regime.

Unlike the observations, all three simulations have nonzero quiescent fractions below M = 109 M* (Figure 4, inset). For both SIMBA and TNG, the quiescent fraction, though nonzero, remains small, with fq  ∼ 0.01 at M = 109 M*. Overproduction of quiescent galaxies in EAGLE is more pronounced.

We examine the distribution of quiescent galaxies as a function of environment in two mass bins just above and below the observed quenching threshold for isolated galaxies (108.5–9.0 and 109.5–10.0 M*) in Figure 5. At intermediate stellar masses (top panel), EAGLE and SIMBA are a good qualitative match to observations, with low quiescent fractions for isolated galaxies and no environmental dependence in the isolated regime. The TNG lies below the SDSS observations, with an average fq  ∼ 0.05 beyond dhost  = 1.5 Mpc. At lower stellar masses (bottom panel), EAGLE's overproduction of quiescent galaxies in all environments is more pronounced, while SIMBA and TNG are in better agreement, with a near absence of quiescent galaxies beyond dhost  = 1.5 Mpc.

Figure 5.

Figure 5. Quiescent fraction of galaxies as a function of dhost in two mass bins. The top panel shows the observed distribution of galaxies with M* = 109.5–10.0 M, and the bottom panel shows M* = 108.5–9.0 M. The gray region indicates dhost< 1.5 Mpc, corresponding to nonisolated galaxies.

Standard image High-resolution image

In general, the rapid increase of fq at dhost < 1.5 Mpc is in qualitative agreement with observations. However, it is notable that all three simulations also overproduce low-mass nonisolated quiescent galaxies at a many-sigma tension given the derived errors. Feedback models for satellite galaxies are often tuned to reproduce galaxies in the Local Group. The SAGA Survey (Mao et al. 2021) recently showed that the Local Group satellites may be overquenched relative to those around more typical Milky Way–mass galaxies, which may be driving this tension.

In Figure 6, we show how the observed isolated quiescent fraction changes with different definitions of "isolation." For each simulation, we compare fq derived using dhost and fq from each simulation's halo finder. In each case, the halo finder produces a larger quiescent fraction (with the exception of M* > 1010.5 M in TNG and EAGLE). This effect varies as a function of stellar mass (as well as halo finder used) and highlights the differences between selecting a sample of central galaxies versus isolated galaxies.

Figure 6.

Figure 6. Effects of using a halo finder (filled symbols) vs. dhost (open symbols) to select isolated galaxies on the resultant quiescent fractions in TNG (top), EAGLE (middle), and SIMBA (bottom). In all three simulations, using the halo finder leads to higher quiescent fractions when compared to dhost, driven by the criterion's more restrictive nature.

Standard image High-resolution image

5.2. Resolution Effects

All three cosmological simulations used in this work have multiple runs at varying resolutions. SIMBA and EAGLE both have smaller boxes with higher resolution, and we use these boxes to test our results specifically with respect to resolution. While we do not use the higher-resolution version of TNG (TNG50), we note that the effect of resolution on the colors and color bimodality of galaxies is described in Nelson et al. (2018). They argue that the main effect of resolution on galaxy colors is in using purely the stellar particles to derive a star formation history. While we avoid this effect by including an instantaneous SFR in the star formation histories used to generate our spectra (to reflect the most recent star formation), future work should explicitly address the resolution convergence properties of TNG for the quiescent fractions of low-mass galaxies.

In Figure 7, we show the quiescent fraction as a function of stellar mass for the SIMBA and EAGLE higher-resolution boxes (lighter colors and dotted and dotted–dashed lines) and compare those to the default box and resolution from Figure 4 (darker colors and dashed lines). Because the high-resolution boxes are much smaller in volume, we determine the effect from cosmic variance on the quiescent fraction. To do so, we recompute the quiescent fraction for subboxes of 30 Mpc on a side from our default simulation boxes and show the full range of recovered quiescent fractions with the gray shaded regions in each panel. Due to the smaller total number of galaxies, we also use wider mass bins to calculate the quiescent fraction.

Figure 7.

Figure 7. Effects of resolution on the isolated quiescent fraction in EAGLE (top panel) and SIMBA (bottom panel). In both panels, dashed lines show the fiducial quiescent fractions from the large volume, while dotted or dotted–dashed lines show the quiescent fractions measured from the higher-resolution (25 Mpc)3 volumes. The gray shaded regions represent the maximum variation due to cosmic variance when considering (30 Mpc)3 subvolumes from each large-volume simulation. In the top panel, data from the two high-resolution simulations are offset horizontally so that the error bars may be distinguished.

Standard image High-resolution image

For the EAGLE simulations, we specifically compare to two boxes of 25 Mpc on a side, run using the "Reference" version of the code (identical to the 100 Mpc box), and using the "Recal" version of the code, which is recalibrated at this higher resolution of m DM  = 1.21 × 106 M and mgas = 2.26 × 105 M to counterbalance the resolution effects on the subgrid physics (see Schaye et al. 2015 and Schaller et al. 2015 for a complete argument and description of the recalibration and comparison of the different versions).

The quiescent fraction of the REF high-resolution box is slightly suppressed compared to that for the large box, falling below what can be explained by cosmic variance alone (gray shading) for the lowest masses. The recalibrated high-resolution box (RECAL) is further suppressed relative to both the fiducial and REF boxes, and this effect is most extreme at the lowest stellar masses. Schaye et al. (2015) argued that recalibrating at different resolutions is most appropriate, as the (un)physical effects of resolution may be hard to trace. Of the three EAGLE runs examined here, the RECAL box produces the fewest number of quiescent galaxies below M* = 109 M and even agrees with a quiescent fraction of fq  = 0 within the uncertainties for three of the four stellar mass bins. The recalibrated high-resolution EAGLE box is therefore in closest agreement with the observations, though fq does not show the same strong dependence on stellar mass in the mass range M* = 109–10 M that is observed in SDSS. Additionally, both high-resolution boxes still show an overabundance of low-mass quiescent galaxies in isolation.

For SIMBA, the quiescent fraction of low-mass galaxies increases slightly at higher resolution, which is run with identical physics to the larger-scale box. We find that the high-resolution box agrees with the fiducial run, with a quiescent fraction that is slightly elevated relative to the fiducial run at M* = 108.5–9.0 M. Improved resolution allows for the measurement of the quiescent fraction down to slightly lower stellar masses, which we find to be fq  ∼ 0.06. This is still in tension with SDSS, where no quiescent galaxies are observed in isolation below M* = 109 M.

We highlight that the effects of resolution, as exemplified by these case studies, are not uniform across simulations and can increase or decrease the quiescent galaxy fraction.

5.3. Stellar Mass

In earlier work, Munshi et al. (2013) found that the stellar masses measured for simulated galaxies with a combination of synthetic Petrosian magnitudes and Bell & de Jong (2001) mass-to-light (M/L) ratios underestimate the true stellar mass by about 50%. This underestimate varied only slightly as a function of mass. Similarly, Leja et al. (2019) showed that stellar masses inferred from SED modeling with nonparametric star formation histories are roughly ∼0.2 dex larger than those obtained under the assumption of an exponentially declining star formation history.

We confirm the approximate magnitude of this effect using the g − r color-based approximation for stellar mass from Mao et al. (2021), which produces stellar masses that underestimate the true stellar masses by ∼0.3 dex. In addition to shifting the quenching threshold to lower stellar masses, accounting for systemic error in the measurement of stellar mass leads to a softening of the slope of the quiescent fraction.

5.4. Aperture Effects

Our fiducial mock galaxy spectra are based on galaxy properties and star formation histories calculated using all of the star and gas elements considered to be part of a simulated galaxy's subhalo. In comparison, observations of galaxies in SDSS are aperture-limited. There are two SDSS apertures that are important for our results: the photometric aperture and the spectroscopic fiber aperture. The first roughly corresponds to what is considered the size of a galaxy and is relevant for the stellar masses used in this work, and the second is for the Dn 4000 and Hα EW measurements.

With respect to "galaxy size"–aperture effects, previous work using the EAGLE and TNG simulations found that these aperture effects are significant for high-mass (M* ≳ 1011 M) galaxies but negligible for lower-mass galaxies (Schaye et al. 2015; Donnari et al. 2021), and we reach similar conclusions when comparing stellar masses in the full subhalo and within a 30 kpc aperture for TNG. Therefore, we conclude that aperture effects are likely insignificant for stellar masses, in particular when compared to other uncertainties, as discussed in the previous section.

The aperture of the SDSS spectroscopic fiber typically covers a few kpc in the central region of a galaxy. This means that aperture effects are again crucial to take into account at higher masses but can be small for low-mass galaxies, as there, the SDSS fiber aperture may cover a significant fraction of the galaxy. We have checked the differences in Dn 4000 and Hα EW when only considering the central r < 2 kpc region of each galaxy. While Hα EW is more variable as a function of aperture (driven by changes in the amount of continuum contained within the aperture), galaxies do not significantly shift in or out of the quiescent region (which is based on both Hα EW and Dn 4000). In particular, we find only small differences in the quiescent fractions for low-mass galaxies, which shift upward by ∼10%.

For some observed galaxies, the SDSS fiber may have been centered on an off-center bright (star-forming) region. With smaller (lower-mass) galaxies, this possibility diminishes purely due to the fact that more of the total galaxy fits within the fiber aperture. Moreover, this is more likely to happen in actively star-forming galaxies, and as we purely use the spectral indices to classify galaxies as star-forming or quiescent, off-center fiber positioning is unlikely to affect our results.

6. Quenching Mechanisms in Simulations

6.1. Comparison to Previous Work

Our work is not the first to compare populations of quiescent galaxies from simulations to observations, though we are the first to compare three large-volume simulations to observations in a fully consistent manner. Here we review past studies evaluating the quiescent populations of each simulation.

6.1.1. EAGLE

The quiescent fraction of galaxies in the EAGLE simulations has been discussed in a multitude of papers using a variety of galaxy subsamples and of definitions and tracers of the quiescent fraction. Schaye et al. (2015) found that the passive fraction for all galaxies defined based on specific SFRs (sSFRs; SFR/M) is in reasonably good agreement with results from GAMA (Bauer et al. 2013) and SDSS (Moustakas et al. 2013), although these observational results use slightly different criteria. The threshold where the galaxy population goes from predominantly blue and star-forming to predominantly red and quiescent is found to be at a slightly higher stellar mass in EAGLE than for observed low-redshift data sets (Schaye et al. 2015; Trayford et al. 2015). Furlong et al. (2015) and Trayford et al. (2015, 2017) showed that the apparent larger quiescent fraction at the low-mass end is largely due to resolution effects, as low-mass, low-SFR galaxies have very few star-forming gas particles, and that the recalibrated higher-resolution box (Recal-25) agrees with observations.

We agree with their results that the recalibrated high-resolution box improves the quiescent fraction compared to the observations (see Figure 7 and Section 5.2). Nevertheless, we still find that EAGLE produces a higher quiescent fraction for low-mass galaxies than is observed in SDSS. It is noteworthy that EAGLE has a particularly high fraction of low-mass galaxies that are classified as quiescent but have current star formation at very low rates. This gives credence to the possibility that the SFR of these galaxies is poorly resolved and that they are misclassified as quiescent. However, our quiescent definition requires both low Hα and high Dn 4000. While the Hα can be affected by these resolution effects, the Dn 4000 is based on the stellar continuum and probes a longer timescale. Our quiescent fraction may therefore be less dependent on resolution than for a purely Hα-based definition.

At higher masses, the quiescent fraction we find for EAGLE seems particularly low compared to earlier work based on SFR or sSFR. However, Trayford et al. (2017) showed similar results for Dn 4000-defined quiescent fractions. In addition, they showed that dust can affect the discrepancy; when including dust, the passive fraction in EAGLE based on a cut of Dn 4000 > 1.8 is 35% reduced compared to results from SDSS for galaxies in the mass range 1010 M < M* < 1011 M, while without including dust, this discrepancy increases to 70% (Trayford et al. 2017). These discrepancies are larger than what we find in this work and may be because our stellar mass–dependent Dn 4000-based definition of quiescence, which is identical to the one used for the SDSS sample, reaches a lower value at similar masses.

6.1.2. Illustris-TNG

Donnari et al. (2020, 2021) provided an in-depth exploration of the quenched galaxy fraction in TNG, exploring the effects of aperture, quenching definition, SFR timescale, environmental misidentification, and incompleteness on the quenched fraction. While they focused on galaxies with M* > 109 M and used SFR- or SPS-based definitions of quiescence, we nonetheless find their fiducial quiescent fraction to be in excellent qualitative agreement with our results, specifically the existence of a quenching threshold just above M* = 1010 M, a small population of isolated (central) galaxies that are quiescent below M* = 109.5 M, and a rapid increase in the number of quiescent galaxies from M* = 1010–11 M.

Donnari et al. (2019) described the galaxy star-forming sequence and the quenched fraction using different definitions and tracers for TNG100. They found that both a UVJ-selected quenched fraction and a distance from the star-forming sequence-selected quenched fraction agree reasonably well with observations, although they used a different UVJ selection for TNG than for the observations they compare to (from Muzzin et al. 2013).

In a by-eye comparison, our TNG quiescent fractions (defined based on Dn 4000 and Hα) are overall lower than the UVJ- and SFS-based quenched fraction from Donnari et al. (2019), in particular at the lower-mass end. A similar difference can be found between the UVJ-selected (Muzzin et al. 2013) and Dn 4000- and Hα-selected (Geha et al. 2012) observed quiescent fractions at M < 1010 M. As it is unlikely that low-mass star-forming galaxies move into the UVJ selection region due to dust, this suggests that these low-mass systems are predominantly red but still exhibit some low (relatively recent) star formation that can be traced by their Dn 4000 and/or Hα EW.

6.1.3. SIMBA

Davé et al. (2019) compared the sSFRs of SIMBA galaxies to observations from the GALEX-SDSS-WISE Legacy Catalog (Salim et al. 2016, 2018) and found good agreement, which is also consistent with our findings of very good agreement between the quiescent fraction from SIMBA and that observed in SDSS.

Davé et al. (2019) also suggested that the few low-mass quiescent galaxies are in fact satellites misclassified by the FOF-based halo finder. Additionally, Christiansen et al. (2020) and Borrow et al. (2020) showed that jet feedback from AGNs in SIMBA can influence large regions around the AGN host galaxy and entrain and influence gas in nearby galaxies. The SIMBA 100 Mpc box contains a massive cluster (M* ∼ 1015 M), and the high-resolution 25 Mpc box also has an anomalously large central halo. However, we confirm that in both boxes, the majority of the quiescent galaxies we select are at least 2 Mpc away from the cluster or most massive halo center, and most are more than 5 Mpc away. Our results suggest that there may be additional effects driving the nonzero quiescent fraction at low mass (see Figure 5, bottom panel).

Low-mass quiescent galaxies in SIMBA are also found to be oversized compared to their star-forming counterparts (Davé et al. 2019). However, this should lead to a higher fraction of quiescent galaxies removed from the survey due to the SDSS surface brightness limits. Fixing this issue would therefore only increase the population of low-mass quiescent galaxies in SIMBA.

6.2. Quenching Mechanisms of Low-mass Isolated Galaxies in Simulations

Given the observed overabundance of quiescent galaxies at low mass in the simulations, we consider what feedback mechanisms might drive quenching in these systems.

Black holes. Some of the low-mass galaxies in our sample will contain black holes (almost all galaxies above M* = 109 M in TNG and EAGLE and M* = 109.5 M in SIMBA). Feedback from central supermassive black holes is often thought to be able to quench galaxies and is possibly important in keeping galaxies quiescent (e.g., Somerville & Davé 2015; recent discussions include Terrazas et al. 2016, 2020; Chen et al. 2020). However, in the models, the black hole feedback often becomes effective only at certain black hole masses (Bower et al. 2017; McAlpine et al. 2017; Weinberger et al. 2018; Thomas et al. 2019; Habouzit et al. 2021; Terrazas et al. 2020), which are not reached in low-mass galaxies. Moreover, not all of the low-mass galaxies in these simulations will host black holes (especially at M* < 109 M), so while black hole feedback could play some role, it is unlikely to be the exclusive driver of the quenching of the lowest-mass galaxies in each simulation.

Star formation feedback. Feedback from star formation is generally thought to reduce the efficiency of galaxy formation at the lower-mass end (Schaller et al. 2015; Schaye et al. 2015; Somerville & Davé 2015; Pillepich et al. 2018b; Davé et al. 2019). Outflows from star formation feedback can temporarily expel large amounts of gas and decrease the gas density of galaxies. This is commonly seen in higher-resolution simulations of lower-mass galaxies (halo masses ≲1010 M; Wheeler et al. 2015, 2019; Wright et al. 2019; Rey et al. 2020). These effects may result in temporary quiescence, but it is unclear whether supernovae are capable of removing enough gas to fully quench galaxies at M* ∼ 109 M. Temporarily quiescent low-mass galaxies may, however, compose part of our quiescent low-mass galaxy samples. In the large-scale simulations used in this work, the effectiveness of stellar feedback is likely also affected by resolution at the lowest-mass end, and the effects of feedback and resolution may be hard to disentangle.

Splashbacks. Galaxies moving through larger halos can have their gas stripped away, reducing star formation. Splashback galaxies can be found up to  ∼3Rvir (see, e.g., Diemer 2021) from their true host galaxy, making misidentification possible. However, our dhost isolation criterion is more restrictive than the halo finders (see Figure 6), making it unlikely that more than a small fraction of the isolated quenched galaxies observed in each simulation are splashbacks.

Outflows from nearby massive galaxies. Borrow et al. (2020) and Wright et al. (2019) showed that the gas of low-mass galaxies (both satellites and low-mass centrals) can be stripped or entrained by jets and AGN outflows from more massive halos. While many of the observed isolated galaxies lie many virial radii from any massive halos, this gas removal effect could contribute to the elevated quenched fractions seen in the simulations.

7. Star Formation Histories of Quiescent Galaxies

By forward modeling simulated galaxies into observational space, we gain the ability to look back at the "true" properties of observationally selected galaxies. One of the most informative properties of a galaxy is its star formation history. For each simulated galaxy, we are able to compare the observed properties back to the biases inherent in the recovery of star formation histories from real data.

In Figure 8 (top row), we show the median star formation histories of low-mass (M* = 109.0–9.5 M) isolated galaxies observed as star-forming (blue dashed lines) and quiescent (orange solid lines) at z = 0 in each simulation. The shape of the median star formation histories in each simulation varies significantly. Low-mass star-forming galaxies in SIMBA show late-time rising star formation histories, while the same galaxies in TNG have approximately flat star formation histories, both in clear contrast with their quiescent counterparts.

Figure 8.

Figure 8. Top row: median star formation histories of low-mass (M* = 109.0–9.5 M) isolated galaxies in each simulation, split into star-forming (blue dashed) and quiescent (orange solid) using the Dn 4000–Hα EW criteria at z = 0. Shaded regions encompass the 20th–80th percentiles of each distribution. Bottom row: median cumulative stellar mass as a function of look-back time for the star-forming and quiescent samples of isolated galaxies at low mass in each simulation. Vertical dotted–dashed lines show the average time at which 90% of a galaxy's stellar mass has formed (τ90). Though the star formation histories of quiescent galaxies vary across the three simulations, the average quenching timescales derived for the three populations are effectively identical (τ90 ∼ 5 Gyr).

Standard image High-resolution image

In EAGLE, the star formation histories are declining at late times for both the star-forming and quiescent populations, and in fact, a significant fraction of the low-mass quiescent galaxies appears to be forming stars at very low rates at late times. This may be connected to the relatively low scatter in the EAGLE galaxy star-forming sequence (see Hahn et al. 2019) and star formation histories (Iyer et al. 2020). In particular, the declining and low SFRs for EAGLE's low-mass galaxies may be connected to the relatively high quiescent fraction that we find using the Dn 4000 and Hα EW selection. Figure 2 illustrates that, of the three simulations in our sample, EAGLE has a notably large fraction of low-mass, low-SFR galaxies that are borderline star-forming when considering their sSFRs but have such low SFRs that they are classified as quiescent when using the Dn 4000 and Hα EW measurements.

In the bottom row of Figure 8, we show the median cumulative stellar mass as a function of look-back time for the same star-forming and quiescent samples. The vertical dotted–dashed lines show the average time at which quiescent and star-forming galaxies formed 90% of their stellar mass (Skillman et al. 2017; Weisz et al. 2015). Despite the apparent variation in the late-time star formation histories of the quiescent populations observed in each simulation, τ90 is nearly identical: τ90, TNG = 4. 9 ± 1.6 Gyr, τ90, EAGLE = 4. 9 ± 1.4 Gyr, and τ90, SIMBA  = 4.6 ± 2.0 Gyr. The measurement of τ90 ∼ 5 Gyr for low-mass quiescent galaxies in isolation is a testable prediction for the observational sample and should provide insight into the timescale of feedback mechanisms that drive quenching in low-mass galaxies.

Of the galaxies observed as quiescent, a subset in each simulation has nonzero instantaneous SFRs at z = 0 (approximately 50% of low-mass galaxies in TNG and EAGLE and 20% in SIMBA). The empirical definition of quiescence used throughout this work (based on Dn 4000 and Hα EW) selects galaxies with a range of evolutionary histories and z = 0 properties. We find that the SFR = 0 galaxies do not separate out cleanly in Dn 4000–Hα EW space (highlighted in Figure 2). In order to select a totally quiescent sample of galaxies, additional probes of SFR would be required. Characterizing exactly how much ongoing star formation a galaxy can experience while still being selected by a given definition of quiescence will help inform the selection of appropriate quenching criteria (e.g., UVJ versus SFR versus Dn 4000) going forward.

In future work, we will constrain the average quenching time for observed quiescent galaxies in isolation via SED fitting and apply the same process to the synthetic spectra. Applying the same SED fitting methods to the synthetic data is critical, as the observed differences in the late-time star formation history may only be recoverable above a certain threshold of data quality (e.g., spectrum S/N, photometric wavelength coverage).

8. Conclusions

In this paper, we have produced mock SDSS-like surveys from three large-volume hydrodynamic simulations (Illustris-TNG, EAGLE, and SIMBA) from which we measured the quiescent fraction of isolated galaxies and compared back to extant constraints from the local universe (Geha et al. 2012). We find the following.

  • 1.  
    The three simulations examined in this work, when transformed into observational space using an identical methodology, produce three different dependencies of the quiescent fraction on stellar mass. Above M* = 109.5 M, all three simulations qualitatively reproduce the declining quiescent fraction with decreasing stellar mass observed in SDSS.
  • 2.  
    All three simulations have nonzero quiescent fractions below M* = 109.0 M, in contrast to observations of the SDSS volume. This suggests that current models of feedback in the low stellar mass regime are too efficient.
  • 3.  
    Our empirical definition of quiescence selects low-mass galaxies with a range of star formation histories when viewed over longer (many Gyr) timescales. However, these populations all show similar quenching timescales (τ90 ∼ 5 Gyr), which can be compared back to observations. Understanding how sensitive a particular definition of quenching is to formation history can inform future population studies.
  • 4.  
    Measuring the quiescent fraction in a higher-resolution box does not fully resolve the overabundance of quiescent galaxies below M* = 109.0 M. In fact, improved resolution can lead to either a significant decrease or slight increase in the measured quiescent fraction, depending on the simulation. While the low-mass quiescent fraction is not converged in the large-scale simulations, the discrepancy with observations persists at higher resolutions.

The low-mass quenching threshold of isolated galaxies represents an observational boundary condition: a stellar mass regime where internal feedback mechanisms become ineffective. Observations of the isolated galaxy quiescent fraction provide a unique probe into the delicate balance of internal feedback mechanisms in low-mass galaxies. Understanding how well or poorly modern simulations of galaxy evolution can reproduce this feature is a novel test of feedback prescriptions but requires the creation of mock observations and surveys to enable appropriate comparisons between the observed universe and simulations.

Our method for producing mock surveys from large-volume hydrodynamic simulations can also be applied to zoom-in simulations and semianalytic models and adapted to match other surveys and observations. Future work will explore the star formation histories recovered from synthetic observations as compared to those derived from observations, as well as the observed gas properties of simulated galaxies.

The constraining power of the observations we compare to are set by the size of the SDSS volume and the limiting magnitude and surface brightness of the survey. Future wide-field surveys, such as the Vera Rubin Observatory's Large Synoptic Survey Telescope (LSST; Ivezić et al. 2019) and the Dark Energy Survey Instrument (DESI; DESI Collaboration et al. 2016), have the potential to substantially improve our census of the local universe, providing new constraints on the population of low-mass quiescent galaxies in the field.

Making direct comparisons between observations and simulations requires careful translation from the simulation to observational frame (or vice versa). In doing so, we gain novel insights into the ways that feedback shapes the evolution of galaxies.

We thank the Illustris collaboration and the Virgo Consortium for making their simulation data publicly available and the SIMBA collaboration for sharing their data with us. The EAGLE and SIMBA simulations were performed using the DiRAC-2 facility at Durham, managed by the ICC, and the PRACE facility Curie, based in France at TGCC, CEA, Bruyères-le-Châtel.

This research was supported in part through the computational resources and staff contributions provided by the Quest high-performance computing facility at Northwestern University, which is jointly supported by the Office of the Provost, the Office for Research, and Northwestern University Information Technology. The data used in this work were, in part, hosted on facilities supported by the Scientific Computing Core at the Flatiron Institute, a division of the Simons Foundation, and the analysis was largely done using those facilities.

The Isolated & Quiescent (IQ) Collaboratory thanks the Flatiron Institute for hosting the collaboratory and its meetings. The Flatiron Institute is supported by the Simons Foundation. C.M.D. is supported by a Professor's Grant from the Howard Hughes Medical Institute (PI: Geha). N.M. acknowledges support from the Klauss Tschira Foundation through the HITS Yale Program in Astrophysics (HYPA).

Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy Office of Science, and the Participating Institutions.

The SDSS-IV acknowledges support and resources from the Center for High Performance Computing at the University of Utah. The SDSS website is www.sdss.org. The SDSS-IV is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS Collaboration, including the Brazilian Participation Group, the Carnegie Institution for Science, Carnegie Mellon University, the Center for Astrophysics ∣ Harvard & Smithsonian, the Chilean Participation Group, the French Participation Group, Instituto de Astrofísica de Canarias, The Johns Hopkins University, the Kavli Institute for the Physics and Mathematics of the Universe (IPMU)/University of Tokyo, the Korean Participation Group, Lawrence Berkeley National Laboratory, Leibniz Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Astrophysik (MPA Garching), Max-Planck-Institut für Extraterrestrische Physik (MPE), the National Astronomical Observatories of China, New Mexico State University, New York University, the University of Notre Dame, Observatário Nacional/MCTI, The Ohio State University, Pennsylvania State University, Shanghai Astronomical Observatory, the United Kingdom Participation Group, Universidad Nacional Autónoma de México, the University of Arizona, the University of Colorado Boulder, the University of Oxford, the University of Portsmouth, the University of Utah, the University of Virginia, the University of Washington, the University of Wisconsin, Vanderbilt University, and Yale University.

Software: Astropy (Astropy Collaboration et al. 2013, 2018), FSPS (Conroy et al. 2009; Conroy & Gunn 2010), IPython (Pérez & Granger 2007), Matplotlib (Hunter 2007), NumPy (Harris et al. 2020), Python (Van Rossum & Drake 2009), Python-FSPS (Foreman-Mackey et al. 2014), SciPy (Virtanen et al. 2020).

Footnotes

Please wait… references are loading.
10.3847/1538-4357/abc014