<figure><img src="image_1.png" /></figure>

# Quenching as a Contest between Galaxy Halos and Their Central Black Holes

Zhu Chen¹, S. M. Faber², David C. Koo², Rachel S. Somerville³,⁴, Joel R. Primack⁵, Avishai Dekel⁶,
Aldo Rodríguez-Puebla⁷, Yicheng Guo⁸, Guillermo Barro⁹, Dale D. Kocevski¹⁰, A. van der Well¹,¹², Joanna Woo¹³,¹⁴,
Eric F. Bell¹⁵, Jerome J. Fang¹⁶, Henry C. Ferguson¹⁷, Mauro Giavalisco¹⁸, Marc Huertas-Company¹⁹,²⁰,²¹,²²,
Fangzhou Jiang⁶, Susan Kassin¹⁷, Lin Lin²³, F. S. Liu²⁴,²⁵, Yifei Luo², Zhijian Luo¹, Camilla Pacifici²⁶, Viraj Pandya²,
Samir Salim²⁷, Chenggang Shu¹, Sandro Tacchella²⁸, Bryan A. Terrazas¹⁵, and Hassen M. Yesuf²⁹

¹Shanghai Key Lab for Astrophysics, Shanghai Normal University, 200234, Shanghai, People's Republic of China

²Department of Astronomy and Astrophysics, University of California, Santa Cruz, CA 95064, USA

³Center for Computational Astrophysics, Flatiron Institute, New York, NY 10010, USA

⁴Department of Physics and Astronomy, Rutgers University, Piscataway, NJ 08854, USA

⁵Physics Department, University of California, Santa Cruz, CA 95064, USA

⁶Racah Institute of Physics, The Hebrew University, Jerusalem 91904, Israel

⁷Instituto de Astronomía, Universidad Nacional Autónoma de México, México, D.F., Mexico

⁸Department of Physics and Astronomy, University of Missouri, Columbia, MO 65211, USA

⁹Department of Physics, University of the Pacific, Stockton, CA 95211, USA

¹⁰Department of Physics and Astronomy, Colby College, Waterville, ME 04901, USA

¹¹Max-Planck-Institut für Astronomie, D-69117 Heidelberg, Germany

¹²Sterrenkundig Observatorium, Universiteit Gent, Krijgslaan 281 S9, B-9000 Gent, Belgium

¹³Department of Physics & Astronomy, University of Victoria, Victoria BC V8P 1A1, Canada

¹⁴Department of Physics, Simon Fraser University, 8888 University Dr, Burnaby BC, V5A 1S6, Canada

¹⁵Department of Astronomy, University of Michigan, Ann Arbor, MI 48109-1107, USA

¹⁶Astronomy Department, Orange Coast College, Costa Mesa, CA 92626, USA

¹⁷Space Telescope Science Institute, 3400 San Martin Road, Baltimore, MD 21218, USA

¹⁸Department of Astronomy, University of Massachusetts, Amherst, MA 01003, USA

¹⁹LERMA, Observatoire de Paris, PSL Research University, CNRS, Sorbonne Universités, UPMC Univ. Paris 06, F-75014 Paris, France

²⁰Univeristé de Paris, 5 Rue Thomas Mann, F-75013 Paris, France

²¹Departamento de Astrofísica, Universidad de La Laguna, E-38206 La Laguna, Tenerife, Spain

²²Instituto de Astrofísica de Canarias, E-38200 La Laguna, Tenerife, Spain

²³Shanghai Astronomical Observatory of CAS, Shanghai 200030, People's Republic of China

²⁴College of Physical Science and Technology, Shenyang Normal University, Shenyang 110034, People's Republic of China

²⁵University of California Observatories and the Department of Astronomy and Astrophysics, University of California, Santa Cruz, CA 95064, USA

²⁶Goddard Space Flight Center, Code 665, Greenbelt, MD 20771, USA

²⁷Department of Astronomy, University of Indiana, Bloomington, IN 47404, USA

²⁸Harvard-Smithsonian Center for Astrophysics, Cambridge, MA 02138, USA

²⁹Kavli Institute for Astronomy and Astrophysics, Peking University, Beijing 100871, People's Republic of China

Received 2019 September 25; revised 2020 May 2; accepted 2020 May 23; published 2020 July 7

## Abstract

Existing models of galaxy formation have not yet explained striking correlations between structure and star formation activity in galaxies, notably the sloped and moving boundaries that divide star-forming from quenched galaxies in key structural diagrams. This paper uses these and other relations to “reverse engineer” the quenching process for central galaxies. The basic idea is that star-forming galaxies with larger radii (at a given stellar mass) have lower black hole (BH) masses due to lower central densities. Galaxies cross into the green valley when the cumulative effective energy radiated by their BH equals ~4× their halo gas-binding energy. Because larger-radii galaxies have smaller BHs, one finds that they must evolve to higher stellar masses in order to meet this halo energy criterion, which explains the sloping boundaries. A possible cause of radii differences among star-forming galaxies is halo concentration. The evolutionary tracks of star-forming galaxies are nearly parallel to the green-valley boundaries, and it is mainly the sideways motions of these boundaries with cosmic time that cause galaxies to quench. BH scaling laws for star-forming, quenched, and green-valley galaxies are different, and most BH mass growth takes place in the green valley. Implications include the radii of star-forming galaxies are an important second parameter in shaping their BHs; BHs are connected to their halos but in different ways for star-forming, quenched, and green-valley galaxies; and the same BH-halo quenching mechanism has been in place since z ~ 3. We conclude with a discussion of BH-galaxy coevolution and the origin and interpretation of BH scaling laws.

*Unified Astronomy Thesaurus concepts: Galaxy evolution (594); Galaxy quenching (2040); Galaxy structure (622); Galaxy physics (612); Star formation (1569); Supermassive black holes (1663)*

## 1. Introduction

How and why star formation comes to an end (“quenches”) in massive galaxies is one of the outstanding puzzles of galaxy evolution. Though gaseous infall onto halos is decreasing in the current universe, there is consensus that this effect by itself is not

sufficient to curtail star formation to the low levels that are observed in quenched galaxies (see review by Somerville & Davé 2015), and in any case, quenching began at $z \gtrsim 3$, when the universe was still quite gas rich. The number of massive quenched galaxies at the knee of the quenched mass function has
approximately tripled since $z = 1$ (Bell et al. 2004; Faber et al. 2007; Huang et al. 2013; Ilbert et al. 2013; Muzzin et al. 2013; Tomczak et al. 2014), and the vast majority of that growth is in central galaxies, not satellites (Tinker et al. 2013). It is therefore clear that quenching of central galaxies is continuing today, but whether by the same or different mechanisms as at earlier times is not known.<fnref n="30" />

Many studies have modeled galaxy quenching, from semianalytic models (SAMs) using parametric prescriptions, to hydrodynamic simulations with more detailed descriptions of gas and feedback physics, and basic data like the mass functions of quenched and star-forming galaxies have been more or less successfully matched. However, in our opinion, model makers have not yet focused fully enough on the full panoply of physical data exhibited by quenching galaxies versus mass and time. In this paper, we assemble a more complete list of relevant observational data and use it to construct a simple picture in which galaxies evolve through a structural parameter space while star-forming (and building black holes, BHs) and then encounter a “quenching boundary,” which is defined in this paper as the entrance to the green valley. As shown below, this boundary is visible as the dividing line between quenched and star-forming galaxies in various structural diagrams. Key properties of the boundary (see Section 1.2 and Figure 1) are its considerable length in stellar mass (more than 1 dex), its sharpness, smooth motion in the diagrams with time, and its slope or tilt whereby more massive galaxies along the boundary have larger radii, higher central densities, and larger central BHs than smaller galaxies.

If quenching depended only on stellar mass, the boundary lines in structural diagrams would be concentrated at a given stellar mass. The fact that they are extended and sloping implies that at least one more parameter in addition to mass is involved in quenching. We show in this paper that taking this second parameter to be the effective radius and tying the BH mass to it can naturally explain most data. This is demonstrated by showing that the observed boundaries are very close to those obtained if galaxies enter the green valley when the total effective energy emitted by their central BHs equals a certain (small) multiple of the binding energy of the hot gas in their dark halos.

Before quenching, the model starts by placing star-forming galaxies in dark halos according to abundance matching, and their stellar masses, radii, and BH masses are computed from their halo properties using simple rules. As galaxies enter the green valley, we find that BH growth accelerates, which leads ultimately to full quenching. The full history of BH growth from star-forming to quenched is sketched in the “cartoon” shown in Figure 4.

As of this writing, the exact nature of the mechanism that quenches galaxies is unknown. Two broad classes have been suggested. One is so-called “halo quenching,” which involves the fact that halo gas becomes shock-heated as halos develop and is therefore less able to cool and fall onto galaxies (e.g., Ostriker & Rees 1977; Blumenthal et al. 1984; Birnboim & Dekel 2003; Cattaneo et al. 2006; Dekel & Birnboim 2006; Keres et al. 2009). The other is “active galactic nucleus (AGN)

quenching” by feedback from growing BHs, which either ejects gas from galaxies (“ejective feedback”; e.g., Silk & Rees 1998; Hopkins et al. 2006; Springel et al. 2006) or heats or disturbs the halo gas in such a way as to prevent it from falling in (“preventive feedback,” “strangulation”; e.g., Bower et al. 2006; Croton et al. 2006; Somerville et al. 2008, 2012; Henriques et al. 2015; Peng et al. 2015).<fnref n="31" />

Both halo quenching and AGN quenching have considerable observational support. On the halo side, N-body simulations provide reliable models of halo growth (e.g., Behroozi et al. 2013), the basic physics of halo gas density and temperature within those halos is reasonably well known (see references above), and X-ray data confirm that massive halos are filled with hot gas at the predicted temperatures and densities (e.g., Arnaud et al. 2010). On the AGN side, nearly all quenched galaxies have massive spheroids (Bell 2008; Bluck et al. 2014), and these spheroids harbor massive BHs (Magorrian et al. 1998; Ferrarese & Merritt 2000; Gebhardt et al. 2000). AGN activity peaks just as galaxies are entering quiescence (e.g., Nandra et al. 2007; Schawinski et al. 2007, 2010; Aird et al. 2019). Tight correlations exist between BH mass and host spheroid properties, such as bulge mass ($M_{BH}-M_{bulge}$; Häring & Rix 2004; Kormendy & Ho 2013) and effective velocity dispersion ($M_{BH}-\sigma_e$; Ferrarese & Merritt 2000; Gebhardt et al. 2000; Kormendy & Ho 2013; Ho & Kim 2014; van den Bosch 2016), signaling some sort of close interplay between spheroid building and BH growth. Finally, matter falling onto a BH may release of order 10% of its rest-mass energy $E = mc^2$, which is comparable in amount to the binding energy of gas in the halo (e.g., Silk & Rees 1998; Ostriker et al. 2005; Oppenheimer 2018), enough to either heat, blow away, or otherwise significantly disturb the gas and induce quenching. For reviews of these topics, see Cattaneo et al. (2009) and Kormendy & Ho (2013, hereafter KH13).

Despite these successes, neither class of model as currently formulated naturally explains the sloping nature of the boundaries that divide quenched from star-forming galaxies above $M_* \sim 10^{10} M_\odot$. Halo quenching depends on only the temperature and density of the halo gas and thus tends to set in at a characteristic halo mass near $M_{crit} \sim 10^{11-12} M_\odot$ (e.g., Birnboim & Dekel 2003; Keres et al. 2005; Dekel & Birnboim 2006), and therefore at a corresponding stellar mass because the stellar-mass/halo-mass relation is tight (Wechsler & Tinker 2018). Halo quenching therefore tends to predict a too-sharp boundary in stellar mass near $M_* \sim 3 \times 10^{10} M_\odot$ and also does not naturally explain the systematic variation (sloping boundaries) in other galaxy properties near that mass.

The predicted observational signature of AGN quenching depends on the mode of feeding mass onto the BH. In “radiative mode” (Somerville et al. 2018; also called “quasar mode,” Croton et al. 2006; “bright mode,” Somerville et al. 2008; or “high mode,” Pillepich et al. 2018), AGN mass accretion is rapid and is often triggered by major mergers (Hopkins et al. 2006) and/or disk-gas instabilities (Cattaneo et al. 2006). Because these are random events that can occur over a broad range of time in a galaxy’s life, it is hard to see why the quenching boundary should be sharp, let alone particularly sloped. In “jet mode” (Somerville et al. 2018; also

<footnote n="30">³⁰ In this paper, quenched galaxies are defined as objects lying more than −1.0 dex below the local star-forming main-sequence ridgeline at their mass, and star-forming galaxies are galaxies that lie above the line lying −0.45 dex below the ridgeline. Green-valley galaxies lie between these populations from −0.45 dex to −1.0 dex.</footnote>

<footnote n="31">³¹ Feedback is also produced by young stars and supernovae, but the characteristic energy is too small to affect gas in massive galaxies. The major consequence is to modulate the gas content of lower-mass galaxies while they are still star-forming (Somerville & Davé 2015).</footnote>
called “radio mode,” Croton et al. 2006; “kinetic mode,”
Weinberger et al. 2017; or “low mode,” Pillepich et al. 2018),
gas drains onto the BH more slowly, creating feedback energy
and/or momentum that ejects gas and/or heats the halo.
However, highly collimated jets cannot couple efficiently to the
halo gas while the gas is still cool and dense (e.g., Bower et al.
2006; Cattaneo et al. 2006; Croton et al. 2006), and so in
practice, jet-like feedback also needs to wait until the halo mass
exceeds Mcrit. Prospects may be better in kinetic mode
simulations that are more isotropic (e.g., IllustrisTNG; Pillepich
et al. 2018), but sloping boundaries have not yet been
demonstrated.

The above is a brief summary of a rich and complex literature, and variations in points of view and terminology are common. For example, Peng et al. (2010) also identified two quenching mechanisms, which they called “mass quenching” and “environment quenching.” However, their mass quenching is a more general feedback-driven process that includes both stellar feedback as well as AGN feedback, while their environment quenching refers to satellite galaxies only. We ignore satellites in this paper owing to their different quenching physics (e.g., Woo et al. 2017), and halo quenching, when it is mentioned, always refers to central galaxies.

The upshot of the preceding is that current models have not
yet explained the extent of the quenching boundary in stellar
mass or why the radii, central densities, and BH masses of
boundary galaxies vary systematically along it. In this paper,
we assume that the quenching process has imprinted these
trends, and we attempt to “reverse engineer” the physical nature
of the quenching mechanism from them. This leads us to a
picture in which galaxy quenching is dominated by AGN
quenching and becomes perceptible when total effective BH
energy output exceeds a certain multiple of the gravitational
binding energy of the halo gas. The point of perceptibility is
defined to be −0.45 dex below the star-forming main sequence,
where the number of green-valley galaxies begins to exceed the
lower tail of the star-forming main sequence in CANDELS
(Fang et al. 2018). For convenience, we will call this the
“quenching boundary” even though galaxies have presumably
started to quench before passing this point.

The next section reviews previous works that have also tried
to use the structural properties of galaxies as a guide to
understanding quenching. It is followed by a section that
collects together all of the basic properties of quenching and
star-forming galaxies that a satisfactory theory of quenching
should explain. The final section in this Introduction sketches
the nature of the model and outlines the rest of the paper.

1.1. Previous Work

Many studies have remarked on the striking structural
parameter correlations of both quenched and quenching
galaxies. Let R_e be the optical effective radius containing half
the light and Σ_e be the effective mass surface density given by
Σ_e = 0.5M_*/πR_e^2. Kauffmann et al. (2003) showed that
quenched SDSS galaxies have higher Σ_e than star-forming
galaxies at the same mass, and Franx et al. (2008) showed the
same trend in distant quenched galaxies. The corresponding
offset ridgeline of quenched galaxies to smaller R_e at fixed M*
was shown for SDSS galaxies by Shen et al. (2003), Omand
et al. (2014), and Haines et al. (2017) and for distant galaxies
by van der Wel et al. (2009, hereafter vdW09), van der Wel
et al. (2014, hereafter vdW14), and Haines et al. (2017).

Kauffmann et al. (2003, 2006, 2012) further noted that the
boundary between quenched and star-forming galaxies corre-
sponds closely to a line of constant effective surface density:
R<sub>e</sub> ∼ M<sub>*</sub><sup>0.5</sup>. Separately, Bell (2008), Kauffmann et al. (2012),
and Bluck et al. (2014) showed that local quenched galaxies
have high Sérsic indices characteristic of bulges, and the same
was shown for distant galaxies by Bell et al. (2012) and Wuyts
et al. (2012). Collectively, these studies established that
quenched galaxies are both smaller (more “compact”) than
star-forming galaxies at the same stellar mass and also that
quenched galaxies have prominent bulges.

Cheung et al. (2012) introduced the use of central properties to predict quenching. After testing numerous quantities, they concluded that high stellar-mass surface density within a projected radius of 1 kpc (Σ₁) is the best predictor of quenching at z ~ 0.7. Fang et al. (2013, hereafter F13) found a similar correlation between Σ₁ and stellar mass for SDSS galaxies, with quenched and green-valley galaxies occupying a narrow sloping ridgeline at the top of the Σ₁–M* distribution. With an rms scatter of only 0.16 dex (F13), the Σ₁–M* ridgeline is among the tightest scaling relations known for quenched galaxies. F13 further showed that SSFR declines sharply near this ridgeline and suggested that the residual ΔΣ₁ could be a marker for quenching, possibly signaling a threshold in BH feedback. Wake et al. (2012) (also Bluck et al. 2016; Teimoorinia et al. 2016) found a similar result for SDSS galaxies using central velocity dispersion σ₁ within 1 kpc (which F13 showed correlates closely with Σ₁). Finally, Luo et al. (2020) showed that high ΔΣ₁ correlates with high bulge-to-total ratio (B/T), consistent with the known quenching-spheroid connection.

van Dokkum et al. (2014) demonstrated that the threshold
relationship between Σ₁ and the specific star formation rate
(SSFR) found by F13 holds throughout cosmic time for high-
mass galaxies. Lower-mass galaxies were added by Tacchella
et al. (2015), Barro et al. (2017, hereafter B17), Whitaker et al.
(2017), Lee et al. (2018), and Mosleh et al. (2018), who found
that the Σ₁–M* relations for quenched and star-forming
galaxies at z = 0 are replicated almost unchanged out to
z = 3—both slopes are constant with time, and both zero
points march downward in unison by ~0.3 dex from z = 3 to
now.<fnref n="32" />

Various quenching models have tried to reproduce the
structural differences between quenched and star-forming
galaxies. Van der Wel et al. (2009, hereafter vdW09)
introduced a model in which galaxies evolve in R<sub>e</sub>–M* and
quench when they reach a boundary set by their effective
velocity dispersion, σ<sub>e</sub>. Because quenched galaxies with
higher σ<sub>e</sub> today are older at fixed M* (Bernardi et al. 2005;

<footnote n="32">³² We note parenthetically that, while Σ₁ and Σ_e are both surface densities, their behaviors versus mass and time are different (B17). The outer envelopes of quenched galaxies grow due to dry mergers after quenching, which increase R_e and decrease Σ_e while leaving Σ₁ roughly unchanged (e.g., Daddi et al. 2005; Trujillo et al. 2006; van Dokkum et al. 2008, 2010; Bezanson et al. 2009; Naab et al. 2009; Oser et al. 2012; Szomoru et al. 2012). The quenched-galaxy relation is therefore broader in Σ_e–M* than in Σ₁–M*, and its zero point declines with time, gradually encroaching on the domain of star-forming galaxies. The division between star-forming and quenched galaxies is therefore not as clean in Σ_e–M* as in Σ₁–M* (B17). Σ₁ and Σ_e are also different for star-forming galaxies, and the difference makes Σ₁ more suitable for predicting BH mass. We show below that Σ₁ ~ M*/R_e (Section 3.1), whereas Σ_e ~ M*/R_e². Because BH mass scales approximately as σ_e⁴ (e.g., KH13, van den Bosch 2016) and σ_e² ~ (M*/R_e)² ~ Σᵢ, a close power-law relation is predicted between M_BH and Σᵢ, as borne out below. In contrast, no simple power-law relation is predicted between M_BH and Σ_e.</footnote>
Graves et al. 2009), vdW09 reasoned that the threshold $\sigma_e$ for quenching must have been higher in the past. These facts motivated them to adopt a moving boundary that shifts upwards with time in $R_e-M_*$ (and downwards in $\sigma_e-M_*$). This moving boundary was actually detected by Haines et al. (2017) in $R_e-M_*$. However, neither of these works specified how galaxies approach the boundary before quenching. This was added by van Dokkum et al. (2015, hereafter vD15), who envisioned that galaxies evolve on parallel tracks through $R_e-M_*$ while star-forming and then quench when meeting a moving boundary in $\sigma_e$. A similar picture for the motion of galaxies in $R_e-M_*$ and the origin of quenched galaxies was presented by Cappellari (2013, 2016).

B17 and Tacchella et al. (2015) considered the evolution of galaxies in $\Sigma_1-M_*$. They concluded from the narrowness of the star-forming locus that galaxies must evolve approximately along it while star-forming. In their picture, galaxies reach a threshold value of $\Sigma_1$ for their mass, quench, and then pile up along the $\Sigma_1-M_*$ ridgeline. Because the ridgeline in $\Sigma_1-M_*$ moves down with time, this quenching boundary also moves down. Thus, like vdW09 and vD15, B17 envisioned a moving boundary but modeled it as a power law in $\Sigma_1-M_*$, not as a line at constant $\sigma_e$. None of these papers offered any physical explanation for their boundaries—they were simply set empirically to match the data.

Omand et al. (2014) introduced yet another aspect of quenching by focusing on $R_e-M_*$. They observed that star-forming SDSS galaxies populate a shallow ridgeline in $R_e-M_*$ with rms scatter $\sigma_{R_e} \sim 0.20$ dex at fixed mass, while quenched galaxies occupy a narrower and more steeply sloped ridgeline at the bottom edge of this distribution (Omand et al. 2014; see also Shen et al. 2003; vdW09; vdW14; Lange et al. 2016). Contours of equal quenched fraction are sloping stripes running parallel to the quenched ridgeline, and quenched galaxies therefore coexist with star-forming galaxies at the same $M_*$. This fact had been noticed before (e.g., Kauffmann et al. 2003), but Omand et al. seized on it as evidence that quenching is controlled by at least one more parameter other than mass and suggested that galaxies are a two-parameter family whose quenching is governed by $R_e$ as well as $M_*$.

Further evidence that both $R_e$ and $M_*$ play a role came from the work of van den Bosch (2016) on local BHs. Starting with the $M_{BH}-\sigma_e$ relation, he noted that $\sigma_e^2 \sim M_*/R_e$ from the virial theorem and proposed that BHs populate a “fundamental plane” in $M_{BH}-R_e-M_*$. In words, $R_e$ can be thought of as cleaning up the scatter in $M_{BH}$ versus $M_*$ that appears when star-forming galaxies are added to quenched galaxies (KH13) in the sense that smaller-radii galaxies must have larger BH. Similar arguments were made by Krajnović et al. (2018). Both works stressed that $R_e$ versus $M_*$ is a map of BH mass, an insight that will prove basic to our model.

Finally, Lilly & Carollo (2016, hereafter L16) presented yet another quenching model. Galaxies evolve along parallel log-log tracks in $R_e-M_*$, as in vdW14 and vD15, but quenching depends on $M_*$ probabilistically, increasing with $M_*$ until all massive galaxies eventually quench. Given that earlier-quenching galaxies tend to be smaller and that effective radii tend to shrink a lot on quenching in their model, L16 showed that the contours of equal quenched fractions would actually be tilted despite the fact that quenching formally depends on $M_*$ only. This model is discussed further in Section 7.9.

## 1.2. Collected Properties of Quenched versus Star-forming Galaxies

The previous section has introduced several different types of data that a successful quenching model must explain. For convenience, these and other data yet to be mentioned are collected together with basic references in Table 1. The first two lines of this table describe the green valley, which in this paper is taken to be the main channel through which high-mass central galaxies evolve from star-forming to quenched. Data on the mass function of quenched galaxies indicate that the green valley consists of two separate components: a high-mass channel that appears early at $z \ge 2.5$ with a roughly log-normal mass function, and a later channel after $z \sim 1.5$ that appears as a gradually rising low-mass tail (e.g., Ilbert et al. 2010; Huang et al. 2013; Muzzin et al. 2013).<fnref n="33" /> This tail is generally due to satellite quenching, not AGN quenching (Tinker et al. 2013; Fossati et al. 2017; Guo et al. 2017) and is therefore not treated in this paper. The high-mass channel is the object of our model, and its two striking features are the near constancy of its peak mass near log $M_*/M_\odot \simeq 10.7$ from $z=3$ to now<fnref n="34" /> and its substantial width, approaching 1 full dex in stellar mass. Line 3 of Table 1 also notes the continuing steady flow of galaxies though the high-mass channel down to recent times (e.g., Straatman et al. 2015, F18).

Lines 4–7 in Table 1 describe the properties of quenched and star-forming ridgelines in the $R_e-M_*$ plane. Star-forming galaxies need accurate treatment in the model because they feed the quenching process and their properties determine what quenched galaxies ultimately will look like. Important facts about the $R_e-M_*$ plane are the constancy of the ridgeline slopes for both quenched and star-forming galaxies over time and the increase in both ridgeline zero points since $z=3$. Lines 8–11 summarize similar properties of quenched and star-forming ridgelines in $\Sigma_1-M_*$. Important facts are again the constancy of both ridgeline slopes over time, accompanied (in this case) by decreasing zero points.

Line 12 of Table 1 says that the ridgeline zero point of quenched galaxies in $M_{BH}-M_*$ also decreases with time. Because less is known about the history of this ridgeline, we do not cite data on its slope back in time, but its general behavior seems similar to the motion of the quenched ridgeline in $\Sigma_1-M_*$. The distribution in $M_{BH}-M_*$ of star-forming AGN galaxies has been studied locally by Ho & Kim (2014) and Reines & Volonteri (2015) and in distant galaxies by Ding et al. (2017, 2020). This ridgeline seems to move down by a similar amount to the quenched ridgeline, which keeps the two populations in lockstep. To the extent that AGNs can be taken to represent all star-forming galaxies, this suggests that the entire $M_{BH}-M_*$ distribution for all galaxies maintains the same shape and moves down in bulk, just like that in $\Sigma_1-M_*$. Line 13 points out that there seem to exist various scaling laws linking $M_{BH}$ to both $M_*/R_e$ and $\sigma_e$ and that the slopes of these laws differ depending on whether all galaxies or only quenched ridgeline galaxies are included (KH13; van den Bosch 2016).

<footnote n="33">³³ Note that the earlier-cited increase in the number of massive quenched galaxies with time referred to the high-mass channel near the knee of the quenched mass function, not the low-mass tail.</footnote>

<footnote n="34">³⁴ Bundy et al. (2006) and Haines et al. (2017) both claim to see a decrease in the characteristic quenching mass with time. However, their definition of quenching mass is sensitive to the relative numbers of quenched and star-forming galaxies and tends to fall as galaxies quench. We focus on the typical green-valley mass itself, which is constant, as seen in, e.g., B17, Straatman et al. (2016), Haines et al. (2017), Pan et al. (2016), and Fang et al. (2012).</footnote>
<table><caption>**Table 1**
Observations to be Explained<sup>a</sup></caption><thead><tr><th>Item No.<sup>b</sup></th><th>Observation</th><th>References</th></tr></thead><tbody><tr><td>1</td><td>The stellar-mass distribution in the high-mass green valley peaks at approximately the same value ~10<sup>10.7</sup> M<sub>☉</sub> from z = 3 to now.</td><td>Ilbert et al. (2010)<br>Huang et al. (2013)<br>Barro et al. (2017)</td></tr><tr><td>2</td><td>The high-mass green valley is ~1 dex wide in stellar mass from z = 3 to now.</td><td>Ilbert et al. (2010)<br>Huang et al. (2013)<br>Barro et al. (2017)</td></tr><tr><td>3</td><td>There is a steady flow of high-mass central galaxies across the quenching boundary from z = 3 to now.</td><td>Bell et al. (2004)<br>Faber et al. (2007)<br>Tinker et al. (2013)</td></tr><tr><td>4</td><td>The ridgeline of star-forming galaxies in R<sub>e</sub> vs. M<sub>*</sub> has constant slope R<sub>e</sub> ~ M<sub>*</sub><sup>0.21</sup> from z = 3 to now.</td><td>Tomczak et al. (2014)</td></tr><tr><td>5</td><td>The zero point of the star-forming ridgeline in R<sub>e</sub> vs. M<sub>*</sub> has increased by ~0.4 dex from z = 3 to now.</td><td>van der Wel et al. (2014)</td></tr><tr><td>[6]</td><td>The ridgeline of quenched galaxies in R<sub>e</sub> vs. M<sub>*</sub> has constant slope R<sub>e</sub> ~ M<sub>*</sub><sup>0.75</sup> from z = 3 to now.</td><td>van der Wel et al. (2014)</td></tr><tr><td>[7]</td><td>The zero point of the quenched ridgeline in R<sub>e</sub> vs. M<sub>*</sub> has increased by ~0.75 dex from z = 3 to now.</td><td>van der Wel et al. (2014)</td></tr><tr><td>8</td><td>The ridgeline of star-forming galaxies in Σ<sub>1</sub> vs. M<sub>*</sub> has constant slope Σ<sub>1</sub> ~ M<sub>*</sub><sup>0.88</sup> from z = 2.5 to 0.5.</td><td>Barro et al. (2017)</td></tr><tr><td>9</td><td>The zero point of the star-forming ridgeline in Σ<sub>1</sub> vs. M<sub>*</sub> has decreased by ~0.30 dex from z = 2.5 to 0.5.</td><td>Barro et al. (2017)</td></tr><tr><td>[10]</td><td>The ridgeline of quenched galaxies in Σ<sub>1</sub> vs. M<sub>*</sub> has constant slope Σ<sub>1</sub> ~ M<sub>*</sub><sup>0.66</sup> from z = 2.5 to 0.5.</td><td>Barro et al. (2017)</td></tr><tr><td>[11]</td><td>The zero point of the quenched ridgeline in Σ<sub>1</sub> vs. M<sub>*</sub> has decreased by ~0.30 dex from z = 2.5 to 0.5.</td><td>Barro et al. (2017)</td></tr><tr><td>[12]</td><td>The zero points of the quenched and AGN galaxy ridgelines in M<sub>BH</sub> vs. M<sub>*</sub> appear to be decreasing with time by similar amounts.</td><td>Bennert et al. (2011)<br>Ding et al. (2017)</td></tr><tr><td>13</td><td>M<sub>BH</sub> can be expressed as a power law of either M<sub>*</sub>/R<sub>e</sub> or σ<sub>e</sub>, but the slopes of quenched galaxies are shallower than the slopes of all galaxies taken together.</td><td>Kormendy & Ho (2013)<br>van den Bosch (2016)</td></tr><tr><td>14</td><td>The ridgelines of quenched galaxies are narrow in Σ<sub>1</sub> vs. M<sub>*</sub> and M<sub>BH</sub> vs. M<sub>*</sub> but not so narrow in R<sub>e</sub> vs. M<sub>*</sub>.</td><td>Fang et al. (2013)<br>Barro et al. (2017)<br>Kormendy & Ho (2013)</td></tr><tr><td>15</td><td>The quenched ridgelines in Σ<sub>1</sub> vs. M<sub>*</sub> and M<sub>BH</sub> vs. M<sub>*</sub> show age stratification, with older galaxies having higher values of Σ<sub>1</sub> and M<sub>BH</sub> at fixed M<sub>*</sub>.</td><td>van der Wel et al. (2014)<br>Tacchella et al. (2017)</td></tr><tr><td>[16]</td><td>The zero points of quenched and star-forming galaxies have evolved together in lockstep in Σ<sub>1</sub> vs. M<sub>*</sub> since z = 3.</td><td>Luo et al. (2020)<br>Barro et al. (2017)</td></tr><tr><td>[17]</td><td>The zero points of quenched and star-forming galaxies have evolved differently in R<sub>e</sub> vs. M<sub>*</sub> since z = 3.</td><td>van der Wel et al. (2014)</td></tr><tr><td>18</td><td>Star-forming galaxies mix with quenched galaxies on the quenched ridgeline in Σ<sub>1</sub> vs. M<sub>*</sub> today but lie well below the quenched ridgeline in M<sub>BH</sub> vs. M<sub>*</sub>.</td><td>Fang et al. (2013)<br>Barro et al. (2017)<br>Terrazas et al. (2017)</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Referring to central galaxies only.

<sup>b</sup> Items without brackets are thought to be fully explained by the model; items with brackets are only partly explained. The model being scored is the binding-energy model (blue lines in Figure 1), for which there is only one adjustable parameter, $f_B = 4$ at entry into the green valley, which has been used to match the zero point of the blue lines to that of the red lines. All of the other successes listed here are new. See further discussion in Section 7.1.

Line 14 refers to the width of the quenched ridgelines, which are quite narrow in Σ<sub>1</sub>–M<sub>*</sub> and M<sub>BH</sub>–M<sub>*</sub> at fixed M<sub>*</sub> but not so narrow in R<sub>e</sub>–M<sub>*</sub>, especially at early times. Our model says that these quenched ridgelines are not excesses due to sighting along a plane edge on (like the fundamental plane for ellipticals or the Tully–Fisher relation for spirals) but are actual pileups of galaxies that stop evolving and come to a halt along particular loci. The slopes and zero points of these loci therefore contain important information about quenching, which is basic data for the model. However, the quenched ridgeline in R<sub>e</sub>–M<sub>*</sub> is somewhat broader than the other two due to variable amounts of galaxywide dissipation during quenching (e.g., Zolotov et al. 2015; Tacchella et al. 2016a) and to variable amounts of dry merging after quenching (e.g., Naab et al. 2009; vdW09), both of which tend to decouple global R<sub>e</sub> from central Σ<sub>1</sub> and M<sub>BH</sub>. The narrowness of the quenched ridgelines in Σ<sub>1</sub>–M<sub>*</sub> and M<sub>BH</sub>–M<sub>*</sub> therefore suggests that those relations more closely preserve the memory of quenching; significantly, both are properties that relate to the centers of galaxies, not to their outer

parts. However, line 15 points out that recent data are beginning to reveal age stratification even within both of these ridgelines, with the quenched stellar populations above the ridgelines being older (Tacchella et al. 2017; Martín-Navarro et al. 2018; Luo et al. 2020). The authors attribute this to “progenitor bias” (van Dokkum & Franx 1996; Carollo et al. 2013; Poggianti et al. 2013) caused by the moving boundary, which evidently leaves early-quenched objects stranded above the boundary at late times.

Lines 16 and 17 highlight the different zero-point behaviors of quenched versus star-forming galaxies in the different spaces. In Σ<sub>1</sub>–M<sub>*</sub>, the ridgelines of both populations move down together in lockstep so that the relative positions of quenched versus star-forming galaxies are preserved and the Σ<sub>1</sub>–M<sub>*</sub> diagram looks the same at all times except for this bulk downward shift (see Figure 2 in B17). However, in R<sub>e</sub>–M<sub>*</sub>, the ridgeline of quenched galaxies moves up faster than the ridgeline of star-forming galaxies and the two populations begin to overlap, as in vdW14.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Data from all five CANDELS fields divided by redshift, patterned after similar figures in vdW09 and B17. The data and sample are described in Section 2; only face-on galaxies with $b/a \ge 0.5$ are used. Satellites and centrals are included even though the model technically refers to centrals only. Panel (a) (top row): effective radius $R_e$ along the semimajor axis (SMA) vs. stellar mass, colored by $\Sigma_1$ for star-forming galaxies. Panel (b) (second row): the same as panel (a) but for quenched galaxies. Panel (c) (third row): $\Sigma_1$ vs. stellar mass colored by $R_e$ for star-forming galaxies. Panel (d) (bottom row): same as panel (c) but for quenched galaxies. Star-forming galaxies are defined as lying above the line running -0.45 dex below the main-sequence ridgeline, and quenched galaxies lie more than -1.0 dex below it. Green-valley galaxies (between these two limits) are not shown. The figure illustrates the dramatic difference in structure between star-forming and quenched galaxies and the sharp boundaries between them. The red lines are the empirical quenching boundaries (tilted boundaries) that are developed from the data in Section 3. The blue lines are predicted by the halo gas-binding-energy model in Section 5, which says that entry to the green valley occurs when the total emitted BH energy ($E_{BH}$) equals $4 \times$ the binding energy of hot halo gas. The agreement between the red and blue lines is the major result of this paper.</figcaption></figure>

The first column in Table 1 numbers the items for future reference. Item numbers without brackets indicate phenomena that we think are fairly well explained by the new model, while item numbers with brackets indicate phenomena that are only partially explained. We stress again that the model only predicts $R_e$ for star-forming galaxies before they reach the quenching boundary and that $R_e$ is not modeled later, when dissipation and dry merging become important. Because lines 6, 7, and 17 all involve quenched galaxies, their values of $R_e$ are not accurately known, which is why they are in brackets. Other brackets are described in Section 7.1.

Figure 1 illustrates many of the properties listed in Table 1. This figure borrows heavily from similar figures in vdW14

and B17 using CANDELS data. The upper two panels compare star-forming and quenched galaxies in the space of $R_e$ versus $M_*$, and the lower two panels compare the same samples in the space of $\Sigma_1$ versus $M_*$ (in order to sharpen differences, green-valley galaxies are not plotted). Looking ahead, the red lines represent empirically determined quenching boundaries that are derived in Section 3, and the blue lines represent the quenching boundaries of the “binding-energy” model of this paper, in which the threshold quenching BH energy is equal to a multiple of the binding energy of halo gas.

Focusing first on the lower panels with $\Sigma_1$, we see that virtually all quenched galaxies lie on a ridgeline above the red quenching boundary and their numbers are growing with time,
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Schematic illustration of how star-forming galaxies evolve in the empirical power-law model. Galaxies evolve along approximately power-law tracks in a four-dimensional space comprising $M_*$, $R_e$, $\Sigma_1$, and $M_{BH}$. The panels show projections onto $M_*$ of each of the other three quantities. The arrows show the evolutionary tracks of individual galaxies, and the sloping (tilted) dashed lines show the quenching boundaries at which galaxies enter the green valley. The fact that the dashed lines are not vertical indicates that at least one other parameter beyond $M_*$ is involved in quenching. Panel (a): $R_e$ vs. $M_*$. The points are face-on CANDELS galaxies with $b/a > 0.5$ in the redshift range $z = 0.6–0.9$. Blue points are star-forming galaxies on the star-forming main sequence with $\Delta\log SSFR$ above –0.45 dex; red points are fully quenched galaxies with $\Delta\log SSFR$ below –1.0 dex. Green valley galaxies are not shown. Panel (b): $\Sigma_1$ vs. $M_*$. The objects are the same as in panel (a). Panel (c): $M_{BH}$ vs. $M_*$. This diagram reproduces Figure 2 from Terrazas et al. (2016), which plots BH masses for 90 nearby galaxies. Blue points are star-forming galaxies with the same selection criteria as in Terrazas et al. (2016); red points are now both quenched and green-valley galaxies. Note the strong bias to high masses in this panel compared to panels (a) and (b). All three panels of the figure are schematic: the slopes of the dashed boundary lines are accurate, but the evolutionary track slopes are illustrative only, and the dashed boundaries move with time.</figcaption></figure>

while the vast majority of star-forming galaxies lie below the same line. Galaxies approach the quenching boundary from the bottom left in this plane, as shown in Figure 2 below. A smattering of star-forming galaxies also lies above the line, which is the so-called “elbow effect” pointed out by F13 and B17 in which star-forming and quenched galaxies coexist at nearly the same values of $\Sigma_1$ (see also Whitaker et al. 2017; Lee et al. 2018). This is mentioned in line 18 of Table 1 and signifies a highly nonlinear relation between $\Sigma_1$ and the decline of star formation in the green valley. Line 18 also says that the downward trend in star formation rate versus $M_{BH}$ is more gradual (straighter; see Figure 3, Terrazas et al. 2017) so that star-forming galaxies never lie on the quenched ridgeline in $M_{BH}-M_*$ but always lie well below it (Figure 2, Terrazas et al. 2016). These different behaviors of $\Sigma_1$ versus $M_*$ and $M_{BH}$ versus $M_*$ near quenching are illustrated in Figures 2 and 5 below, and matching them will prove to be a major clue to BH growth rates in the green valley. Note that the quenched ridgelines, and hence also the red lines, maintain constant slope with time in $\Sigma_1-M_*$ but move down in zero point. The ridgelines of star-forming galaxies also maintain constant slope and move down in lockstep with the quenched ridgeline in the same diagram. Collectively, these are lines 8–11, 14, and 16 in Table 1.

The upper two panels in Figure 1 show the same objects in $R_e-M_*$, where the pattern is reversed—virtually all quenched galaxies now fall below the red quenching boundary while the vast majority of star-forming galaxies fall above it. The picture here is that galaxies approach the quenching boundary from above left, as shown in Figure 2. Note that the quenched ridgelines, and hence also the red lines, again maintain constant slope but now move up in zero point with time. The ridgelines of the star-forming galaxies also maintain constant slope and move up. However, in this case, the zero point of the quenched galaxies moves up faster than the zero point of the star-forming galaxies, and the two populations are not in lockstep the way they are in $\Sigma_1-M_*$. Collectively, these are properties 4–7, 14, and 17 in Table 1.

Three more features of Figure 1 are important. First, one sees that quenched galaxies cover roughly a full dex in mass that is centered at about $10^{10.7} M_\odot$ at all redshifts. These are lines 1 and 2 in Table 1. Second, there is a dearth of quenched galaxies near and below $10^{10} M_\odot$ at $z > 2$ that is gradually filled in at lower redshift. This is the growing tail of low-mass quenched galaxies that was referred to above and that we attribute to satellite quenching; it is not treated in the current paper. The third point is that the trends seen Figure 1 are ancient, going all the way back to $z \sim 3$. The enduring nature of these relations seems to imply that the quenching mechanism for high-mass galaxies above $10^{10} M_\odot$ has remained basically the same over all these billions of years despite major declines in star formation rates, gas contents, and merger rates. To our knowledge, no SAM or hydro model has as yet reproduced the full collection of properties listed in Table 1 and illustrated in Figure 1.

### 1.3. Outline for the Remainder of the Paper

The first step in this paper, in Section 3, is an empirical parameterization of the quenching boundary that is fitted to real data. Specifically, we assume that the effects of quenching become apparent on galaxy star formation rates when the total effective energy emitted by their BH equals some total energy quota needed by the halo to stop cooling: $E_{BH} = E_{quench}$, and we derive an empirical power-law expression for $E_{quench}$ as a function of $M_{vir}$ and redshift from the observed (red) boundaries in Figure 1. Starting with this empirical formulation for $E_{quench}$ rather than a theoretical expression provides a simple form that can be experimentally manipulated and also shows intuitively how the real halo energy quota must scale if quenching is in fact a balance between cumulative BH energy output and halo mass.

The second major step, in Section 5, compares this empirical halo energy quota to the gas-binding energies of real halos. This comparison is shown by the blue lines in Figure 1, which are seen to agree fairly well with the red lines with no
adjustment in slope and a plausible one-time choice of zero
point. The blue lines are constructed by assuming that
quenching sets in when the total effective energy emitted by
the BH equals some constant, $f_B$, times the binding-energy
content of the halo gas. The precise quenching mechanism is
unspecified—it could be gas heating and ejection, a physical
deflection of infalling primordial gas, or actual gas unbinding
—and the cumulative energy input from BHs is in principle
large enough to accomplish any of these things. Equating the
onset of quenching to the tipping point between total BH
energy and halo gas-binding energy has also been proposed by
others (e.g., Silk & Rees 1998; Ostriker et al. 2005; Booth &
Schaye 2010, 2011; Bower et al. 2017; Oppenheimer 2018;
Davies et al. 2019, 2020; Oppenheimer et al. 2020). Our major
new result is that this simple model explains the observed
quenching boundaries in Figure 1 quite well, along with many
other data in Table 1.

Four other aspects of the model deserve mention. First, our
knowledge of halos and the halo–galaxy connection is used to
model how halos and their central galaxies evolve together
before quenching. A second key assumption is the rule for
putting BHs into star-forming galaxies, which assumes that BH
mass scales as a power law of $\Sigma_1$ during the star-forming phase,
which is suggested by empirical data. This rule implies that
BHs grow modestly but steadily (averaging over possible short-
term high-accretion episodes) while galaxies form stars.<fnref n="35" /> The
third feature replaces this early, modest power law with a
steeper power law (versus $\Sigma_1$) in the green valley (see
Figure 4). This also agrees with a variety of data and can
also reconcile the different values of BH scaling law slopes
found in different works. The fourth and final feature,
following Omand et al. (2014), van den Bosch (2016), and
Krajnović et al. (2018), is to model star-forming galaxies as a
two-parameter family labeled by $R_e$ and $M_*$. This says that BH
mass depends on galaxy radius as well as $M_*$—larger galaxies
have smaller BHs because they have lower central $\Sigma_1$—which
means they have to evolve to higher mass before they
overcome their halos and quench. This last is the key factor
that explains the sloped quenching boundaries seen in Figures 1
and 2.

Section 6 explores four variants of the basic model. All fail,
but each one teaches an important lesson, as discussed in
Appendix D. Section 7 is the Discussion. Topics include the
model’s success in matching the data goals in Table 1, the
origin of scatter in effective radii, how quenching happens in
the binding-energy model compared to other models, compac-
tion, the connection between BHs and dark halos, BH scaling
laws, and what it means to say that BHs and galaxies
“coevolve.” Finally, Section 8 reviews the basic findings and
important takeaway messages, and four appendices present
various formulae, summarize important power laws, and
describe the variant models in more detail.

We end this Introduction with a list of limitations and
phenomena that are not covered by the model. All scaling laws
are assumed to be perfectly tight with zero scatter. Mergers and
satellite quenching are ignored. It is assumed that entry into the
green valley is a one-way street and that rejuvenation does not
occur (see Chauke et al. 2019). Because of its use of $\Sigma_1$ to track

BH growth, the model cannot be used before $z \sim 3$, where
many galaxies fit entirely within 1 kpc. Likewise, no descrip-
tion is provided for the microscopic state of gas in the halo or
how the BH energy interacts with it. Finally, the data matched
are mainly global structural relations, and no attempt is made to
study the pattern of star formation or quenching within
galaxies, to predict mass functions or other structure functions,
or to predict AGN frequencies. These topics are left for
future work.

All data used in this paper are consistent with $H_0 = 70 \text{ km s}^{-1} \text{ Mpc}^{-1}$, $\Omega_m = 0.3$, and $\Omega_\Lambda = 0.7$.

## 2. Data

The data used in this paper come from the multiwavelength
and ancillary data sets produced by the Cosmic Assembly
Near-infrared Deep Extragalactic Legacy Survey on the
Hubble Space Telescope (CANDELS; Grogin et al. 2011;
Koekemoer et al. 2011). All five fields are used (GOODS-S,
UDS, COSMOS, EGS, and GOODS-N) with $H \le 24.5$ mag.
AB magnitudes, spectroscopic/photometric redshifts, stellar
masses, star formation rates, and $A_V$ attenuation values are
taken from the following catalogs: Guo et al. (2013), Galametz
et al. (2013), Nayyeri et al. (2017), Stefanon et al. (2017),
and Barro et al. (2019). Methodologies for the stellar masses and
$A_V$ values are described in Santini et al. (2015). $R_e$ values
and axis ratios are measured using GALFIT by van der Wel
et al. (2012). The selection and completeness of the sample are
discussed in Appendix A.

An important parameter in this paper is the central stellar
density, $\Sigma_1$, whose measurement is explained in B17. Briefly,
mass profiles were constructed from multiband Hubble Space
Telescope (HST) images, and raw values of $\Sigma_1$ were measured
based on the projected stellar mass within 1 kpc. A PSF
correction was then applied from a grid of high-resolution
model images having the same values of $R_e$ and Sérsic index.
Details are given in Appendix C of B17. Two separate sets of
values were computed using slightly different techniques by
Barro and Chen, and they agree well. All data including $\Sigma_1$ are
available at the Rainbow database (Barro et al. 2011).<fnref n="36" /> The
values of $\Sigma_1$ used here are the Chen values.

CANDELS galaxies have not yet been divided into centrals
and satellites even though the model is nominally for central
galaxies only. Environments for CANDELS galaxies are being
measured, and this distinction is left for future work.

# 3. Power Laws in the Empirical Model

## 3.1. Input Power Laws

This section and the next construct an empirical model that
describes the evolutionary tracks that galaxies take while star-
forming, and then quenches them as they pass over the
quenching boundary. Several key relations are power laws. The
present section derives their slopes, and the next section derives
their zero points and how they move with time. Full versions of
all laws, including the zero points, are given in Appendix C.

Figure 2 plots three structural diagrams that illustrate the
basic landscape of the model. Star-forming galaxies evolve
along straight, power-law trajectories in the four-dimensional
space comprising $M_*$, $R_e$, $\Sigma_1$, and $M_{BH}$. The three panels in
Figure 2 show the three two-dimensional projections of this

<footnote n="35">³⁵ BH growth as a power of $\Sigma_1$ during the star-forming phase is the default model for most of the paper. However, the Discussion section considers a complementary picture in which BHs remain smaller than this throughout most of the star-forming phase and only begin to grow rapidly later (see Section 7.8).</footnote>

<footnote n="36">³⁶ http://arcoiris.ucsc.edu//Rainbow_navigator_public/</footnote>
four-dimensional space onto $M_*$.
Eventually, galaxies encounter a boundary in this four-dimensional space, as shown by the dashed lines, where they start to quench. Crossing the boundary is equated with entering the green valley, i.e., the point where AGN feedback has a just-perceptible effect on the star formation rate. The zero point of the boundary is defined using a well-defined feature in the $\Sigma_1-M_*$ distribution of SDSS galaxies from Luo et al. (2020) and is then mapped into the other two projections using the power-law mappings developed in this section. The blue points are star-forming galaxies, and the quenching boundaries mark their maximum extent along the arrows. In panels (a) and (b), the red points are fully quenched galaxies only, and in panel (c), they are both quenched and green-valley galaxies. Panels (a) and (b) show CANDELS galaxies while panel (c) shows a selection of local galaxies with measured BHs from Terrazas et al. (2016).

The use of power laws for boundaries and scaling laws is motivated by the fact that the ridgelines of quenched galaxies are approximately power laws. As Section 6 shows, this all but requires that the equations of quenching also be power laws in basic galaxy/halo properties, which places important constraints on quenching physics (see Section 6).

The power laws used in this paper are in several categories:

(1) A structural mapping relates one set of structural variables to another. The main example is the mapping $\Sigma_1 \sim M_*^g R_e^{-p}$, which maps the $R_e-M_*$ plane onto the $\Sigma_1-M_*$ plane. The variables p and q depend on the Sérsic index n (B17), and hence a given p, q mapping applies only to galaxies with similar Sérsic indices. The needed mapping is for star-forming galaxies, which generally have n in the range 1–2 (B17).

(2) A quenched ridgeline is the fit at a given redshift to the locus of galaxies that have quenched. It is narrowest and best defined in $\Sigma_1-M_*$, where it has the slope 0.66 at all redshifts (B17). Quenched galaxies, as we have noted, are defined as objects that lie more than 1.0 dex below the star-forming main sequence at their mass and redshift.

(3) A star-forming ridgeline is the fit at a given redshift to the locus of galaxies that are still star-forming. It is conventionally shown in the $R_e-M_*$ plane (e.g., vdW14), but it also exists in the $\Sigma_1-M_*$ plane, being mapped there by the p and q values for star-forming galaxies. Star-forming galaxies are defined as all objects that lie above the line lying $-0.45$ dex below the star-forming main sequence at their mass and redshift.

(4) A scaling law relates one variable to the power of another, usually with some implied physics. Examples are our assumed scaling laws that BH mass scales as $M_{BH} \sim \Sigma_1^\nu$ in star-forming galaxies and that the total effective energy radiated during the course of BH mass growth varies as $E_{BH} \sim \eta M_{BH} c^2$.

(5) An evolutionary track is a locus that represents the evolution of one variable versus another over time. An example is the stellar-mass/halo-mass (SMHM) relation, which we take to have the form $M_* \sim M_{vir}^s$. Because we assume no change in slope or zero point with redshift, it also is an evolutionary track. Hence, some scaling laws, like the SMHM relation, can also be evolutionary tracks.

(6) A quenching boundary is a power-law locus in the space of some set of structural variables at a given redshift that signals the beginning of the end of star formation when galaxies cross over it. There are three projected boundaries, in

$\Sigma_1-M_*$, $R_e-M_*$, and $M_{BH}-M_*$, which together form a single line in the master 4D space.

The following five power laws are combined to derive expressions for the slopes of the projected quenching boundaries. The first two of these are derived from abundance matching:

1. SMHM relation for star-forming galaxies: This evolutionary track is parameterized as

$$M_* = \text{const.} \times M_{\text{vir}}^s, \quad (1)$$

where $M_*$ is the total stellar mass and $M_{vir}$ is the halo virial mass. $M_*$ is defined as the current stellar mass including stellar remnants but does not include the mass lost in stellar evolution; the definitions of halo virial mass and virial radii are given in Section 4.1. In the basic version of the model that we are describing now, the SMHM relation is assumed to be redshift independent based on the nearly constant zero point obtained by Rodríguez-Puebla et al. (2017, hereafter RP17) from abundance matching (see Figure D1). Other works give slightly different results, and Variant 2 in Section 6 explores the effect of varying the SMHM zero point with redshift. The exponent s is set to 1.75, which is the slope below the knee from RP17. We note that this is very close to the value s ~ 5/3 predicted from the energetics of SN feedback versus halo potential well depth (Dekel & Silk 1986; Dekel & Woo 2003). The SMHM law below the knee applies to star-forming galaxies only, and because no scatter is assumed, it is simply referred to as the SMHM relation. Above the knee, the SMHM or its inverse, the halo-mass to stellar-mass relation (HMSM), is used as appropriate.<fnref n="37" /> The model does not use star formation rates explicitly, but using a particular SMHM relation for star-forming galaxies implicitly assumes a matching star formation rate to make it true.

2. Median galaxy effective radius, $\langle R_e \rangle$, versus halo virial radius, $R_{vir}$, for star-forming galaxies: This evolutionary track assumes that the median $\langle R_e \rangle$ of star-forming galaxies at a given stellar mass is the same fixed fraction of the halo virial radius, $R_{vir}$, at all masses and redshifts ($R_e$ is the observed half-light semimajor axis from van der Wel et al. (2012), called SMA in Figures 1–3). Thus,

$$\langle R_e \rangle = \text{const.} \times R_{\text{vir}}. \quad (2)$$

This assumption is consistent with the fact that the observed median radii of star-forming galaxies appear to equal about the same fixed fraction of model halo virial radii back to $z=3$ to within approximately 20% (Kravtsov 2013; Shibuya et al. 2015; Huang et al. 2017; RP17; Somerville et al. 2018).<fnref n="38" />

The next three power laws feature the central role of $\Sigma_1$:

<footnote n="37">³⁷ The knee occurs at approximately $M_* \sim 10^{10.7} M_\odot$, which is not quite massive enough to encompass all star-forming galaxies. This problem is addressed in Variant 1 in Section 6.</footnote>

<footnote n="38">³⁸ All observed radii in this paper are projected 2D semimajor-axis half-light radii. These can differ significantly and systematically from 3D half-mass stellar radii, as shown in Somerville et al. (2018). Equation (2) is obeyed more accurately by 3D half-mass radii, but the differences are unimportant given the approximate nature of our treatment.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 3.** $R_e$ and $\Sigma_1$ vs. $M_*$ from all five CANDELS fields for star-forming galaxies only. This figure differs from Figure 1 in showing only star-forming galaxies and to a lower mass limit. Galaxies satisfy $H \le 24.5$ mag and axis ratio $b/a \ge 0.5$, and are classed as star-forming using the same criterion as in Figure 1. Upper panel: effective radius $R_e$ along the semimajor axis (SMA) vs. stellar mass, colored by $\Sigma_1$. Lower panel: $\Sigma_1$ vs. stellar mass colored by $R_e$. The loci of constant $\Sigma_1$ (above) and constant $R_e$ (below) are visible as colored stripes. These indicate the presence of a plane in $R_e-\Sigma_1-M_*$ space for star-forming galaxies, and the two plots are different projections of it. $\Sigma_1 \sim M_*^q R_e^{-p} \sim M_*^{1.1} R_e^{-1}$ is the best mapping of this plane from one projection to another (Appendix A). The red points in the plots of SMA at high redshift are compact star-forming galaxies that have been identified as the precursors of compact quenched “red nuggets” at this redshift (Barro et al. 2013, 2014; Nelson et al. 2014; Williams et al. 2014; vD15).</figcaption></figure>

1. $\Sigma_1$ versus $R_e$ and $M_*$ for star-forming galaxies: This relation is an empirical structural mapping of the form

$$ \Sigma_1 = \text{const.} \times M_*^q R_e^{-p}, \quad (3) $$

where $\Sigma_1$ is the projected stellar mass within 1 kpc of the center of the galaxy. Figure 3(a) plots $R_e$ versus $M_*$ for CANDELS star-forming galaxies color coded by $\Sigma_1$, and Figure 3(b) plots $\Sigma_1$ versus $M_*$ color coded by $R_e$. The slanting colored stripes in both panels indicate the existence of a plane in $\Sigma_1-R_e-M_*$ space (a similar plane was shown for VIPERS galaxies at $z = 0.5-1.0$ by Haines et al. 2017). Fitted values for the exponents of this plane are $p = 1.00 \pm 0.01$, $q = 1.10 \pm 0.01$, which are essentially constant with redshift, as is the zero point. Comparisons of data to predictions from this relation are shown in Figure A1 in Appendix A, where rms residuals in $R_e$ and $\Sigma_1$ are both ~0.2 dex. This $p, q$ mapping applies to star-forming galaxies only and is valid because their Sérsic indices are all in the range $n = 1-2$ (B17). (The analogous mapping for quenched galaxies has $p = 0.50$, $q = 1.10$ and is less tight owing to the more variable Sérsic indices of quenched galaxies.)

2. $\Sigma_1$ versus central velocity dispersion, $\sigma_1$: F13 empirically demonstrated the scaling law that

$$ \Sigma_1 = \text{const.} \times \sigma_1^{2.0 \pm 0.2} \quad (4) $$

for SDSS galaxies, where $\sigma_1$ is the average velocity dispersion along the line of sight within 1 kpc. Because SDSS does not resolve rotational motion close to the centers of galaxies, $\sigma_1^2$ is a measure of the total kinetic energy at the center of the galaxy, including both rotation and dispersion. The relation emerges from the virial mass estimator $M_* \sim R\sigma^2/G$ applied to the mass inside a

sphere if we note that

$$ M_* \sim \Sigma_1 \times \pi R^2 $$

and evaluate all quantities at $R = 1$ kpc. It then follows that $\Sigma_1 \sim \sigma_1^2$ except for (evidently small) geometric factors due to light and mass profile differences, variable mass projected in front of and behind the central sphere, variable orbital geometries, and different inclinations. Equation (4) is assumed to be valid for star-forming and quenched galaxies with the same constant at all redshifts. The main use of Equation (4) is not for actual computations but to help motivate what comes next.

The last step is to specify rules for BH growth. The basic assumption is that BH mass scales as a power law of $\Sigma_1$, but with potentially different slopes and zero points for star-forming, quenched, and green-valley galaxies. The adopted laws are shown in Figure 4. Star-forming galaxies enter along the blue track at the bottom. At some point, they cross the quenching boundary and then turn off into the green valley along the steep (green) vectors. This happens at different values of $\Sigma_1$ for different-$R_e$ galaxies at different times (see below), which creates multiple tracks through the green-valley regime. Eventually BHs stop growing, and galaxies come to a halt along the fully quenched ridgeline (dashed red line). Further evolution is not included in the model.

The series of broken-power-law tracks in Figure 4 is motivated by the following logic: (1) BHs in local quenched galaxies obey a power-law relation $M_{BH} \sim \sigma_e^u$ (KH13). We assume $\sigma_1 \sim \sigma_e$, in which case $M_{BH} \sim \Sigma_1^v$, where $v = u/2.0$ from Item 4 above. The value of $v$ is the slope of the upper, quenched ridgeline in Figure 4. (2) A steep track in the green valley is then indicated by the fact that star-forming and quenched galaxies can have similar $\Sigma_1$ values today at the
<figure><img src="image_5.png" /><figcaption>**Figure 5.** Cartoon illustrating the different evolution of $\Sigma_1$ and $M_{BH}$ through the green valley. Panel (a): $\Delta \log \Sigma_1$ is the residual of $\log \Sigma_1$ with respect to the quenched ridgeline in $\Sigma_1-M_*$ and $\Delta \log \text{SSFR}$ is the residual of $\text{logSSFR}$ with respect to the star-forming main sequence. $\Sigma_1$ shows a characteristic “elbow” shape whereby SSFR falls sharply but $\Sigma_1$ grows only slightly (F13; B17) through the green valley. Panel (b): The analogous plot for $\Delta \log M_{BH}$ based on data from Terrazas et al. (2017, Figure 3). In contrast to $\Sigma_1$ in panel (a), BH mass grows considerably in the green valley. The delta quantities above the figures are used to plot the green-valley vectors in Figure 4.</figcaption></figure>

evidence is actually growing for two different laws with different zero points, one for quenched galaxies (KH13) and one for star-forming galaxies (Ho & Kim 2014; Reines & Volonteri 2015; Savorgnan et al. 2016), separated by the green valley. The slope of the quenched law is better determined, which is why we have used it as our starting point, but the basic idea of two laws is consistent with our scheme in Figure 4. More discussion of the adopted laws and their derivation is given in Section 7.

The preceding paragraphs have motivated the law that $M_{BH} = \text{const.} \times \Sigma_1^\nu$ for star-forming galaxies and that BHs grow more steeply relative to $\Sigma_1$ in the green valley. We now extrapolate those rules and apply them to star-forming and green-valley galaxies at all redshifts. The main argument in favor of this assumption is the similar appearance of the $R_e-M_*$, $\Sigma_1-M_*$, and $M_{BH}-M_*$ scaling laws at all redshifts, which suggests that the same basic physics has operated at all epochs since $z \sim 3$. Thus, if bulges build BHs today, we assume that they did so in the same way at early times, and we model the behavior of early galaxies on current galaxies.

Before leaving this topic, we note that the cartoon in Figure 4 interprets the blue BH scaling law for star-forming galaxies as an evolutionary track that galaxies evolve along. Indeed, that is how it was introduced above. Considerable evidence suggests that star-forming galaxies in fact spend most of their time below the blue relation and build their BHs only later as they enter the green valley. This means that they would cross the blue line from below when passing into the green valley. This alternate scenario is discussed in Section 7.8.

Returning to the main thread and combining the above logic, we now have:

1. BH mass versus $\Sigma_1$ for star-forming galaxies:

$$M_{\mathrm{BH}} = \mathrm{const.} \times \Sigma_{1}^{\nu}, \qquad (5)$$

where $u$ is the slope of the quenched ridgeline and $v = u/2.0$. This is the blue line in Figure 4, which is assumed

to hold for all star-forming galaxies at all redshifts $z < 3$. Because $M_{BH}$ is uniquely specified by $\Sigma_1$, this means that $M_{BH}-M_*$ maps onto $\Sigma_1-M_*$ and that $M_{BH}$ is now also joined to the space of $\Sigma_1-R_e-M_*$ for star-forming galaxies.

To this point, star-forming galaxies in the model have only one degree of freedom, which is $M_{vir}$ today (or alternatively $M_*$ today), and all galaxies of a given $M_*$ at a given redshift have the same $R_e = \langle R_e \rangle$. It remains now to introduce the crucial radius variations about $\langle R_e \rangle$ as a second parameter, without which there would be no variations in $M_{BH}$ at a given stellar mass. Measurable scatter is observed in $R_e$ at fixed $M_*$ for star-forming galaxies, with $\sigma_{R_e} \sim 0.20$ dex at all masses (vdW14). At this point, the introduction of scatter is purely empirical; a discussion of possible origins is deferred to Section 7.2.

Adding radius scatter makes star-forming galaxies a two-parameter family labeled by $M_*$ and $\Delta \log R_e$. Looking ahead to the model with zero points in Section 4, we find it convenient to reexpress these numbers in terms of the values for a galaxy and its halo at a starting redshift, $z_s$, when the stellar mass $M_* = 10^9 M_\odot$. Thus, the initial parameters at $z_s$ are $M_{*,s} = 10^9 M_\odot$, $M_{vir} = M_{vir,s}$, $R_e = R_{e,s}$, and $R_{vir} = R_{vir,s}$. How these parameters are set is described in Section 4.2.

To model the scatter, we multiply each starting radius by a scatter factor $R_{e,s}/\langle R_e \rangle$, where $\langle R_e \rangle$ is the mean value of $R_e$ for $10^9 M_\odot$ galaxies at $z_s$. This scatter factor is assumed to be constant throughout the star-forming phase, and its logarithm is assumed to be Gaussiany distributed with $\sigma_{R_e} = 0.2$ dex. Every galaxy is now labeled by two numbers, the redshift $z_s$ at which $M_* = 10^9 M_\odot$ and its value $\Delta \log R_{e,s}$ at that redshift. The latter is also equivalent to the scatter in $M_{BH}$ and $\Sigma_1$ because the previous equations imply

$$\Delta \log \Sigma_{1,s} = -p \Delta \log R_{e,s}, \qquad (6)$$

$$\Delta \log M_{BH,s} = v \Delta \log \Sigma_{1,s}, \qquad (7)$$
and

$$ \Delta \log M_{BH,s} = -vp \Delta \log R_{e,s}. \quad (8) $$

Note that galaxies with larger $R_e$ at fixed $M_*$ have smaller BHs.

The next step is to evolve galaxies forward in time. To do this, another crucial assumption is made, namely, that the star formation rate at a given mass is independent of $R_e$. This fact is now established from several studies on both the SDSS (Omand et al. 2014; Lin et al. 2019) and CANDELS/3D-HST (Whitaker et al. 2017; F18; Lin et al. 2019). Galaxies with the same starting $M_*$ but different starting radii therefore grow in $M_*$ at the same rate, and each starting mass cohort stays together. Because all relationships are power laws, it follows that $\Delta \log R_{e,s}$, $\Delta \log \Sigma_{1,s}$, and $\Delta \log M_{BH,s}$ are constant through time. Hence, galaxies at a given $M_*$ that start out with large radii always have smaller-than-average BHs, and vice versa. If quenching requires reaching some threshold BH mass that does not vary too rapidly with $M_*$, large galaxies will tend to do this later and attain a larger stellar mass before quenching. The fact that the starting radius confers an $M_{BH}$ enhancement to some galaxies and a penalty to others is what creates the sloped quenching boundaries in $R_e-M_*$, $\Sigma_1-M_*$, and $M_{BH}-M_*$ in Figure 2.

The final step is to specify the equation of the quenching boundary. For this, our first-stage, empirical power-law model assumes that the total effective energy emitted by a BH during its lifetime is equal to some total energy quota that must be transferred to the halo gas in order to inhibit further infall onto the galaxy. We call the latter the halo quenching-energy quota, and in the empirical model, we parameterize it as a simple power of halo mass. Thus,

$$ E_{BH} = E_{\text{quench}} \quad (9) $$

or

$$ \eta M_{BH} c^2 = \text{const.} \times M_{\text{vir}}^t, \quad (10) $$

where $\eta$ is an efficiency factor for converting BH accretion energy into useful feedback energy. We reiterate that crossing over the boundary corresponds to entry into the green valley, and thus that BHs are just starting to affect their halos and galaxies at that point. The use of $E_{\text{quench}}$ to denote the critical halo energy is meant to be generic and does not mean that the energy transferred must be thermal; it could be mechanical, cosmic rays, or any other mode of energy/momentum transfer that can stop halo gas from falling onto the galaxy. This expression for necessary BH feedback as a power law in $M_{\text{vir}}$ is needed only during the star-forming phase in order to calculate the location of the quenching boundary and is not needed (or assumed to be valid) after galaxies pass over the quenching boundary. The value of $t$ and the (evolving) constant are adjusted to match the adopted quenching boundary, which was taken to be offset by 0.2 dex in $\Sigma_1$ from the quenched ridgeline above. We will find in the next section that the empirical value of $t$ must be close to 2.0 in order to match the observed slopes of the quenched ridgelines in $\Sigma_1-M_*$ and $M_{BH}-M_*$, and the constant must vary with time in order to match the moving boundaries in those same spaces.

## 3.2. Quenching Boundary Power Laws

We now show that the above assumptions suffice to compute the slopes of the quenching boundaries in our various spaces as a function of $s$, $t$, $v$, $p$, and $q$. Galaxies on the boundary at the same epoch obey

$$ M_{BH} \sim M_{\text{vir}}^t $$

by Equation (10). But

$$ M_{BH} \sim \Sigma_1^v $$

from Equation (5), and

$$ M_{\text{vir}} \sim M_*^{1/s} $$

from Equation (1). Substituting, we have

$$ \Sigma_1 \sim M_*^{t/sv}, \quad (11) $$

which gives the slope, $t/sv$, of the instantaneous boundary in $\Sigma_1-M_*$. Now apply the mapping

$$ \Sigma_1 \sim M_*^q R_e^{-p} $$

from Equation (3) to obtain

$$ R_e \sim M_*^{q/p-t/svp}, \quad (12) $$

which gives the slope, $q/p-t/svp$, of the instantaneous boundary in $R_e-M_*$. Finally, $M_{BH} \sim \Sigma_1^v$ gives

$$ M_{BH} \sim M_*^{t/s}, \quad (13) $$

which is the slope, $t/s$, of the instantaneous boundary in $M_{BH}-M_*$.

Of these five exponents, $p$ and $q$ are set by the $p$, $q$ mapping while $s$, $t$, and $v$ are set by three additional observational constraints:

1. Slope of the quenching boundary in $\Sigma_1$ versus $M_*$. Even though the spaces $\Sigma_1-M_*$, $R_e-M_*$, and $M_{BH}-M_*$ are all closely related (for star-forming galaxies), the quenched ridgeline is narrowest and best determined in $\Sigma_1-M_*$, and so we start there. We have further assumed that the boundary is offset by a constant distance below the ridgeline by 0.2 dex at all masses (based on a feature in the local SDSS $\Sigma_1-M_*$ distribution found by Luo et al. 2020) and therefore has the same slope as the ridgeline. B17 found the slope of this ridgeline to be nearly constant at 0.66 from $z=3$ to now. Equation (11) now gives$^{40}$

$$ \frac{\Delta \log \Sigma_1}{\Delta \log M_*} = \frac{t}{vs} = 0.66 \pm 0.02, \quad (14) $$

where the error is based on the scatter in four independent measurements at four different redshifts in B17. This is the first constraint.$^{41}$ As discussed below, we do not use the quenched ridgeline slope in $R_e-M_*$ as an

<footnote>$^{40}$ B17 set the cut for quenched galaxies at -0.7 dex below the star-forming ridgeline, which differs from our definition using -1.0 dex in Figure 1. However, the ridgeline slope is insensitive to this definition because green-valley galaxies and fully quenched galaxies occupy nearly the same regions in $\Sigma_1-M_*$ (F13, B17).</footnote>

<footnote>$^{41}$ For simplicity, the shape of the quenched ridgeline in $\Sigma_1-M_*$ is assumed here to be a pure power law, as found by B17. Recent work by Luo et al. (2020) for SDSS galaxies indicates some downward curvature in $\Sigma_1-M_*$ at high masses, which may actually improve quenching predictions for high-mass galaxies, as discussed under Variant 1 in Section 6.1.</footnote>
observational constraint because of the sensitivity of effective radii to dry merging and other post-quenching processes. However, the slope of the quenching boundary in $R_e-M_*$ will be needed in the future. Equation (12) now gives

$$ \frac{\Delta \log R_e}{\Delta \log M_*} = \frac{q}{p} - \frac{t}{vps} \\ = 1.10 - \frac{t}{vs} = 0.44 \pm 0.02, \qquad (15) $$

where the second line follows from inserting $t/vs = 0.66 \pm 0.02$ and $p = 1.0 \pm 0.01$, $q = 1.1 \pm 0.01$ from Appendix A. Note that although the quenching boundary and quenched ridgeline are assumed to be parallel in $\Sigma_1-M_*$, they are not parallel when mapped into $R_e-M_*$ because quenched galaxies have different Sérsic indices and $R_e$ due to dry merging. This is a major reason why the quenched ridgeline in $R_e-M_*$ has the steep slope of 0.75 (vdW14) whereas the quenching boundary slope in Equation (15) is only 0.44.

2. *Slope of the quenching boundary in $M_{BH}$ versus $M_*$:* As shown in Figure 2(c) here and more clearly in Figure 3 of Terrazas et al. (2017), the separation between fully quenched and star-forming galaxies in $M_{BH}-M_*$ is a fixed ~1.5 dex today for galaxies above $10^{10} M_\odot$, and we assume that this is true at all redshifts. This implies that the quenching boundary is parallel to the quenched ridgeline in $M_{BH}-M_*$, as it is in $\Sigma_1-M_*$, and the quenched slope in $M_{BH}-M_*$ can therefore be used as another observational constraint on the unknown exponents $s$, $t$, and $v$. The local observed $M_{BH}-M_*$ relation consists of a narrow ridgeline of fully quenched galaxies with star-forming galaxies scattering below it (KH13; Terrazas et al. 2016; van den Bosch 2016; Krajnović et al. 2018). We equate the observed quenched ridgeline slope to that of local ellipticals and bulges, which is 1.16 from KH13. Thus,

$$ \frac{\Delta \log M_{BH}}{\Delta \log M_*} = \frac{t}{s} = 1.16 \pm 0.08, \qquad (16) $$

where the error comes from KH13. Because classical bulges and ellipticals follow the same law, this slope is insensitive to the precise definition of the quenched sample. This is the second constraint.

3. Finally, the slope of the SMHM relation, $s$, is observed to be $1.75 \pm 0.17$ based on the abundance-matching results of RP17 and Behroozi et al. (2019). This is the third constraint.

Constraints are now available from Items 1, 2, and 3, and the solution is uniquely determined. The combination $s = 1.75 \pm 0.17$, $t/s = 1.16 \pm 0.08$, and $t/vs = 0.66 \pm 0.02$ implies $t = 1.16 \times 1.75 = 2.03 \pm 0.24$ and $v = 1.16/0.66 = 1.76 \pm 0.13$, and from this last one, we derive $u = 2.0v = 3.52 \pm 0.26$, as was employed in Equation (5) above.

We pause here to consider the uncertainties in these values. The formal errors are small, but the important unknowns are whether the cartoon model in Figure 4 is actually correct and whether the same $M_{BH}-\Sigma_1$ scaling relations are in fact valid at all redshifts. The general constancy of the cartoon relations is supported by the fixed power-law appearance of the $\Sigma_1-M_*$ and $M_{BH}-M_*$ distributions back in time (lines 12 and 16 in Table 1)

and by correct predictions for both evolving zero points when combined with the moving boundary hypothesis, as discussed in Section 7.1. In short, the constant slopes and offsets of quenched and star-forming galaxies in $\Sigma_1-M_*$ and $M_{BH}-M_*$ over time impose strict constraints on the constancy of the quenching process, which is why these relations are so important. If further data confirm these diagrams, something like the cartoon model must be taking place.

Another concern is that the slope of the quenched ridgeline has been used to set the slope of the star-forming galaxies in Figure 4, but the quenched slope may have evolved after quenching due to mergers. For example, to match the scaling laws of quenched galaxies today requires $R_e \sim M_*^2$ along the merging trajectory (e.g., Naab et al. 2009; vdW09) and thus $\sigma_e \sim M_*^{-0.5}$, i.e., $\sigma_e$ actually declines as galaxies merge! Both of these vectors (in $R_e$ versus $M_*$ and $\sigma_e$ versus $M_*$) point at angles to their respective ridgelines, and evolution along these vectors might therefore change the original slopes. Mergers of clusters of galaxies add further orbital energy to central galaxies, puff them up, and reduce $\sigma_e$ even more (Boylan-Kolchin et al. 2006; Lauer et al. 2007). Some evidence of steepening in $M_{BH}$ versus $\sigma_e$ at high masses has been noted (Lauer et al. 2007; Gültekin et al. 2009; Krajnović et al. 2018; KH13). To summarize, effective quantities like $R_e$ and $\sigma_e$ are strongly impacted by dry merging, and their expected evolutionary vectors are at angles to today's ridgelines, with unknown consequences for the ridgeline slopes over time. This makes effective quantities unsuitable as slope constraints, which is why we have not used them.

By contrast, merging is probably not as important for $\Sigma_1-M_*$, because centers are more immune to merging than envelopes and the slope of $\Sigma_1-M_*$ is actually measured to be constant over time (B17). In $M_{BH}-M_*$, if premerger BHs merge along with their host galaxies and if BH mass is conserved, galaxies would tend to move along this ridgeline as its slope is near unity. This slope would then also be invariant. Hence, it is plausible that both our calibrating slopes are stable over time.

A final matter is the low value of $u$ that we have derived, which implies that $M_{BH} \sim \sigma_1^{3.52 \pm 0.26}$. This value stems from $v = 1.76$, because $u$ is $2.0 \times v$. But $v$ is simply the ratio of the slopes of the quenched ridgelines in $\Sigma_1-M_*$ and $M_{BH}-M_*$, both of which are well determined. Alternative values for $u$ are $4.40 \pm 0.30$ for quenched ridgeline galaxies from KH13 and $5.35 \pm 0.23$ for all galaxies together from van den Bosch (2016). At first glance, these seem incompatible with our value of $3.52$ and also with each other. (These laws use $\sigma_e$ while we use $\sigma_1$, but the difference is probably not large.) However, as described in Section 7.7 below, the Figure 4 cartoon model can explain the differences as well as the low value that we get for $v$ (and $u$) for star-forming galaxies alone.

To summarize, the empirical power-law model predicts a series of slopes and relations that can be fitted to observed quantities. The following three exponents are the basic unknowns: $M_{BH} \sim M_{vir}^t$ (the halo energy quota needed to halt gas cooling), $M_* \sim M_{vir}^s$ (the SMHM relation), and $M_{BH} \sim \Sigma_1^\nu$ (the assumed relation between BH mass and central stellar density for star-forming galaxies). These three exponents, $t$, $s$, and $v$, predict the slopes of the two quenching boundaries in $\Sigma_1-M_*$ and $M_{BH}-M_*$, which are observed, and $s$ is separately measured from abundance matching. This makes three constraints, which together determine the three exponents.
<table><caption>**Table 2**
Model Ingredients and Assumptions</caption><thead><tr><th>No.</th><th>Item</th><th>References</th></tr></thead><tbody><tr><td>1</td><td>Halo masses grow according to mean trajectories vs. time derived from N-body simulations.</td><td>Rodríguez-Puebla et al. (2017)</td></tr><tr><td>2</td><td>Central galaxies evolve along a power-law SMHM relation while star-forming, with constant slope s = 1.75, constant zero point, and no scatter.</td><td>Rodríguez-Puebla et al. (2017)</td></tr><tr><td>3</td><td>The median effective radii &lang;R<sub>e</sub>&rang; of star-forming galaxies are 0.02 &times; R<sub>vir</sub> at all masses and redshifts.</td><td>Huang et al. (2017)<br>Somerville et al. (2018)</td></tr><tr><td>4</td><td>There is significant real scatter in both R<sub>e</sub> and &Sigma;<sub>1</sub> among star-forming galaxies at fixed M<sub>*</sub> and R<sub>vir</sub>.</td><td>Omand et al. (2014)<br>van der Wel et al. (2014)<br>Fang et al. (2013)<br>Barro et al. (2017)</td></tr><tr><td>5</td><td>A galaxy's offset in radius &Delta; log R<sub>e</sub> from the median &lang;R<sub>e</sub>&rang; is constant during the star-forming phase. Star-forming galaxies evolve along parallel tracks in R<sub>e</sub>-M<sub>*</sub>.</td><td>Assumed, Section 4.1</td></tr><tr><td>6</td><td>For star-forming galaxies, &Sigma;<sub>1</sub> = const. M<sub>*</sub><sup>q</sup> &times; R<sub>e</sub><sup>-p</sup>, where q = 1.1 and p = 1.0.</td><td>Fitted to CANDELS data in Appendix A.<br>Omand et al. (2014)</td></tr><tr><td>7</td><td>The star formation rate of galaxies on the star-forming main sequence is independent of R<sub>e</sub> and &Sigma;<sub>1</sub> at fixed M<sub>*</sub>.</td><td>Whitaker et al. (2017)<br>Fang et al. (2018)<br>Lin et al. (2019)<br>Luo et al. (2020)</td></tr><tr><td>8</td><td>&Delta; log &Sigma;<sub>1</sub> = 0.2 dex and &Delta; log M<sub>BH</sub> = 1.5 dex across the green valley, independent of stellar mass and redshift.</td><td>Terrazas et al. (2017)<br>Assumed in Section 3.1<br>and discussed in Section 7.5</td></tr><tr><td>9</td><td>BH masses in star-forming galaxies grow as M<sub>BH</sub> = const. &Sigma;<sub>1</sub><sup>1.76</sup> with the same zero point and slope at all redshifts and no scatter.<sup>a</sup></td><td>Assumed, Section 3.1</td></tr><tr><td>10</td><td>The effective energy radiated from the BH is E<sub>BH</sub> = &eta;M<sub>BH</sub>c<sup>2</sup>, where &eta; = 0.01.</td><td>Assumed, Section 4.3</td></tr><tr><td>11</td><td>The energy quota that must be absorbed by the halo gas to initiate entry into the green valley in the empirical power-law model is E<sub>quench</sub> = k<sub>e</sub>(z) (M<sub>vir</sub>/10<sup>12</sup>M<sub>&odot;</sub>)<sup>t</sup> with no scatter.</td><td>Assumed, Section 4.3</td></tr><tr><td>12</td><td>The energy quota that must be absorbed by the halo gas to initiate entry into the green valley in the binding-energy model is 4 &times; E<sub>bind</sub> = 4 &times; 1/2 f<sub>bar</sub> M<sub>vir</sub> V<sub>vir</sub><sup>2</sup> with no scatter.</td><td>Assumed, Section 5</td></tr><tr><td>13</td><td>Mergers during the star-forming and green-valley phases are ignored, and galaxies do not rejuvenate in the green valley.</td><td>Assumed throughout</td></tr></tbody></table>

**Note.**

# 4. The Evolving Empirical Power-law Model

## 4.1. Model for Galaxy Radii: Evolutionary Tracks

The power-law exponents *s*, *t*, and *v* are now specified. The present section adds moving zero points and evolves star-forming galaxies along evolutionary tracks to create an empirical evolving power-law model, so-called because the halo energy quota, *E*<sub>quench</sub>, is represented as a power law in *M*<sub>vir</sub> with the exponent *t* set to 2.0, and a moving boundary is introduced to match the moving (red) quenching boundaries in Figure 1. In what follows, this name is shortened to the “empirical power-law model.”

The model starts by using our knowledge of how halos evolve in $R_{vir}$ versus $M_{vir}$ to derive evolutionary tracks for central star-forming galaxies in $R_e-M_*$, $\Sigma_1-M_*$, and $M_{BH}-M_*$. The ingredients of and assumptions in the model are listed in Table 2. Accurate evolutionary tracks are necessary during the star-forming phase in order to have galaxies cross the quenching boundaries in the right place. To achieve this, a halo-based model is used to predict the evolution of $M_*$ and $R_e$ for star-forming galaxies. The starting point is the formulae for halo-mass growth and related galaxy parameters of central star-forming galaxies from Rodríguez-Puebla et al. (2015) and RP17. Figure 6 plots the resulting stellar and virial parameters for three representative different-mass halos and

their associated galaxies versus time. The adopted constant value for the ratio $\langle R_e \rangle / R_{vir} = 0.02$ comes from the abundance-matching results in RP17 for star-forming galaxies below $M_* = 10^{10.5} M_\odot$ (see Item 2 in Section 3.1 and Figure D1 in Appendix D). See also the discussion of radius evolution in Somerville et al. (2018), which justifies a constant value.

The conventions we are using for $M_{vir}$ and $R_{vir}$ are based on the definition of halo overdensity $\Delta_{vir}$ by Bryan & Norman (1998) for an isothermal sphere mass profile. The formula relating virial mass and radius is

$$M_{vir} = \frac{4\pi}{3} \Delta_{vir} \rho_{m} R_{vir}^{3}, \qquad (17)$$

where $\rho_m$ is the mean matter density $\Omega_m\rho_c$ (Rodríguez-Puebla et al. 2016a) and $\Delta_{vir}$ is given by Bryan & Norman (1998). The quantity $V_{vir}$ is defined such that $M_{vir} \equiv V_{vir}^2 R_{vir}/G$.

The tracks in Figure 6 are used to evolve central star-forming galaxies in the $R_e-M_*$ plane through time. In order to check them and to develop handy power-law approximations, Figure 7 replots them versus $M_{vir}$ and $M_*$, where the relations in Equations (3) and (5) are used to derive $\langle \Sigma_1 \rangle$ and $\langle M_{BH} \rangle$ from $\langle R_e \rangle$ and $M_*$. The solid lines are for the standard empirical power-law model, which has a fixed zero point in the SMHM relation. The dashed lines refer to the Variant 2 model, which puts the evolving zero point into the SMHM relation instead of

<footnote>a This is the default assumption in both the empirical power-law and binding-energy models. An alternative version in which BHs are smaller than this in star-forming galaxies and cross over this scaling law as they enter the green valley is discussed in Section 7.8.</footnote>
<figure><img src="image_6.png" /><figcaption>**Figure 6.** Evolutionary tracks for $M_*$, $M_{vir}$, $\langle R_e \rangle$, $R_{vir}$, $T_{vir}$, and $V_{vir}$ for three sample star-forming central galaxies vs. redshift. $M_{vir}$ and $R_{vir}$ are based on the evolutionary model for dark halos summarized in RP17. $M_*$ is derived from $M_{vir}$ using the SMHM relation $M_* \sim M_{vir}^s$ with $s = 1.75$ from Section 3.1. The constant value $\langle R_e \rangle / R_{vir} = 0.02$ is taken from RP17 (see Item 2 in Section 3.1 and Figure D1). Tracks start at the starting redshift, $z_s$, which is the redshift where stellar mass $M_* = 10^9 M_\odot$. The quantity $V_{vir}$ is used in Section 5, which discusses the halo binding-energy interpretation of the halo energy quota. The quantity $T_{vir}$ is used to compute heating and cooling rates in Variants 3 and 4 in Appendix D.3.</figcaption></figure>

in the halo quenching condition as in Equation (19) and is discussed in Section 6 and Appendix D.

The predicted tracks for $\langle R_e \rangle$ in Figure 7 show the motions of single galaxies through time, whereas plots of observed $R_e$ versus $M_*$ are snapshots of different galaxies at the same time. Such snapshots are seen to populate a broad range in $R_e$, have a shallower slope than the individual evolutionary tracks, and have a ridgeline zero point that evolves upwards in time (vdW14). By definition, the peak of the ridgeline is the locus of galaxies with $R_e = \langle R_e \rangle$. The properties of a ridgeline can therefore be calculated from the formula for $\langle R_e \rangle$ versus time from RP17. From Equations (1), (2), and (17), we get

$$ \begin{aligned} \langle R_e \rangle &\sim R_{vir} \\ &\sim M_{vir}^{1/3} \Delta_{vir}^{-1/3} \rho_m^{-1/3} \\ &\sim M_*^{1/(3s)} \Delta_{vir}^{-1/3} \rho_m^{-1/3}. \end{aligned} \qquad (18) $$

Thus, the predicted ridgeline slope $\Delta \log \langle R_e \rangle / \Delta \log M_* = 1/(3 \times 1.75) = 0.19$, and the zero point evolves as $\Delta_{vir}^{-1/3} \rho_m^{-1/3}$, which is approximately equal to $h(z)^{-0.74}$ from $z = 0.5$ to $z = 3$, where $h(z) \equiv H(z)/H(0)$ and the Bryan & Norman (1998) prescription for $\Delta_{vir}$ is used. Both trends agree

well with observed data as noted by vdW14, who found an average ridgeline slope of 0.22 and a zero-point evolution of $h(z)^{-0.66}$. Further support for this radius model comes from vD15, who reviewed observational and theoretical data on evolutionary tracks for star-forming galaxies in $R_e-M_*$. They found an average evolutionary track slope $\Delta \log R_e / \Delta \log M_* = 0.3$, in reasonable agreement with our average track slope of 0.40 in Figure 7.

A final consistency check can be made by comparing to slopes in $\Sigma_1-M_*$ using the $p, q$ mapping. A given log slope $m$ in $R_e-M_*$ corresponds to a slope $q - pm$ in $\Sigma_1-M_*$. For $p = 1$, $q = 1.1$, the predicted ridgeline slope for star-forming galaxies is then $1.1 - 1.0 \times 0.19 = 0.91$. This compares well to the fitted ridgeline slope in B17, which averages 0.88 at all redshifts. The corresponding evolutionary track slope is $1.1 - 1.0 \times 0.40 = 0.70$, just slightly steeper than the quenched ridgeline slope of 0.66. An additional check on the evolutionary track slope can be made by identifying a median star-forming galaxy at $10^{8.9} M_\odot$ at $z = 2.5$ in B17, computing a stellar mass increase to $10^{10.2} M_\odot$ at $z = 0.75$ (RP17), and placing it in the middle of the star-forming cloud in $\Sigma_1-M_*$ of B17 at the new redshift and mass. The resulting changes are
<figure><img src="image_7.png" /><figcaption>**Figure 7.** The same evolutionary tracks from Figure 6 for central star-forming galaxies are replotted vs. mass, with $\langle\Sigma_1\rangle$ and $\langle M_{BH}\rangle$ added. The black solid lines are for the empirical power-law model with fixed SMHM zero point. The black dashed lines are analogous tracks for the Variant 2 model, which has an evolving zero point in the SMHM relation (Section 6.2). The corresponding red lines show approximate power-law slopes (with arbitrary zero points) that are valid for galaxies in the range $M_{vir} = 10^{11.5-12.5} M_{\odot}$ and $z = 3$ to 0.5. The first exponent refers to the empirical power-law model with a fixed SMHM zero point; the second exponent refers to the Variant 2 model with an evolving SMHM zero point. These slopes are useful for rough estimates; the real curves are used when computing model predictions.</figcaption></figure>

1.30 dex in $M_*$ and 0.90 dex in $\Sigma_1$, for a net evolutionary track slope of 0.69, in excellent agreement with the predicted slope (inferred from $R_e-M_*$) of 0.70.

Finally, star-forming galaxies are assumed to evolve along the power-law SMHM relation (Equation (1)) until they enter the green valley.

In summary, a simple halo-based model for $M_*$ and $R_e$ appears to match fairly well the observed average radii of star-forming galaxies, their evolution in time after $z=3$, and estimates of how individual galaxies move in $R_e-M_*$ and $\Sigma_1-M_*$.

## 4.2. Starting Conditions for Each Galaxy

In order to ensure that galaxies cross the quenching boundary in the right place, star-forming galaxies must be properly positioned in $\Sigma_1-M_*$ and $R_e-M_*$ to start with. Power-law evolutionary tracks are now in place, but each model galaxy needs to be provided with the appropriate starting conditions for entry onto its track. We start by choosing a galaxy stellar mass today and use the observed HMSM relation at $z=0$ from RP17 to calculate a corresponding halo mass.42 This halo mass is then propagated backwards in time using

the RP17 halo-mass evolution formula together with the power-law SMHM relation, $M_* = k_m M_{vir}^s$, until $M_* = 10^9 M_\odot$, which sets the starting redshift, $z_s$. This halo mass, stellar mass, and virial radius are then evolved forward using the halo evolution model of RP17 and the same SMHM relation. For the empirical power-law model, which we are describing here, the zero point of the SMHM relation is kept fixed, and $k_m$ is a constant.

It remains to scatter the starting effective radii. A mean starting $\langle R_{e,s} \rangle$ at $z_s$ for galaxies with $M_* = 10^9 M_\odot$ is computed from the $R_{vir}$ formula in Equation (2), and galaxies are scattered above and below this value with a Gaussian distribution having $\sigma \log(R_e/\langle R_e \rangle) = 0.20$ dex, which is the observed rms scatter in the radii of star-forming galaxies (vdW14; Somerville et al. 2018). An important simplification, noted above, is that the star formation rates of ridgeline galaxies on the star-forming main sequence are observed to be independent of $R_e$ at all redshifts (Omand et al. 2014; Whitaker et al. 2017; F18; Lin et al. 2019) so that all galaxies in a given $M_*$ cohort stay together in mass regardless of their starting $R_e$.

## 4.3. Evolutionary Tracks in Relation to the Quenching Boundary

When the above evolutionary tracks were computed, a surprising fact emerged. We had anticipated that galaxies would evolve through the $R_e-M_*$ and $\Sigma_1-M_*$ planes toward the sloping boundary, where SSFRs fall -0.45 dex below the

<footnote n="42">42 The HMSM relation is the inverse of the SMHM relation, being the mean value of $M_{vir}$ observed at fixed $M_*$ today. It is taken from RP17 and is not a power law.</footnote>
star-forming main-sequence value. Galaxies would move toward this boundary because their evolutionary tracks were angled relative to it, and eventually the two would intersect. However, we saw above that the observed quenching boundary in $\Sigma_{1}-M_{*}$ has slope 0.66 (B17) whereas the estimated evolutionary track has slope 0.70, i.e., nearly the same. Likewise, the quenching boundary in $R_{e}-M_{*}$ has slope 0.44 whereas the evolutionary track has slope 0.40, again nearly the same. In other words, the evolutionary tracks of star-forming galaxies make very acute “angles of attack” relative to the quenching boundaries. The concern then arises that many galaxies would never quench, i.e., that galaxies initially located far from the boundary would never cross it.

The solution to this problem turned out to be the moving boundary. A downward motion of -0.3 dex in the $\Sigma_{1}-M_{*}$ quenched ridgeline zero point was measured by B17 from z = 2.5 to 0.5, and the corresponding upward motion in $R_{e}$ is +0.3 dex. Our quenching boundaries are tied to these ridgelines. This motion, which seems small at first sight, is actually the key ingredient in the empirical power-law model that causes galaxies to quench. Individual galaxies evolve roughly parallel to the boundary, but meanwhile, the boundary sweeps over them, from above in $\Sigma_{1}$ and from below in $R_{e}$.

To make further progress, this motion needs to be inserted into the expression for the quenching boundary. It was shown above in Equation (10) that equating the total energy emitted by the BH to a power of halo mass matches the slope of the quenching boundary if t = 2.0. We now take that power-law relation and add a moving zero point:

$$E_{\text{BH}} = E_{\text{quench}} \qquad (19)$$

or

$$\eta \frac{M_{\text{BH}}}{10^8 M_{\odot}} c^2 = k_{\text{e}}(z) \left( \frac{M_{\text{vir}}}{10^{12} M_{\odot}} \right)^t, \qquad (20)$$

where $k_{\text{e}}(z)$ varies with time and we have taken $\eta = 0.01$ (see Oppenheimer et al. 2020). The left side is the effective emitted BH energy, $E_{\text{BH}}$, and the right side is the halo energy quota for quenching, $E_{\text{quench}}$. From $M_{\text{BH}} \sim \Sigma_{1}^{\nu}$ (Equation (7)), we find that the change $\Delta \log k_{\text{e}}(z)$ with time is equal to $\nu$ times the change in the zero point of the quiescent $\Sigma_{1}$ ridgeline in B17 at fixed $M_{*}$. A fit of this zero point versus redshift, called $B^{Q}$ in B17, is shown in Figure B1 in Appendix B. It yields $\log B^{Q} = 9.44 + 0.66 \log(H(z)/H(0))$, and thus, $k_{\text{e}}(z) = k_{\text{e}}(0) \times h(z)^{0.66\nu}$.

The constant $k_{\text{e}}(0)$ at z = 0 is used as a free parameter to match the absolute location of the quenching boundary at z = 0, assuming $\eta = 0.01$. Rather than match to the quenched ridgeline itself, however, we have offset entry into the green valley by 0.2 dex below the quenched ridgeline in $\Sigma_{1}-M_{*}$ based on the data on SDSS green valley galaxies of Luo et al. (2020), and we apply that same shift at all redshifts. The value of $k_{\text{e}}(0)$ is then found to be $10^{17.66} \text{ km}^2 \text{ s}^{-2}$. Finally, the moving boundary in $\Sigma_{1}-M_{*}$ is mapped onto $R_{e}-M_{*}$ using the p, q mapping. This prescription for $k_{\text{e}}(z)$ (and t) perfectly matches the red lines in Figure 1.

We pause here to make two points. The first is that the moving zero point in Equation (20) is needed to match the empirically determined quenching boundaries (red lines in Figure 1), which move with time. However, a more physical model for quenching based on real halo properties is coming up in Section 5, which also naturally produces a moving

<figure><img src="image_8.png" /><figcaption>Figure 8. Sample evolutionary paths for a selection of galaxies of a given stellar mass in $R_{e}$ vs. $M_{*}$ in the empirical power-law model. The black line is for an average galaxy with radius $R_{e} = \langle R_{e} \rangle$ and $M_{*} = 10^{10.6} M_{\odot}$ today (if it does not quench first). The blue tracks above and below that line show galaxies of the same starting mass but scattered up and down in $R_{e}$ by $\pm 0.75 \sigma_{\text{Re}}$ and $\pm 1.5 \sigma_{\text{Re}}$ (two extra values are added on the low side). Tracks with lower values of $R_{e}$ have higher BH masses. The observed distribution of galaxies is shown for comparison at z = 0.5–1.0, colored by $\Sigma_{1}$. The diagonal red lines show the location of the quenching boundary, which moves upward with time. Entry to the green valley occurs when a galaxy crosses the instantaneous boundary; such locations are shown by the red/white dots, which correspond to z = 1.76, 1.36, 0.90, and 0.43, from the bottom. The evolution of galaxies in $R_{e}-M_{*}$ after entering the green valley is not covered by the model and is not shown.</figcaption></figure>

boundary. The moving boundary thus has two independent justifications: the need to match observational data on the one hand and the changing physical properties of gaseous halos on the other. The second point is that the presentation and labeling of quantities in Equation (19) explicitly evoke the notion of a BH–halo energy balance, with BH energy on the left side and a halo quenching-energy quota on the right side. We think that this presentation is helpful because it points the way to the more physically based model in Section 5. On the other hand, one could omit the labels and view Equation (20) simply as an empirical relation between $M_{\text{BH}}$ and a power of $M_{\text{vir}}$ (plus a moving zero point) and save the physical interpretation until later. That is equally valid and, depending on one’s preferences, might be regarded as a “cleaner” presentation.

Evolutionary paths and quenching points of typical Milky Way–mass galaxies in $R_{e}$ versus $M_{*}$ are illustrated in Figure 8. These tracks use the actual trajectories of halos, not the approximate power laws. The black line is for an average galaxy with radius $R_{e} = \langle R_{e} \rangle$ and $M_{*} = 10^{10.6} M_{\odot}$ today. The blue tracks above and below this line are galaxies of the same starting mass but scattered up and down in $R_{e}$ by $\pm 0.75 \sigma_{\text{Re}}$ and $\pm 1.5 \sigma_{\text{Re}}$, where $\sigma_{\text{Re}} = 0.2$ dex (two extra values are added on the low side). Tracks with lower values of $R_{e}$ have higher BH masses and vice versa. The observed distribution of galaxies at z = 0.5–1.0 is shown for comparison in the background. The diagonal red lines show the location of the quenching boundary at various times. Entry into the green valley occurs when a galaxy crosses the instantaneous boundary (red and white dots). It is seen that the tracks are inherently quite parallel to the boundary, and it is really the upward motion of the boundary zero point that causes galaxies to cross it.

To summarize, a galaxy is identified as quenching perceptibly (i.e., entering the green valley) when the quenching condition in Equation (20) is satisfied, and the redshift when
<figure><img src="image_9.png" /></figure>

**Figure 9.** The BH energy output ($E_{BH}$, rainbow bands) and the halo quenching-energy quota ($E_{quench}$, red dashed lines) are shown for star-forming galaxies vs. redshift in the evolving power-law model. Three representative $z=0$ stellar mass/halo-mass combinations are shown. $E_{BH}$ is the total effective energy, $\eta M_{BH}c^2$, emitted by the BH to that epoch, and $E_{quench}$ is the total energy quota that needs to be absorbed by the halo as galaxies enter the green valley (Equation (20)). Galaxies below the red dashed lines have large radii, small BHs, and are still star-forming, while galaxies above the red dashed lines have small radii, large BHs, and have entered the green valley. (Note: because galaxies above the lines are in the green valley, the real values of their $M_{BH}$ and $E_{BH}$ would be different than those shown. Star-forming galaxies do not actually populate these locations; see text.) The width of the rainbow bands corresponds to the scattering of the initial values of $R_{e,s}$ by $\pm 1.5\sigma_{Re} = 0.30$ dex, which scatters $M_{BH}$ by $\pm 0.52$ dex from Equation (8). Galaxies with smaller $R_e$'s (and bigger BHs) at a given $M_*$ quench sooner. Notice the acute “angle of attack” between $E_{BH}$ and $E_{quench}$ and the tendency for a larger fraction of massive galaxies to quench. These trends fall out naturally from setting $t=2.0$ and the evolving zero point of $E_{quench}$ to match the observed quenching boundary. The blue solid lines show the halo quenching-energy quota, $4 \times E_{bind}$, from the binding-energy model in Section 5. It is seen that $E_{quench}$ follows $4 \times E_{bind}$ very closely, showing that the empirical power-law energy quota we have used to quench a halo (Equation (20)) is closely equal to the actual gas-binding energies of halos. The red dashed lines and blue solid lines here correspond to the red and blue quenching boundaries in Figure 1.

this happens and its location are marked on the tracks. If the model is working, the locus of perceptibly quenching galaxies at a given redshift should match the adopted quenching boundary in $\Sigma_1-M_*$ at that redshift, which is taken to be 0.2 dex below the observed quenched ridgeline.

All of the power-law expressions used in the model and their slopes and zero points are collected and explained for reference in Appendix C.

## 4.4. Results for the Empirical Power-law Model

Results for the empirical power-law model are shown in Figures 9 and 10. Figure 9 plots the BH heating-energy, $E_{BH}$, and the halo heating-energy quota, $E_{quench}$, versus redshift for three different halo–galaxy combinations. The curves are labeled with ending stellar and halo masses today. The red dashed lines show $E_{quench}$ from Equation (20), which depends only on halo mass and is therefore shown as a single line for each $M_*$. The rainbow-colored bands show $E_{BH}$; their spread reflects the range of BH masses coming from a $\pm 1.5\sigma$ scatter in $R_e$ (smaller $R_e$ and larger BHs are at the upper edges). A given galaxy quenches when its $E_{BH}$ line crosses above $E_{quench}$. We stress that $E_{quench}$ is not necessarily thermal energy but rather expresses some sort of generic energy (or momentum) transfer from the BH to the halo gas that stops cooling.

We note that the $E_{BH}$ trajectories have been plotted ignoring the fact that galaxies perceptibly quench when they cross over the $E_{quench}$ lines and their subsequent $M_*$ and $M_{BH}$ evolution switches to the green-valley vectors in Figure 4. Star-forming galaxies therefore do not actually exist in these locations, but we have left them in place in order to better convey the fraction of galaxies at each mass and redshift that have quenched.

Several features are apparent. First, the $E_{quench}$ curves bend over with time for a given stellar mass despite the fact that the halo mass is always growing. In the empirical power-law model, this bend comes from the decreasing value of $k_e(z)$ in

the $E_{quench}$ zero point (Equation (20)). Without this bend, $E_{quench}$ would continue to rise with $M_{vir}$, and galaxies would quench very late or not at all. Note that small-$R_e$ galaxies with bigger BHs quench sooner (upper envelopes). Finally, one sees that a critical parameter is the “angle of attack” that a given $E_{BH}$ trajectory makes with $E_{quench}$. We mentioned this angle before in the $\Sigma_1-M_*$ and $R_e-M_*$ diagrams; the concept here is similar —if the angle of attack is large, all galaxies of a given starting $M_*$ will quench in a short space of time, and thus at nearly the same mass. This means that the length of the quenching locus (and the width of the green valley) for that mass galaxy will be small. On the other hand, if the angle is acute, quenching takes a long time from smallest to largest $R_e$, and the extent of the quenching locus in $M_*$ will be large. Matching the observed large width of the green valley requires acute angles of attack, as discussed further below.<fnref n="43" />

Figure 10 shows corresponding snapshots of galaxies in the $R_e-M_*$ and $\Sigma_1-M_*$ planes at four redshifts. Each dot is a model star-forming galaxy of a given mass, and attached to each point is a short line showing its recent evolutionary track. Five model galaxies are plotted in each mass cohort, corresponding to galaxies scattered in radius by $\pm 0.75\sigma_{Re}$ and $\pm 1.5\sigma_{Re}$. They are superimposed on a shaded region denoting the observed $\pm 1.5\sigma_{Re}$ range of observed star-forming galaxies at that redshift from vdW14. The model galaxy points populate the observed range of galaxies fairly well, which validates our scheme for starting galaxies at $10^9 M_\odot$ with $R_e = 0.02 \times R_{vir}$ and growing $R_e$ in proportion to $R_{vir}$ thereafter.

Galaxies that are entering the green valley at each redshift are shown by the pale green points, which are populated more densely in starting mass in order to show their locus clearly.

<footnote n="43">43 The need for acute angles of attack to match the wide width of the green valley is exaggerated here by the fact that all of the assumed power-law relations in the model are scatter free. Scatter around the scaling relations would help to broaden the green valley naturally and should be included in future versions of the model.</footnote>
<figure><img src="image_10.png" /><figcaption>**Figure 10.** Model star-forming galaxies from the empirical power-law model in $R_e$ vs. $M_*$ (panel (a)) and in $\Sigma_1$ vs. $M_*$ (panel (b)). Seven masses and four redshift ranges are shown. Each dot is a galaxy, and the short line attached to it is its recent evolutionary track. The pale green circles show the loci of galaxies that are entering the quenching boundary at that redshift, and the red lines underneath them are the observed quenching boundaries from Figure 1 (red lines). A perfect match has been attained by choosing the right values of s, t, and v and by allowing the zero point in Equation (20) to evolve appropriately. The model galaxy points have been scattered by a total of $\pm 1.5\sigma_{Re}$ in $R_{e,s}$, and the shaded regions outline the corresponding $\pm 1.5\sigma_{Re}$ areas that are covered by star-forming galaxies in vdW14. The similarity of the two distributions indicates that our method for choosing starting values of $R_{e,s}$ (Section 4.2) has achieved a reasonable match between the regions populated by the model star-forming galaxies vs. observed galaxies. Notice that the individual evolutionary tracks are almost parallel to the quenching boundary and that galaxies cross the boundary because the boundary moves down in $\Sigma_1$ (or up in $R_e$) to intercept them. Notice also the “tongue” of star-forming galaxies at late times that extends to higher masses than observed. This is discussed in Variant 1 in Section 6.</figcaption></figure>

The straight red lines underneath them are the quenched ridgelines of B17 displaced downward by 0.2 dex in $\Sigma_1$ to represent the adopted quenching boundary (entry into the green valley). They are the same as the red lines shown in Figure 1. The good match between the red lines and the green points verifies that the zero point and redshift evolution of the quenching expression in Equation (20) have been set correctly. The predicted slopes also match, indicating that the adopted values of t/vs, p, and q are also proper. They are taken from Section 3 and are t = 2.03, v = 1.76, s = 1.75, p = 1.0, and q = 1.1.

**5. Interpretation of $E_{quench}$ in Terms of Halo Gas-binding Energy**

Galaxies enter the green valley in the empirical power-law model when the total effective energy radiated by their accreting BHs equals the halo-heating-energy quota $E_{quench} = k_e(z) (M_{vir}/10^{12} M_\odot)^t$ in Equation (20).$^{44}$ The formulation of $E_{quench}$ as a power law of $M_{vir}$ was adopted to match the approximate

power-law nature of the quenching boundary. However, there was little physical motivation, and it does not quite work—an ad hoc evolving zero-point correction had to be added.

An additional obvious energy to be considered is the binding energy of the halo gas. It is already known that $E_{BH}$ is of the same order of magnitude as halo gas-binding energy and thus that actual BHs emit enough power to influence their halos (Silk & Rees 1998; Ostriker et al. 2005; Bower et al. 2017; Oppenheimer 2018; Davies et al. 2019; Nelson et al. 2019; Oppenheimer et al. 2020). If in addition we can show that halo gas-binding energy varies similarly to $E_{quench}$ as a function of $M_{vir}$ and time, we will have additional evidence that BH energy beating against halo gas-binding energy is the gatekeeper for quenching.

Our adopted expression for halo gas-binding energy is

$$E_{\text{bind}} = \frac{1}{2} f_{\text{bary}} M_{\text{vir}} V_{\text{vir}}^2, \qquad (21)$$

where $f_{bary}$ is the fraction of $M_{vir}$ in the form of hot gas, set equal to 0.1 throughout. Inserting this into Equation (20), we now have a new quenching condition:

$$E_{\text{BH}} = f_{\text{B}} \times E_{\text{bind}} \qquad (22)$$

or

$$\eta \frac{M_{\text{BH}}}{10^8 M_{\odot}} c^2 = f_{\text{B}} \times \frac{1}{2} f_{\text{bary}} M_{\text{vir}} V_{\text{vir}}^2, \qquad (23)$$

where $f_B$ is a multiple of halo gas-binding energy that needs to be supplied by the BH at the quenching boundary.

The quantity $f_B \times E_{bind}$ is plotted as solid blue lines for the three typical halo masses in Figure 9, and the equivalent

<footnote>$^{44}$ In adopting this criterion, we have ignored both radiative losses in the AGN feedback and possible energy transfer to the halo by supernovae (e.g., Dekel & Silk 1986) and by gravitational heating (e.g., Dekel & Birnboim 2008; Johansson et al. 2009). Recent hydro simulations give some sense of the importance of the first two factors. In both EAGLE (Bower et al. 2017) and TNG50 (Nelson et al. 2019), successful AGN winds tend to conserve energy, which justifies the use of cumulative BH energy in Equations (22) and (23). In contrast, supernova winds are prone to large radiative losses (e.g., Somerville & Davé 2015), and adding them to the current framework would require a major energy-loss correction factor. In TNG50, the relative roles of supernovae versus BHs are hard to separate near $M_*=10^{10.5} M_\odot$, but AGN winds clearly dominate above this (Nelson et al. 2019).</footnote>
quenching boundaries for this new model are the blue lines in Figure 1. Comparing these with the observed quenching boundaries (red lines), we see fairly good agreement in slope, average zero point, and zero-point evolution. The agreement in average zero point has been optimized by choosing $f_B = 4$, which is the best value given $\eta = 0.01$.

In contrast to the average zero point, the agreements in slope and zero-point evolution are unforced, as can be seen by considering the ratio of $E_{\text{quench}}$ to $E_{\text{bind}}$. At fixed $z$, $E_{\text{quench}}$ varies as $M_{\text{vir}}^2$, while $E_{\text{bind}}$ varies as $M_{\text{vir}}V_{\text{vir}}^2 \sim M_{\text{vir}}^{1.67}$, which is similar but slightly shallower. This slight difference accounts for the shallower slope of the blue lines in Figure 1 and also explains why the blue line lies above the red dashed line for small galaxies but below it for massive galaxies in Figure 9.<fnref n="45" />

For a given galaxy over time, it may also be shown using the above equations plus $M_{\text{vir}} \sim V_{\text{vir}}^3/h(z)^{1.11}$ that $E_{\text{bind}}/E_{\text{quench}}$ evolves as $M_{\text{vir}}^{-0.33}h(z)^{0.74-0.66v} \sim M_{\text{vir}}^{-0.33}h(z)^{-0.42}$. But $M_{\text{vir}}$ grows approximately as $h(z)^{-1.50}$ over this redshift range, and so $E_{\text{bind}}/E_{\text{quench}} \sim h(z)^{0.08}$. In other words, $M_{\text{vir}}$ is rising while $h(z)$ is falling, and the two almost cancel along a given evolutionary track. This is why the red dashed lines and the blue solid lines are parallel in Figure 9. An additional insight comes from noting that $E_{\text{quench}}$ varies as $k_e(z)M_{\text{vir}}^{2.0}$ while $E_{\text{bind}}$ varies as $M_{\text{vir}}^{2.0}R_{\text{vir}}^{-1}$. Thus, the evolving term $k_e(z)$ in the empirical power-law model is needed to mimic the evolution of $R_{\text{vir}}^{-1}$ over time. The question of the zero-point evolution of the boundary over time is different from the evolution of an individual galaxy and is considered in Section 7.1.

This is a good place to make a point about the uniformity of quenching. If $f_B$ varies from halo to halo, say, with rms scatter $\delta_{\text{fB}}$ dex at a given epoch, this will induce the same amount of scatter in $M_{\text{BH}}$ and a scatter of $\delta_{\text{fB}}/v$ in $\Sigma_1$. The $\Sigma_1$ ridgeline has rms scatter equal to 0.16 dex today (F13), which implies that $\delta_{\text{fB}} = 0.28$ dex. However, this is an upper limit because age stratification has been detected (Tacchella et al. 2017; Luo et al. 2020), so the real scatter in $\delta_{\text{fB}}$ must be less. The point is that the halo binding-energy criterion in Equation (23), represented by $f_B$, must be highly uniform from halo to halo in order to create the tight ridgelines in $\Sigma_1-M_*$ that are seen.

To summarize, in agreement with previous work, we have shown that the energy available in the course of BH growth is the same order of magnitude as or greater than the binding energy of the halo gas, and thus that BHs potentially produce enough energy to modify the state of their halo gas and retard gas infall. In addition, we have shown that halo gas-binding energy has the right general behavior versus $M_{\text{vir}}$ and time to match the observed sloping quenching boundaries in plots of galaxy structural parameters and thus to explain important structural scaling laws of quenching galaxies. The notion that halo gas-binding energy controls BH growth, and through that, quenching, was modeled by Booth & Schaye (2010, 2011) and has been seen and analyzed in the EAGLE simulations (Davies et al. 2019; Oppenheimer et al. 2020) and in IllustrisTNG (Terrazas et al. 2020). Comparisons to these and other models are discussed in Section 7.

# 6. Variants of the Empirical Power-law Model

This section describes four variations on the empirical power-law model. Summaries are presented here, and details are given in Appendix D. Each variant fails in some respect, but the reasons are instructive, yielding important insights into the sensitivities and basic physics of the problem. Readers interested in final conclusions can skip directly to the Discussion in Section 7.

## 6.1. Variant 1: A Broken Power Law for the SMHM Relation

A feature of the empirical power-law model that does not match observed data is the long “tongue” of star-forming galaxies extending up to $M_* > 10^{12} M_\odot$ at late redshifts in Figure 10. This tongue is due to the fact that we have used the curved HMSM relation from RP17 at $z=0$ to find the starting redshifts, $z_s$, but have calculated stellar masses going forward using the power-law SMHM relation of Equation (1), $M_* = k_m M_{\text{vir}}^3$. Because the latter keeps climbing with $M_{\text{vir}}$, massive galaxies that start near the knee of the real SMHM relation today can wind up with a higher stellar mass than their assumed initial mass. Galaxies with $M_* > 10^{10.7} M_\odot$ and large $R_e$ (and small BHs) are the objects that populate the massive tongue.

To attempt to cure this problem, the Variant 1 broken-power-law model uses the same power-law SMHM relation for star-forming galaxies up to $10^{10.7} M_\odot$ but switches to a flatter power law for SMHM above that. Results are shown in Figures D1 and D2 in Appendix D.1. With the new relation, massive galaxies now largely stop growing but are permanently stuck in a star-forming state because their BHs remain small. An excess of massive star-forming galaxies should therefore accumulate over time, which is not seen.

However, there may be another solution, based on recent information about the shape of the quenching boundary from a large sample of SDSS galaxies (Luo et al. 2020). The quenched ridgeline in $\Sigma_1-M_*$ in this new sample is seen to bend over to lower $\Sigma_1$ at high stellar mass, which would create a matching bend in the quenching boundary. Our assumed linear power law lies about 0.2 dex above the new SDSS ridgeline at $10^{11.4} M_\odot$, and it is seen from Figure 10 that bending the ridgeline down by this amount would remove the problematic massive galaxies.

The lesson of Variant 1 is that the fraction of quenched versus star-forming galaxies at high mass is very sensitive to the exact shape of the quenching boundary, and a deep physical understanding of the boundary will be needed to predict the fraction accurately. A corollary is that the shape and amplitude of the quenched mass function of quenched galaxies at high values of $M_*$ is difficult to predict.

## 6.2. Variant 2: Move the Evolving Zero Point from $E_{\text{quench}}$ to the SMHM Relation

We have noted that the quenching boundary evolves in $\Sigma_1-M_*$, and the empirical power-law model matches this by evolving the zero point, $k_e(z)$, of $E_{\text{quench}}$ in Equation (20). An alternative approach is to move the evolving zero point to the SMHM relation in Equation (1) and remove it from $E_{\text{quench}}$. The expression that we have used for the SMHM relation in this case is derived in Appendix D.2. It makes the zero point of $M_*$ versus $M_{\text{vir}}$ larger by about 0.5 dex at $z=0$ compared to $z=2$. This grows BHs more quickly, which allows them to

<footnote n="45">⁴⁵ We note that the slope of the blue lines could be increased in Figure 1 to match that of the red lines if the slope of the SMHM relation, $s$, were reduced from 1.75 to 1.44. This change is larger than the quoted error of 0.17 but is perhaps not outside the realm of possibility.</footnote>
overcome halos without resorting to $k_e(z)$. Changes in the SMHM relation back in time that bear on this are described in Appendix D.2. There is some evidence for a change in the right direction, but the amount is smaller than required.

This model achieves the same match to the quenching boundaries as the empirical power-law model, but the growth of galaxy properties versus galaxy mass is different, as shown by the dashed lines in Figure 7. BHs grow more rapidly at late times, so quenching starts too late, i.e., not until around $z = 1$ (see Figure D3 in Appendix D.2). Trying to move the onset of quenching earlier by increasing $\eta$ in Equation (20) produces overquenching at late times. Also, the more rapid growth in $E_{\text{quench}}$ with time no longer parallels the slower growth in halo gas-energy content, $E_{\text{bind}}$, as shown by the differing red and blue lines in Figure D3 in Appendix D.2, and so the appealing physical argument equating $E_{\text{quench}}$ with halo gas-binding energy is therefore also broken.

The important lesson from Variant 2 is that changes in the zero point of the real SMHM relation back in time would require compensating changes in $E_{\text{quench}}$, which could threaten its interpretation as halo gas-binding energy and also disturb the agreement on time of quenching. Accurate measurements of the SMHM relation back in time are therefore needed to finalize the model.

### 6.3. Variants 3 and 4: Mass Rates or Energy Totals?

A major difference between the empirical power-law model and other halo-heating models is that the latter often use heating and cooling rates as opposed to cumulative energy totals as in Equation (20). Use of rates is especially common in SAMs (e.g., Bower et al. 2006; Croton et al. 2006; Somerville et al. 2008; Henriques et al. 2015). Variants 3 and 4 explore whether Equation (20) can be recast in terms of rates. In both variants, the halo energy quota, $E_{\text{quench}}$, on the right side of Equation (20) is replaced by the halo mass-cooling rate, $\dot{m}_{\text{cool}}$, from Croton et al. (2006; also White & Frenk 1991; Cole et al. 2000; Somerville et al. 2008; Guo et al. 2011; Henriques et al. 2015).

To replace the total BH energy on the left side of Equation (20) by a BH heating rate, the BH mass-accretion rate must also be known. In Variant 3, this rate is obtained by differentiating the instantaneous $M_{\text{BH}}$ values from the empirical power-law model. In Variant 4, the BH accretion expression introduced by Croton et al. (2006) for radio mode is used, but as this has no dependence on $R_e$, this model is only a one-parameter family. Equations are given in Appendix D.3, and the accretion rates versus time are compared in Figure D4. It is seen that the Variant 3 BH accretion rates fall strongly at late times, the Variant 4 accretion rates grow super-exponentially at late times, while total BH masses in the empirical power-law model grow moderately at all times.

The consequences of these different rates are illustrated in Figures D5–D8. Variant 3 fails because the late turndown in BH accretion rates means that galaxies stop quenching at late times. Variant 4 fails because the late turnup in BH accretion rates makes the characteristic green-valley mass fall by 1.0 dex from $z = 1.4$ to now, and because it is only a one-parameter family.

Collectively, Variants 3 and 4 illustrate two points: (1) BH and halo effects need to march closely together in time in order to feed galaxies smoothly into the green valley at constant rate and mass, and this is achieved more easily with cumulative

energies than with rates; and (2) a successful quenching model needs two parameters working together to match the observed 2D structural diagrams of galaxies, which automatically excludes one-parameter models that depend only on mass, like Variant 4.

## 7. Discussion

The Discussion has nine subsections. The first reviews and explains the successes that were ascribed to the model in Table 1. From there, we consider possible explanations for the scatter in the effective radii of star-forming galaxies, which is a necessary ingredient in the model. Subsequent sections discuss quenching, compaction, BH growth through the green valley, and the different BH scaling laws found by different authors. The last opens the way to thinking about the BH-halo connection, what is meant by BH-galaxy “coevolution,” and the implications of having three separate yet matching quenching boundaries in $\Sigma_1-M_*$, $\Sigma_1-M_{\text{BH}}$, and $M_{\text{BH}}-M_*$. We close with a comparison of our picture to the L16 quenching model.

### 7.1. Scorecard for the Binding-energy Model

This section reviews the ability of the model to match the observational data listed in Table 1. By “model” we henceforth mean the empirical power-law model from Section 4 that includes the cartoon evolutionary tracks in Figure 4 but with the halo energy quota, $E_{\text{quench}}$, replaced by $4 \times$ halo gas-binding energy, $E_{\text{bind}}$—i.e., the blue lines in Figure 1. This model will henceforth be called the “binding-energy model”, and it is the major result of this paper. The basic ingredients of this model are (1) the known theoretical evolution of the halo properties $M_{\text{vir}}$, $R_{\text{vir}}$, and $V_{\text{vir}}$; (2) the assumed links between halo properties and galaxy properties $M_*$, $\langle R_e \rangle$, and $\langle \Sigma_1 \rangle$; (3) the scattering of the ratio $R_e/\langle R_e \rangle$; (4) the criterion for entry into the green valley $E_{\text{BH}} = f_B \times E_{\text{bind}}$; and (5) the relation $M_{\text{BH}} = \text{const.} \times \Sigma_1^{1.76}$ for star-forming galaxies.

Column 1 in Table 1 indicates the successes of this model by using the line number without brackets to indicate observations that we think are well explained by the model and the line number with brackets to indicate partial success. We note that the binding-energy model has only one free parameter, and that is the constant $f_B$ in Equation (23), which has been set equal to 4 to make the zero point of the blue lines match that of the red lines in Figure 1 on average.<fnref n="46" /> Aside from the general insistence on power laws, none of the observations in Table 1 were used to determine the five basic ingredients above, and hence all of the successes in Table 1 are new.

As noted, the starting point for the binding-energy model places star-forming central galaxies at the centers of dark halos with $M_*$ and $\langle R_e \rangle$ derived from the halo $M_{\text{vir}}$ and $R_{\text{vir}}$. As pointed out by vdW14 and noted in Section 4.1, this halo-based approach matches rather well the ridgeline slope, zero-point evolution, and evolutionary track slopes of star-forming galaxies in $R_e-M_*$ and, through the $p, q$ mapping, also in $\Sigma_1-M_*$. These are lines 4, 5, 8, and 9 in Table 1, which are indicated as full successes. Moving down, the slopes of the quenching boundaries are constant with time in the binding-energy model, as observed, but the predicted slope is a bit too

<footnote n="46">⁴⁶ The adopted value of $f_B$ depends on the assumed value of the BH efficiency factor $\eta$, which we have set to 0.01, and we are actually adjusting the ratio of these two quantities.</footnote>
shallow in $\Sigma_1-M_*$ and a bit too steep in $R_e-M_*$, and the zero-point motions with time are a bit too small (i.e., the blue lines do not quite match the red lines in Figure 1). Because of these mismatches, lines 10, 11, and 12 for the quenching boundaries—and their associated quenched ridgelines—have been recorded as partial successes.

We pause here to comment on line 12, which says that the $M_{BH}-M_*$ relation for quenched galaxies is offset to higher values at early times (Peng et al. 2006; Treu et al. 2007; Woo et al. 2008; Gu et al. 2009; Bennert et al. 2010, 2011; Decarli et al. 2010; Ding et al. 2017, 2020). This offset has often been interpreted to mean that BHs form before the stars of their host galaxies (e.g., Peng et al. 2006; Treu et al. 2007; Bennert et al. 2011; Ding et al. 2017), making BHs appear larger back in time for their host stellar mass. This interpretation is contrary to the binding-energy model, which says that BHs experience their major growth in the green valley, and thus that BHs form after most of the stars. The model has a different explanation for the offset, namely, that the quenching boundary in $M_{BH}-M_*$—and therefore its associated ridgeline—was higher in the past. The predicted motion can be calculated from the observed motion of the $\Sigma_1-M_*$ ridgeline and is approximately $v \times 0.3 = 0.5$ dex back to $z=3$, which is in good agreement with the observed shift in $M_{BH}-M_*$ (Bennert et al. 2011; Ding et al. 2017, 2020). The point is that offsets of $M_{BH}$ versus $M_*$ with time are not necessarily due to growth-time differences between stars and BHs but may instead reflect motions of the quenching boundary.47 The same point was made by Caplar et al. (2018), who explored a model similar to ours in which BHs grow rapidly just as galaxies enter the green valley.

Moving on, line 16 in Table 1 says that the observed zero point of the quenched ridgeline in $\Sigma_1-M_*$ evolves downward as $h(z)^{0.66}$ (Appendix B), while the observed zero point of the star-forming ridgeline in the same plane evolves downward by the same amount (B17). This is the key equality that maintains star-forming galaxies in the same relationship to quenched galaxies in $\Sigma_1-M_*$ and is the reason why the diagram looks the same at all redshifts (except for the bulk vertical shift). This synchronized motion means that the quenching boundary cuts through the star-forming distribution over the same range of masses at all redshifts. This is seen in the predicted model quenching loci of Figure 10, where the smallest galaxies always enter quenching near $10^{9.5} M_\odot$ and are fully quenched by $10^{10} M_\odot$, which is the observed lower edge of the high-mass quenching channel at all redshifts (B17). Finally, this same constant vertical offset between the quenching boundary and the star-forming cloud ensures a smooth, steady flow of galaxies across the boundary at all times (line 3 in Table 1). The delicacy of achieving this is illustrated by Variant 3, which fails to do this.

The synchronized downward motion of the quenched and star-forming clouds in $\Sigma_1-M_*$ is thus a fundamental feature of the data that demands an explanation, and the model provides estimates for the motions of both zero points separately. The evolution of the star-forming ridgeline at fixed mass in $\Sigma_1-M_*$ goes as $R_e^{-p} \sim R_e^{-1} \sim R_{vir}^{-1} \sim h(z)^{0.74}$ (this last is from Equation (18) in Section 4.1). The zero-point motion of the quenching boundary (and thus the quenched ridgeline) varies as $E_{bind}^{1/v} \sim (M_{vir}^2/R_{vir})^{1/v}$, which at fixed $M_*$ goes as

$R_{vir}^{-1/v} \sim h(z)^{0.74/v} \sim h(z)^{0.42}$. The two exponents of $h(z)$ are therefore similar but not quite the same, differing by the factor $v=1.76$. This explains why the zero point of the blue line moves down more slowly than the red line in Figure 1, leaving a residual in $\Sigma_1$ of about 0.1 dex. Based on this partial success, line 16 is given a bracket.

It is helpful to rephrase the above mathematical language into something more intuitive. Simply put, BH masses (at fixed $M_*$) grow smaller at late times because galaxies have larger $R_e$. However, halo binding energies (at fixed $M_{vir}$) also become smaller with time because halos have larger $R_{vir}$. Together, these two factors cause the downward boundary motion. This point is not new; it was made by Booth & Schaye (2011), who predicted that the quenching zero point of models in which the BH energy output is balanced against halo gas-binding energy would decrease with time. Our new insight is that the parallel evolution of BHs and halos is why the green-valley mass remains constant with time. In more detail, the decline in halo binding energy is not quite fast enough to keep halos and BHs in perfect lockstep, but real halos might have additional factors, such as late-infall dry-up or the disappearance of cold streams, that might help them quench more easily at late times (see Section 7.3).

A moving quenching boundary has another consequence. If the evolution of individual galaxies after quenching does not precisely follow the moving boundary, early-quenched galaxies will be systematically offset from late-quenched ones, leading to age stratification in the quenched ridgeline due to “progenitor bias” (e.g., Carollo et al. 2013; Poggianti et al. 2013).

B17 listed three factors that help quenched galaxies stay close to the downward-evolving ridgeline in $\Sigma_1-M_*$. They are (1) stellar mass-loss coupled with adiabatic expansion, (2) core scouring by BH mergers in merging galaxies, and (3) dry merging. With three mechanisms helping, we expect that age stratification in $\Sigma_1-M_*$ will be only moderate.

In contrast, post-quenching evolution in $M_{BH}-M_*$ ought not to track the moving boundary as closely. Stellar mass loss moves galaxies to the left, above the relation. Also, if BHs merge along with their host galaxies, this would tend to evolve objects along the original locus because the log slope of the relation is near unity, at best leaving their vertical residuals unchanged. Finally, core scouring does not alter either $M_{BH}$ or $M_*$. The net result is that quenched objects are not expected to follow the moving ridgeline as closely in $M_{BH}-M_*$ as in $\Sigma_1-M_*$, and thus progenitor bias and age stratification should be larger in this ridgeline.

Age stratification in the quenched ridgeline has been measured in $M_{BH}-M_*$ by Martín-Navarro et al. (2018) and in $\Sigma_1-M_*$ by Tacchella et al. (2017) and Y. Luo et al. (2020, in preparation). Early-quenched galaxies fail to fully keep up with the moving boundary in both spaces, but the mismatch in $\Sigma_1-M_*$ is significantly smaller, as the above logic predicts. The magnitudes of both effects are generally consistent with the model, and line 15 has been marked as a success.

Line 2 notes the wide extent of the green valley in $M_*$. As described in Section 4.4, this is explained by the highly acute “angles of attack” of the evolutionary tracks relative to the quenching boundary in $\Sigma_1$ versus $M_*$, $R_e$ versus $M_*$, and $E_{BH}$ versus $E_{quench}$ (and also $E_{bind}$), which enable star-forming galaxies with the same $M_*$ but only slightly different values of $R_e$ and $\Sigma_1$ to evolve to very different final quenching masses.

<footnote n="47">47 The fact of a moving zero point of $M_{BH}$ versus $M_*$ does not contradict our fundamental assumption that the zero point of $M_{BH} = \text{const} \times \Sigma_1^{-1.76}$ remains the same because $M_*$ and $\Sigma_1$ are not the same quantity.</footnote>
These acute angles are a natural feature of the model, and line 2 is accordingly marked as a success.

Line 14 in Table 1 notes that the ridgelines of quenched galaxies are narrow in $\Sigma_1-M_*$ and $M_{BH}-M_*$ but broader in $R_e-M_*$, especially at high z. The basic narrowness in $\Sigma_1$ and $M_{BH}$ is due to the fact that galaxies stop evolving and pile up along these loci. The broader width of $R_e-M_*$ reflects the greater impacts on $R_e$ of compaction (see Section 7.4) and dry minor merging (e.g., Daddi et al. 2005; Trujillo et al. 2006; van Dokkum et al. 2008; Naab et al. 2009; van Dokkum et al. 2010; Oser et al. 2012; Szomoru et al. 2012). The evolution of the latter goes as $R_e \sim M_*^2$, which is steeper than the quenched ridgeline itself ($R_e \sim M_*^{0.75}$), thus broadening (and steepening) it significantly. By invoking a little help from compaction and dry merging, we consider line 14 a success.

Finally, line 13 in Table 1 says that $M_{BH}$ varies roughly as some power $M_*^w R_e^{-w}$ (or as $\sigma_e^{2w}$) for quenched and star-forming galaxies together but as a smaller power than w for quenched galaxies alone (compare slopes from van den Bosch 2016 and KH13). This difference turns out to be the same problem as the low value of $u = 3.52$ in Section 3.1. Both issues are discussed in Section 7.7, where we show that steep BH growth in the green valley can naturally explain both. The same steep BH growth also explains the greater separation between quenched and star-forming galaxies in $M_{BH}-M_*$ than in $\Sigma_1-M_*$, and both lines 13 and 18 are marked as successes.

The above discussion covers all the items without brackets in Table 1. Lines 6 and 7 (with brackets) involve the slope and zero point of the quenched ridgeline in $R_e-M_*$. As these depend on how $R_e$ evolves after quenching, which is not included in the model, only partial credit is given for getting the locus in roughly the right place. The different zero-point evolution of quenched and star-forming galaxies in $R_e-M_*$ versus $\Sigma_1-M_*$ (line 17) is attributable to the larger impact of post-quenching dry merging on $R_e$ than on $\Sigma_1$, and so partial credit is given to that line as well.

This completes the review of all the items in Table 1. Yet another test would use X-ray luminosities to test whether star-forming galaxies with higher $\Sigma_1$ (and lower $R_e$) at fixed $M_*$ have larger BHs. The findings of Ni et al. (2019) suggest this, but their star-forming sample also includes green-valley galaxies and thus is not quite pure enough. A reanalysis of a pure main-sequence sample is underway (D. D. Kocevski et al. 2020, in preparation). X-ray data for green-valley galaxies are discussed in Section 7.5.

## 7.2. The Origin of Scatter in Effective Radii

A basic ingredient of the binding-energy model is that star-forming galaxies at a fixed stellar mass have different radii. Radius differences correlate with central $\Sigma_1$ (Figure 1), and these induce differences in BH mass and BH feedback. Where do differences in radii come from? Because the model assumes that all star-forming galaxies of the same stellar mass have the same dark halo mass, $R_e$ differences must come from some parameter other than $M_{vir}$.

The connection between dark halos and galaxy radii was reviewed by Somerville et al. (2018) and Wechsler & Tinker (2018). Mo et al. (1998) introduced the formula for the exponential disk scale length $r_d = 1/\sqrt{2}\lambda R_{vir}$, where $\lambda$ is the dimensionless halo spin parameter. This suggested that $\lambda$ variations might cause the scatter in $R_e$. However, the observed rms scatter in $R_e$ is only 0.2 dex (vdW14) compared to 0.25 dex

for $\lambda$ in the models, which leaves no room for any other sources of scatter or measurement errors (Desmond & Wechsler 2015). Hydro simulations also show that cooling gas gains or loses angular momentum unpredictably while falling onto galaxies (Pichon et al. 2011; Stewart et al. 2013; Danovich et al. 2015; Stevens et al. 2017; Jiang et al. 2019), blurring any connection between the original and final angular momenta. The net result is that galaxy radius is found to be almost uncorrelated with halo $\lambda$ in EAGLE (Desmond et al. 2017), and Jiang et al. (2019) find at most poor correlations for VELA (Ceverino et al. 2014) and for NIHAO (Wang et al. 2015).

An interesting alternative to $\lambda$ is halo concentration, c. Jiang et al. (2019) find that the radii of VELA and NIHAO galaxies obey $R_e = 0.02 R_{vir} (c/10)^{-0.7}$, with galaxies being smaller in higher-concentration halos. This is plausible because a greater fraction of the matter in a high-concentration halo falls in at earlier times, when the universe was denser. The binding-energy model then says that such galaxies should have relatively larger BHs owing to their higher stellar densities. For the purposes of this paper, the origin of radius scatter needs no explanation—the $\Delta \log R_e$ offset of a galaxy may simply be a random variable that is assigned at birth and is unconnected to halo properties—it is only required that $\Delta \log R_e$ stay constant during the star-forming phase so that galaxies evolve on parallel tracks in $R_e-M_*$. However, a picture in which the two leading structural parameters of star-forming galaxies, $M_*$ and $R_e$, arise directly from two leading parameters of dark halos, $M_{vir}$ and c, would be both simple and elegant.

Recent analyses of EAGLE and IllustrisTNG see the predicted effect in which more concentrated halos form stars faster, have larger BHs, larger values of $E_{BH}/E_{bind}$, and are more quenched today (Davies et al. 2019, 2020; Oppenheimer et al. 2020). The effect in EAGLE is large—halos in the top quartile of $E_{BH}/E_{bind}$ were quenched several billion years ago while halos in the bottom quartile are still making stars today. Furthermore, Davies et al. (2019) showed that it is the integral of BH feedback, not the instantaneous BH feedback rate, that ultimately ejects baryons, and the best barometer for the degree of quenching is the ratio $E_{BH}/E_{bind}$.

These properties support our assumption in the binding-energy model that quenching is determined by the ratio $E_{BH}/E_{bind}$, not by instantaneous heating or cooling rates. EAGLE quenching becomes perceptible when $E_{BH}/E_{bind} \sim 2$ and is complete when $E_{BH}/E_{bind} \sim 20$ whereas we predict the range of 4–20 for this ratio (see Section 7.5 below). This is a promising level of agreement given the simplicity of our model.<fnref n="48" />

Guided by EAGLE and IllustrisTNG, it is tempting to create a revised binding-energy model by adding halo concentration to $M_{vir}$ as a secondary parameter. This addition of a second halo parameter in addition to $M_{halo}$ would be an example of assembly bias (Wechsler & Tinker 2018). However, doing so would raise important questions. For one thing, N-body simulations show that mass-accretion rates onto halos at fixed time and mass vary with halo concentration (Lin et al. 2019),

<footnote n="48">⁴⁸ It is important to note that the evacuation of gas from halos during quenching that is predicted by the binding-energy model cannot persist when halos later grow to cluster size. Large clusters are seen to be full of gas at approximately the cosmic baryon-to-dark-matter ratio (e.g., Tumlinson et al. 2017), and so they must refill. EAGLE and IllustrisTNG clusters in fact do refill (Davies et al. 2020; Oppenheimer et al. 2020). A different form of BH feedback must take over to maintain quenching at the centers of gas-rich groups and clusters. This later evolutionary stage is not covered by the model.</footnote>
which, if unopposed, might disturb the observed absence of any trend in star formation rate versus galaxy size at fixed $M_*$ (Omand et al. 2014; Whitaker et al. 2017; F18; Lin et al. 2019). Indeed, star formation rates do vary with size in EAGLE (Oppenheimer et al. 2020), in disagreement with the observations. A dodge around this could posit weaker galactic winds in small galaxies on account of their higher escape velocities (Dutton et al. 2010; Lin et al. 2019). But even if current star formation can be evened out in this way, high-concentration halos still need to form stars earlier in order to give their galaxies smaller radii. This means that their star formation rates had to be different in the past, i.e., that mass cohorts cannot stay together over long times, which violates another of our key assumptions. Finally, it is not clear that halos in a given mass cohort would maintain their relative concentration rankings at different epochs, which would be needed to keep their evolutionary tracks parallel.

In summary, adding halo concentration (and perhaps also time of formation) as a second parameter to the binding-energy model is an intriguing suggestion, but important issues remain unresolved.

### 7.3. Quenching in the Binding-energy Model Compared to Other Models

It is instructive to compare the BH accretion and quenching processes in the binding-energy model to analogous mechanisms in current SAM and hydro models. The following summarizes the kinds of feedback mechanisms in these models.

**Stellar feedback.** This is the most basic feedback mechanism, present in all models. Energy/momentum from young massive stars and supernovae drives a wind that expels gas from a galaxy during the star-forming phase.

**Quasar mode** (Croton et al. 2006; Somerville et al. 2008; Henriques et al. 2015), also called “bright mode” or “high mode” (Somerville et al. 2008; Pillepich et al. 2018) or “radiative mode” (Somerville & Davé 2015). This is the principal mode of BH growth in many models, triggered by mergers (Hopkins et al. 2006), disk instabilities (Cattaneo et al. 2006), or other unspecified secular processes (Martín-Navarro et al. 2018). Key features of quasar mode (in most models) are a high BH accretion rate, formation of a thin accretion disk, and high radiative efficiency from the accretion disk (see review by King & Pounds 2015). It is generally quasar-mode winds that quench star formation in SAMs and in many hydro models.

**Radio mode** (Croton et al. 2006; Somerville et al. 2008; Henriques et al. 2015), also called “kinetic mode” or “low mode” (Weinberger et al. 2017; Pillepich et al. 2018), or “jet mode” (Somerville & Davé 2015). Hot gas is accreted slowly onto a central BH, producing high-velocity winds and/or radio jets that escape into the halo. Key features of radio mode (in most models) are a low BH accretion rate, formation of a thick accretion disk, and low radiative efficiency (see review by Yuan & Narayan 2014). Radio mode may quench, as in IllustrisTNG (Weinberger et al. 2017; Nelson et al. 2019; Davies et al. 2020; Terrazas et al. 2020), or it may maintain a low star formation rate in gas-rich halos after other mechanisms have quenched.

**Halo quenching** (Rees & Ostriker 1977; Blumenthal et al. 1984; Silk et al. 1986; Birnboim & Dekel 2003; Keres et al. 2005; Dekel & Birnboim 2006; Dekel et al. 2019). As halos grow in mass, their gas becomes hotter and more dilute, which retards cooling. A critical halo mass has been proposed,

$M_{crit} \sim 10^{11.5-12} M_{\odot}$ (Dekel & Birnboim 2006; Dekel et al. 2019), above which halos fill with hot gas that cannot easily cool, and star formation shuts down. This transition is also accompanied by a change from “cold streams” to “hot flows” (Keres et al. 2005; Dekel & Birnboim 2006; Keres et al. 2009), which are less able to penetrate to the central galaxy.

**Late-infall dry-up.** The specific accretion rate of gas onto halos declines after $z \sim 1$ in all cosmologies, regardless of $\Omega_m$ or $\Omega_{\Lambda}$ (Dekel & Mandelker 2014). Though not a feedback process per se, it lowers the amount of gas available for star formation (Bouché et al. 2010; Davé et al. 2011, 2012; Krumholz & Dekel 2012; Dekel et al. 2013; Lilly et al. 2013; Dekel & Mandelker 2014; Forbes et al. 2014; Feldmann 2015; Mitra et al. 2015; Rodríguez-Puebla et al. 2016b) and retards halo cooling by keeping gas densities low at late times. Along with the transition from cold streams to hot streams (Keres et al. 2005; Dekel & Birnboim 2006), this makes late halos easier to quench and would help to lower the zero point in the quenching boundary relation (red lines) in $\Sigma_I - M_*$, as seen in Figure 1.

Three features of the binding-energy model collectively make it unique compared to all other models. The first is that the power of AGN feedback to quench is based on the rate of BH gas accretion versus the rate of halo gas cooling in some models (e.g., Bower et al. 2006; Croton et al. 2006; Somerville et al. 2008; Henriques et al. 2015) whereas the binding-energy model integrates the total BH energy over time and compares it to the total binding energy of the halo gas. The impact of the difference was explored in Variants 3 and 4 in Section 6, where it was shown that rate formulations do not feed gas to BHs at the right rate. As noted above, two recent hydro simulations also seem to exhibit the total-energy criterion (EAGLE: Davies et al. 2019, Oppenheimer et al. 2020; and IllustrisTNG: Davies et al. 2020, Terrazas et al. 2020).

The second difference is that multiple mechanisms typically operate together in SAM and hydro models to initiate quenching, and the decisive mechanism in any one galaxy can be hard to identify. Such is the case, for example, in the Somerville et al. (2008) SAM (Pandya et al. 2017) and in TNG50 (Nelson et al. 2019). In contrast, entry into the green valley in the binding-energy model is the well-defined point when total AGN output surpasses $f_B \times$ halo gas-binding energy, where $f_B \approx 4$.

Finally, the third major difference between the binding-energy model and most other models is the mechanism that maintains the constant value of the mean green-valley mass, $\langle M_{*,GV} \rangle$, over time. A plausible value for this mass, cited in many works, is the value of $M_*$ at the critical halo cooling mass $M_{crit}$ (e.g., Birnboim & Dekel 2003; Keres et al. 2005; Dekel & Birnboim 2006; Dekel et al. 2019). A second characteristic mass has also emerged in recent hydro models from balancing supernova wind energy against galaxy or halo potential well depth (Dubois et al. 2015; Anglés-Alcázar et al. 2017; Bower et al. 2017; Habouzit et al. 2017; Davies et al. 2019; Dekel et al. 2019; Henriques et al. 2019; Oppenheimer et al. 2020). Supernova wind failure allows central gas to accumulate (e.g., Dekel & Silk 1986), which fuels rapid BH growth and causes quenching. Dekel et al. (2019) even argue that these two masses converge to create a single “golden mass” at early times.

In contrast, the constancy of $\langle M_{*,GV} \rangle$ in our model has nothing to do with $M_{crit}$ or supernova winds but comes instead
from the fact that halos and galaxies evolve so as to keep the quenching boundary located in the same place relative to the blue cloud over time. This means that quenching galaxies always exit the blue cloud with the same characteristic mass. Granted, it is fortunate that $M_{crit}$ for hot halos is near $10^{11.5-12.0} M_{\odot}$ and not larger, because having a hot halo may facilitate BH quenching. Central supernova wind failure might also do this. However, we equate quenching with cross over into the green valley, and, in the binding-energy model, this point is not set by the ability of halo gas to cool, or by the failure of supernova winds to keep galaxy centers gas free, but by the cumulative effect of BH feedback on the halo gas. The characteristic masses set by halo cooling and by supernova winds are not fundamental—they merely need to be small enough to set the stage for the halo binding-energy criterion, which is determinative.

It may be asked if there is a characteristic halo mass in the binding-energy model, and the answer is yes. This is the halo mass, $M_{vir,cross}$, where the growing energy from the BH in a typical galaxy (with $R_e = \langle R_e \rangle$) during the star-forming phase crosses the binding energy of halo gas. Scaling laws in Appendix C may be combined with the quenching condition in Equation (23) to yield

$$ \left( \frac{M_{\text{vir,cross}}}{10^{12} M_{\odot}} \right)^{1.60v-1.67} = 9.5 k_{\text{Re}}^{v} k_{\text{BH}}^{-1} k_{\Sigma_1}^{-v} k_{\text{Rvir}}^{-1} k_{\text{m}}^{-1.1v} \frac{f_{\text{B}}}{4} \times \frac{f_{\text{hot}}}{0.1} \times \left( \frac{\eta}{0.01} \right)^{-1} \times h(z)^{0.74(v-1)}, \quad (24) $$

where the constants $k$ are the zero points of the five relevant scaling laws from Appendix C and we have inserted our best values for $p$, $q$, $s$, and $t$. $M_{vir,cross}$ will be a fundamental mass if the constants $k$ are fundamental numbers. The first two of the constants are already known, and it will be possible to write down the values of the other three once galaxy radii and BH building are fully understood. Two more free parameters at the present time are $f_B$ and $\eta$. These also should be expressible (someday) using fundamental constants.

For reference, we insert $k$-values from Appendix C, use $v = 1.76$, and substitute $M_*$ for $M_{vir}$ to get the crossing mass in terms of $M_*$, which should equal to $M_*$ in the middle of the green valley. We obtain

$$ \frac{M_{*,\text{cross}}}{10^{10} M_{\odot}} = 4.0 \left(\frac{f_{\text{B}}}{4}\right)^{1.54} \times \left(\frac{f_{\text{hot}}}{0.1}\right)^{1.54} \times \left(\frac{\eta}{0.01}\right)^{-1.54} \times \left(\frac{h(z)}{2.3}\right)^{-0.86}, \quad (25) $$

where $h(z)$ has been normalized by its value at $z = 1.5$, where $k_{Rvir}$ is calculated. $M_{*,cross}$ therefore equals $10^{10.6} M_{\odot}$ at $z = 1.5$, which is a good match to the green-valley stellar mass. However, note that $M_{*,cross}$ is not exactly constant with time due to the presence of the term $h(z)^{-0.86}$. This term arises from the same zero-point mismatch between the blue lines and the red lines in Figure 1 that was pointed out in Section 7.1. However, the residuals relative to $z = 1.5$ in mass are -0.2 dex at $z = 3$ and +0.2 dex at $z = 0.5$, which are not large.

To summarize, the binding-energy model, like other models, depends on supernova feedback to keep galaxies on the SMHM

relation while star-forming. For quenching, some form of AGN feedback is required, which could be either high mode or low mode. Because BH energy overcomes halo gas-binding energy in this model, preventive feedback is probably the dominant mode though ejective feedback may also play a role. Halo quenching may be a prerequisite for the binding-energy criterion to work, but it does not set the fundamental green-valley mass. Late-infall dry-up (and the transition from cold streams to hot streams) may make late halos easier to quench than early halos, thereby contributing to the faster lowering of the quenching boundary zero point (compared to the model prediction) that is observed in $\Sigma_1-M_*$ in Figure 1.

### 7.4. Compaction in the Binding-energy Model

Implicit in the binding-energy picture is that galaxies, halos, and BHs evolve smoothly and steadily together during the star-forming phase and that crossing into the green valley is a quiet event with no notable fireworks. That seems a good model for quenching today, particularly for S0 galaxies, which have not experienced recent major mergers because their disks are intact. However, it will be noticed that many dense star-forming galaxies are visible on the wrong side of the boundary at high redshift but then disappear at low redshift. This is especially true using $R_e$ as a measure of density as opposed to $\Sigma_1$ (see red points below the boundary in the top row of Figure 1(a)). The topic of this section is how these dense, discrepant galaxies might be accommodated in the binding-energy model.

These small, dense galaxies are the so-called “blue nuggets” at $z = 2$ (Barro et al. 2013, 2014; Nelson et al. 2014; Williams et al. 2014; van Dokkum et al. 2015; B17). They are thought to be the star-forming precursors of the similarly dense, quenched “red nugget” galaxies (Daddi et al. 2005; Trujillo et al. 2006; van Dokkum et al. 2008; Damjanov et al. 2011), which are visible in Figure 1 at the same redshifts. VELA simulations (Ceverino et al. 2014) suggest that both populations originated from a “compaction” event in which gas was driven strongly to the center, forming a burst of stars that consumed the central gas and drove gas from the rest of the galaxy (Dekel & Mandelker 2014; Zolotov et al. 2015; Tacchella et al. 2016b). This scenario says that star-forming galaxies compactified to form blue nuggets, which then faded to become red nuggets. Simulations by Choi et al. (2018) and by the IllustrisTNG team (Habouzit et al. 2019) produce similar objects.

Actually, VELA galaxies never fully quench owing to their lack of central BHs. But X-ray data show a large population of X-ray sources near the quenching boundary at $z \sim 2$ (Kocevski et al. 2017), right where maximum compaction occurs (Dekel et al. 2019). If these accreting BHs were added to VELA, their feedback might achieve quenching. In contrast to VELA, BHs are present in IllustrisTNG and in the Choi et al. (2018) models, and in these simulations compact blue nuggets quench.

A variety of factors can trigger compaction events in central galaxies, including mergers (Barnes & Hernquist 1991; Hernquist 1991; Hopkins et al. 2006), violent-disk-instability-driven inflows (Dekel & Mandelker 2014), counterrotating streams (Danovich et al. 2015), a triaxial halo core (Tomassetti et al. 2016), and the return of recycled low-angular-momentum gas (Elmegreen et al. 2014). However, most of these events are stochastic in nature and uncorrelated with the underlying galaxy structural parameters. Galaxies therefore tend to evolve along their native evolutionary trajectories in these scenarios, their parameters bobbing up and down due to occasional
compaction events until a final, major drop in star formation
occurs. Such behavior is seen in VELA, where the final drop is
driven by late-infall dry-up (Tacchella et al. 2016a).

We take this as the basic picture for compaction in the
binding-energy model, to which we add BHs. As blue nuggets
compactify, their central densities rise rapidly and their BHs
grow according to the green-valley vectors in Figure 4.
Compaction thus sweeps galaxies rapidly over the quenching
boundary at high redshift, at which point their growing BHs
begin to impede star formation according to Equation (23).

Adding BHs in this way solves four problems. First, the inclusion of BHs makes galaxies quench. Second, quenching according to the binding-energy equation produces sloping quenching boundaries, which are seen at high redshift as well as at low redshift (Figure 1). Third, the model sharpens the rather fuzzy boundary that might otherwise result from compaction alone due to its stochastic nature. Finally, the disappearance of blue nuggets at low redshift is natural since compaction is strong in gas-rich galaxies but should tend to taper off as galaxies become gas poor.

The one remaining issue is timing: why does star formation
continue well after the BH mass has crossed the quenching
boundary? We theorize that this occurs because compaction
delivers an unusually large amount of fresh gas to the center of
a galaxy in a short space of time. In order for the galaxy to
quench, this gas must either be consumed or ejected, which
imposes a delay. During that delay, the galaxy sits on the
wrong side of the boundary, being compact yet still star-
forming. After the central gas has cleared, the recently enlarged
BH can transmit energy and/or momentum all the way to the
halo, as required by the binding-energy picture.

These last remarks echo the “fast track–slow track”
dichotomy for galaxy quenching that has been advanced by
several authors (e.g., Yesuf et al. 2014; Woo et al. 2015; Barro
et al. 2017; Wu et al. 2018; Woo & Ellison 2019). Fast-track
galaxies, sometimes called “post-starburst galaxies” (Dressler
& Gunn 1983; Balogh et al. 1999; Dressler et al. 1999), exhibit
strong Balmer absorption lines and small radii, which indicate
rapid quenching and strong gaseous dissipation (Whitaker et al.
2012; Yesuf et al. 2014; Yano et al. 2016; Almaini et al. 2017;
Wu et al. 2018). The binding-energy model attributes both
characteristics to the rapid arrival of gas at the center. A
starburst ensues while radii are small, which places them
temporarily on the wrong side of the quenching boundary until
BHs can grow and quench. Slow-track galaxies in contrast slide
over the boundary gradually, and their star formation rates
never get strongly out of sync with their structural properties
and BH masses. Because the fast-track scenario requires lots of
gas, the fraction of fast-track to slow-track galaxies should
decrease with time, as observed (Whitaker et al. 2012; Yano
et al. 2016; B17). In conclusion, although fast-track galaxies
are not explicitly included in the binding-energy model, modest
extensions to it seem able to accommodate them.

7.5. *Evolution through the Green Valley*

The manner in which galaxies evolve through the green
valley is important in the model because the slope of fully
quenched galaxies is assumed for BH growth in star-forming
galaxies (Section 3.1), and these two populations are connected
by the green valley. This section seeks to understand why
green-valley BH growth vectors might be both steep and
uniform over time (see Figure 4).

Evidence to support the steepness and uniformity of green
valley tracks was briefly described in Section 3.1, which we
expand on here. The first fact is the sharp “elbow” seen in
ΔlogSSFR versus Δlog Σ₁ by F13 and again in Figure 6(b)
of B17, in which SSFR declines steeply through the green
valley while Σ₁ hardly changes (see cartoon version of this in
Figure 5(a) here). It is this pileup at near-constant Σ₁ that
creates the narrow ridgeline of quenched galaxies in Σ₁–M*
(see Figure 1). The constant morphology of the elbow
(F13; B17; Lee et al. 2018; Mosleh et al. 2018) sets the
growth of log Σ₁ through the green valley at roughly 0.2 dex at
all masses and redshifts z ≤ 3. Consistent with this, the green-
valley entry locus also lies a constant 0.2 dex below the
quenched ridgeline for SDSS galaxies at all masses today (Luo
et al. 2020). In summary, there seems to be a universal Σ₁
“elbow diagram” such that the growth of Σ₁ through the green
valley is small, only about ~0.2 dex, while the SSFR falls by
much more than this (F13; B17; Luo et al. 2020).

The evolution of BH mass through the green valley is
different. Terrazas et al. (2017, Figure 3) have shown that BH
mass continues to grow through the green valley even as star
formation declines. That is, ΔlogSSFR versus Δlog M<sub>BH</sub> does
not show a sharp elbow but rather a smooth, monotonic
decline, with M<sub>BH</sub> increasing across the green valley by about
1.5 dex (see Figure 5(b) here). The same difference is also seen
in Figure 2(c), where star-forming and quenched galaxies are
well separated, in contrast to Σ<sub>1</sub> in Figure 2(b), where the two
populations overlap. Further evidence for large BH growth in
the green valley comes from the enhanced X-ray luminosities
of galaxies as they evolve off the star-forming ridgeline and
into the green valley; these show a remarkably similar pattern
at all redshifts in which L<sub>X</sub> rises by a factor of 2–3 compared to
the star-forming main sequence (Aird et al. 2019, Figure 10(c)).
If X-ray luminosity accurately tracks BH accretion (but see
below), these data suggest high BH growth rates in the green
valley for all quenching galaxies since z = 4. Rapid BH growth
in the green valley based on X-ray fluxes has also been
deduced by Delvecchio et al. (2019), and a similar picture in
which BHs accrete most of their mass just before quenching
has been modeled by Caplar et al. (2018).

To summarize, in the binding-energy model, Σ₁ virtually
stops growing as galaxies enter the green valley whereas M_BH
grows considerably. This is the behavior captured by the steep
green tracks in Figure 4. The different tracks from left to right
represent differences in R_c during the star-forming phase:
galaxies with smaller initial R_c enter the green valley sooner at
lower values of Σ₁ and M_BH and vice versa. As the quenching
boundary moves down with time, the range of green vectors
being populated at any redshift also moves down, as shown by
the moving green parallelogram.

The late rapid BH growth depicted in the green valley in
Figure 4 resembles BH growth in the following hydro models:
RAMSES/Horizon-AGN (Dubois et al. 2015; Habouzit et al.
2017; Dekel et al. 2019), EAGLE (Bower et al. 2017), FIRE
(Anglés-Alcázar et al. 2017), IllustrisTNG (Habouzit et al.
2019), and MARVEL (Bellovary et al. 2019). BH growth in
these models is retarded during star formation due to supernova
feedback, wandering motions of BHs in shallow central
potential wells, and low central gas densities. Retardation ends
when supernovae are no longer able to evacuate the central
regions (for various reasons) and when the growing central
stellar surface density captures the BH and concentrates the gas
around it. These conditions will also occur in the binding-energy model, because declining star formation at entry to the green valley means less supernova feedback, and higher $\Sigma_1$ means a deepening central potential well. In short, the steep BH growth rates in the green valley in the binding-energy model find some support from similarly growing BHs near the end of star formation in hydro models.

In fact, the growth of BHs in the hydro models is so retarded that BH masses lie below the blue line in Figure 4 for most of the star-forming phase, crossing it only later and approaching it from below. Possible tracks are schematically indicated by the steep dashed green lines below the blue line in that figure. This scenario has important implications, as it changes the blue line from an evolutionary track into a boundary that is crossed by galaxies that are moving upward in Figure 4 into the green valley. This would profoundly alter the meaning of the $M_{BH} \sim \Sigma_1^\nu$ scaling law and is discussed further in Section 7.8.

In the above, we have interpreted both loci in Figure 5 as evolutionary tracks that galaxies traverse while evolving from the blue cloud through the green valley and onto the red sequence. Panel (b), showing BH growth, is a re-rendering of data from Terrazas et al. (2017). However, these authors favored a different interpretation in which green-valley galaxies can remain for long periods of time in a quasi-equilibrium state, with BH feedback regulating their star formation at a fairly constant level. We favor the more rapid evolutionary track interpretation because of the observed net flow of galaxies from the blue cloud to the red sequence over cosmic time (Bell et al. 2004; Faber et al. 2007; Huang et al. 2013; Ilbert et al. 2013; Muzzin et al. 2013; Tomczak et al. 2014; Straatman et al. 2015; F18). On the other hand, the Terrazas sample probably does contain galaxies that quenched at different times, and hence this locus is not the track of any single object. If the quenching boundary has been moving down over time, earlier-quenched galaxies would lie too far right in this diagram, and the amount of BH growth would look too large. The total amount of the boundary motion was estimated to be ~0.5 dex in $M_{BH}-M_*$ from $z=3$ to now in Section 7.1. Even if this full correction were applied to the fully quenched galaxies in Figure 5, we would still have a net BH growth factor across the green valley of 1.0 dex, or $\times 10$. We take this as our estimate of the minimum amount of BH growth in the green valley in what follows.

We can now use this growth factor to estimate how much energy BHs emit during their stay in the green valley. Galaxies enter the green valley when their effective BH energy equals $f_B \times$ their halo gas-binding energy, where $f_B$ is set to 4 to match the zero points of the blue lines to the red lines in Figure 1. We assume for the sake of argument that $M_*$ grows like $\Sigma_1$ through the green valley (0.2 dex) and that $M_{vir}$ grows similarly, galaxies being near the knee of the SMHM relation at this point, where the log slope is near unity (Figure D1). Then, $E_{bind} \sim M_{vir}^{5/3}$ grows by $5/3 \times 0.2 = 0.3$ dex (Equation (21)), and the ratio $E_{BH}/E_{bind}$ grows by at least $1.0 - 0.3 = 0.7$ dex, which is a factor of 5. This together with $E_{BH}/E_{bind} = 4$ on entry to the green valley then implies that $E_{BH}/E_{bind} \gtrsim 20$ when galaxies come to a halt on the quenched side. The factor of 20 has margin to accommodate as-yet-unknown radiative losses and other inefficiencies and is plausibly large enough to complete quenching. If green-valley vectors are

indeed universal, this new quantity $E_{BH}/E_{bind} \sim 20$ at arrival on the quenched ridgeline is yet another constant of quenching.

The foregoing discussion shows that postponing major BH growth to the green valley has important consequences for how BHs accrete their gas, generate AGNs, produce feedback, and impact their environment. We therefore need to review the evidence carefully, especially since the Terrazas et al. (2017) sample of BHs is rather small, highly biased toward massive galaxies (as shown in Figure 2(c)), and reflects only BHs today. We start by listing seven different pieces of evidence that do favor uniform and steep BH growth in the green valley over extended cosmic times:

(1) A. Rodríguez-Puebla et al. (2020, in preparation) have assembled a sample of reliable BH masses that is considerably larger than the sample in Terrazas et al. (2017). Using $\Delta\log\text{SSFR}$ to locate a galaxy's position in the green valley, as predicted by Figure 4, they derive a single, unified scaling law that fits all BH masses. The coefficient of $\Delta\log\text{SSFR}$ in the new law is consistent with the slope and length of the green-valley vectors in Figure 4.

(2) The zero points of the $M_{BH}-\sigma_c$ and $M_{BH}-M_*$ relations for star-forming galaxies are roughly 1 dex below those of quenched galaxies (KH13; Ho & Kim 2014; Reines & Volonteri 2015; Savorgnan et al. 2016; Terrazas et al. 2016; Davis et al. 2018). These offsets make sense if it is assumed that BHs grow steeply in the green valley as galaxies pass between the two relations.

(3) A long-standing puzzle has been the very different BH masses of the Milky Way ($4 \times 10^6 M_\odot$) and M31 ($1.5 \times 10^8 M_\odot$; KH13) despite the small difference in their stellar masses of only 0.2–0.3 dex (Sick et al. 2014; Licquia & Newman 2015). No standard BH scaling law has been able to fit both galaxies. The two objects are plotted in Figure 4. M31 lies close to the quenched locus while the Milky Way lies slightly below the star-forming locus, and the slope of the line connecting them is roughly parallel to the adopted growth vector in the green valley. Thus, most of the discrepancy between the two galaxies would be resolved if the Milky Way is just entering the green valley and M31 is just exiting it. This is consistent with the fact that KH13 classify the Milky Way as a pseudo-bulge (i.e., gas-rich and star-forming) but classify M31 as a classical bulge.

(4) At fixed $M_*$, AGN frequency rises strongly in the green valley. This is true for SDSS galaxies (e.g., Schawinski et al. 2007, 2010; Guo et al. 2019) and also for CANDELS galaxies (e.g., Nandra et al. 2007; Gu et al. 2018; Aird et al. 2019).

(5) RP17 used abundance matching to measure evolutionary tracks for $\Sigma_1$ and $M_{BH}$ back in time. A. Rodríguez-Puebla et al. (2020, in preparation) now find, using similar methods, that BHs in quenched galaxies follow evolutionary tracks that look much like the broken power laws in Figure 4.

(6) In Section 7.7, we show that steep BH growth through the green valley (like that in Figure 4) can explain the different BH scaling laws measured by different authors for different galaxy samples, including the very steep slope of the BH “fundamental plane” found by van den Bosch (2016) when fitting all galaxies.

(7) Lastly, there is the theoretical justification from the hydro simulations cited above, which predict that BH growth is unleashed when bulges reach a threshold size, resulting in steep growth trajectories similar to the green vectors in Figure 4.
On the other hand, not all evidence supports sufficiently rapid BH growth in the green valley to satisfy the model. Most troubling are the X-ray AGNs analyzed by Aird et al. (2019), which show only a modest increase in X-ray luminosity as galaxies enter the green valley (2×–3×). If BHs grow by at least ×10 in the green valley, they accrete 10% of their mass on the main sequence and 90% of their mass in the green valley (assuming no further growth after quenching). If main-sequence and green valley lifetimes are similar (which, if anything, probably overestimates time in the green valley), average X-ray fluxes should be roughly an order of magnitude higher in the green valley, not the factor of 2–3 observed. In short, AGNs should not only be prevalent in the green valley, they should also be very luminous X-ray sources, which is not seen.

There appear to be only two ways out. First, the bulk of the X-ray flux might be obscured. Indeed, the sample studied by Aird et al. (2019) has not been corrected for Compton-thick sources with $N_H > 10^{24} \text{ cm}^{-2}$ (Aird et al. 2018). Obscured AGNs are reviewed by Hickox & Alexander (2018), who conclude that their number could be large. The second option is that BH mass accretion in the green valley is radiatively inefficient, i.e., that BHs can accrete much mass without shining in X-rays. This is thought to occur in advection-dominated accretion flows (Yuan & Narayan 2014), but their mass-accretion rates are low whereas our green-valley sources would be accreting fast. The point is that a switch in mode would need to occur as galaxies enter the green valley, whereby the ratio of radiative energy per unit BH mass accreted would need to decline. Whether such a change, or obscuration, or both together can solve the problem of missing X-ray flux in the green valley is unknown. Right now, missing X-ray flux in the green valley appears to be the biggest challenge to the binding-energy model.

Finally, we note that the tight relationship between BH mass and halo mass that exists at the end of quenching, where $E_{BH} \sim 20 \times E_{bind}$, must break down later as the masses of central galaxies and their BHs stagnate while halos continue to grow. The BH-halo relationship in different stages of a galaxy’s lifetime is discussed in the next section.

## 7.6. The Black Hole/Dark Halo Connection

The nature of the link between BHs and halos is controversial. A tight correlation between $M_{BH}$ and $M_{vir}$ was claimed by Ferrarese (2002), Baes et al. (2003), and Bandara et al. (2009) but was disputed by KH13 because neither very large nor very small halos fit their relations. The binding-energy model sheds further light.

The first point is that the different scaling laws in Figure 4 clearly imply different correlations for different kinds of galaxies. The equation $\Sigma_1 \sim M_*^{0.88}$ is the observed ridgeline slope of star-forming galaxies at fixed epoch in $\Sigma_1-M_*$ (B17). Converting this to BH mass using $M_{BH} \sim \Sigma_1^{1.76}$ yields $M_{BH} \sim M_*^{1.55}$. This is steeper than the ridgeline for quenched galaxies ($M_{BH} \sim M_*^{1.16}$), lies below it, and should show more scatter because star-forming galaxies have different $R_e/\langle R_e \rangle$, and therefore different BH masses at fixed $M_*$. These predictions agree tolerably well with observations (KH13; Reines & Volonteri 2015; Savorgnan et al. 2016; Terrazas et al. 2016; Davis et al. 2018; Delvecchio et al. 2019). Finally, using $M_* \sim M_{vir}^{1.75}$ for star-forming galaxies predicts $M_{BH} \sim M_{vir}^{2.71}$, but with large scatter due to $R_e$ variations.

For quenched galaxies, the derived slope depends on whether the quenched galaxies are in halos above or below the knee of the SMHM relation, which is at $M_* = 10^{10.7} M_\odot$ (Figure D1). As most quenched ridgeline galaxies are more massive than this (KH13), we start with the quenched ridgeline slope $M_{BH} \sim M_*^{1.16}$ and use the SMHM relation $M_* \sim M_{vir}^{0.35}$ above the knee to obtain $M_{BH} \sim M_{vir}^{0.40}$ for quenched galaxies, i.e., dramatically shallower than the star-forming slope of 2.17. The main reason for the difference is the flattening of the SMHM relation above the knee at $M_{vir} \sim 10^{12.5} M_\odot$. The net result is that BH mass is predicted to rise strongly in star-forming galaxies with halo mass below $M_{vir} = 10^{12.5} M_\odot$ but much more slowly above that, as BHs stagnate but halos continue to grow. Hence, no single power-law relation is expected between BHs and their halos, and the scatter among star-forming galaxies (and in the green valley) should be large.

In contrast, Ferrarese (2002) found a single tight relation for all galaxies with slope 1.65–1.82 by estimating $M_{vir}$ from the outer rotation curves of spirals in the range $\log M_{vir} = 11.5–12.5$ and the outer kinematic profiles of ellipticals (Kronawitter et al. 2000) in the range $\log M_{vir} = 12.5–13.5$. Baes et al. (2003) found a single slope of 1.27 from a larger sample of spirals but the same ellipticals from Kronawitter et al. (2000). Finally, Bandara et al. (2009) found a slope of 1.55 using strong lensing to estimate halo masses of mainly ellipticals in the range $\log M_{vir} = 13.0–13.8$. None of these results agrees with the model predictions, the slopes being too uniform and equal to some compromise value between the two predicted extremes of 0.40–2.71.

The steep slopes found for the massive end in these papers clearly disagree with elementary estimates. Large spirals like M31 have $M_{BH}$'s a few times $10^8 M_\odot$ and dark halos a few times $10^{12} M_\odot$, while large clusters like Coma have central $M_{BH}$'s a few times $10^9 M_\odot$ but dark halos a few times $10^{15} M_\odot$. The differential log slope is about one-third, in good agreement with the predicted binding-energy slope of 0.40 but vastly smaller than the measured slopes of 1.27–1.82 from these authors. Similar points were made by KH13 and Delvecchio et al. (2019). The problem may be trying to estimate the total mass of huge dark halos from data that do not extend far enough from the center—the resulting values of $M_{vir}$ are too small, and the calculated slopes of $M_{BH}$ versus $M_{vir}$ are too steep.

To summarize, it is not true that there is no relation between the mass of a halo and its central BH. Indeed, it is a fundamental tenet of the binding-energy model that halos build galaxies, which in turn predictably build BHs. However, the resulting slopes and scatter are very different for quenched and star-forming galaxies. This chain of reasoning answers the fundamental question posed by F13, who guessed that a close coordination between halo and BH is needed in order to produce a narrow quenched $\Sigma_1-M_*$ ridgeline but wondered how the halo “knows” what the BH is doing. The answer is that the halo makes the BH during the star-forming phase and then is quenched by it when $E_{BH} \sim 20 \times E_{bind}$ (Section 7.5), closing the loop. BH mass differences due to $\Delta \log R_e$ cause extra scatter about the $M_{BH}-M_{vir}$ relation for star-forming and green valley galaxies but do not fundamentally change this basic picture.
## 7.7. Black Hole Scaling Laws and Fundamental Planes

This section returns to the problem mentioned in Section 3.2 that the value of $v$ in the model seems too small. Specifically, our value $v = 1.76$ implies that $M_{BH} \sim \sigma_e^{3.52 \pm 0.26}$ (if $\sigma_1 \sim \sigma_e$) whereas the observed relation (for classical bulges and ellipticals in KH13) is $M_{BH} \sim \sigma_e^{4.4 \pm 0.3}$. This suggests that $v$ is too small by the factor $4.4/3.52 \simeq 1.25 \pm 0.12$. van den Bosch (2016) fitted all galaxies together and found an even steeper law, $M_{BH} \sim \sigma_e^{5.35 \pm 0.23}$, for a difference with ours of $5.35/3.52 = 1.5 \pm 0.13$. We do not have measurements of $\sigma_e$ in the CANDELS database, but we do have measurements of $M_*$ and $R_e$, which were also fitted by van den Bosch (2016). He found $M_{BH} \sim M_*^{2.91 \pm 0.14} R_e^{-2.77 \pm 0.22}$ from fitting all galaxies, which is close to $M_{BH} \sim \Sigma_1^{2.7 \pm 0.17}$ using the *p*, *q* mapping. This implies that $v$ is too small by the factor $2.7/1.76 \simeq 1.5 \pm 0.14$, i.e., the same as $\sigma_e$. (Note: this last relation also implies that $M_*/R_e \sim \sigma_e^2$.0.) The goal of this section is to understand these slope differences.

We have taken advantage of these rough correspondences between $\Sigma_1$, $\sigma_1$, $\sigma_e$, and $M_*/R_e$ to label the X-axis in Figure 4 with all four quantities, which allows us to crudely compare slope results from different data sets. We stress that $v$ is our slope for star-forming galaxies only (blue line in Figure 4), whereas van den Bosch (2016) lumps all galaxies together and fits them with a single power law. Being an average over both the shallow and steep portions of the tracks, his slope should be steeper than the blue track in Figure 4. A rough estimate is shown by the black dashed line in Figure 4, which is schematically fit to all galaxies. The slope of the black line is 2.6, very close to the van den Bosch (2016) value of 2.7 above. The construction is crude but illustrates clearly how major BH mass growth through the green valley, accompanied by little or no growth in $\Sigma_1$, can produce a big increase in the average slope of $M_{BH}$ versus $\Sigma_1$ (or $\sigma_1$, or $M_*/R_e$) when all galaxies are fitted together.

The difference between our slope and that of KH13 for E's and classical bulges is only 1.25 but should be 1.00 if their sample is pure. However, at the low-mass end, some of the classical bulges in KH13 may actually be in the green valley, and at the high-mass end, Lauer et al. (2007) have argued that $M_{BH}-\sigma$ has been steepened due to galaxy merging. Contamination by both of these effects would steepen the slope in the direction seen. (We note that the derivation of *t* and *v* in Section 3.2 deliberately chose the $M_{BH}-M_*$ relation, not $M_{BH}-\sigma$, for quenched galaxies, precisely because it is less sensitive to merging.)

The discussion above has focused on scaling laws versus one variable, typically $\Sigma_1$ or $\sigma$. However, the close relation of these quantities to $M_*/R_e$ has led some authors to fit using two variables from among $M_*$, $R_e$, and $\sigma_e$. These expanded correlations are often referred to as the BH fundamental plane (BHFP). van den Bosch's planar fit to all galaxies yielding $M_{BH} \sim M_*^{2.91} R_e^{-2.77}$ is an example. However, he noted his opinion that the fundamental variable is really $\sigma_e$ and that $M_*$ and $R_e$ are stand-ins for it. In that sense, the BH "fundamental" plane is not as fundamental as the relationship with $\Sigma_1$ (or with $\sigma$).

Separately, Hopkins et al. (2007a) fitted a plane to quenched ridgeline galaxies only and found $M_{BH} \sim M_*^{1.78 \pm 0.40} R_e^{-1.05 \pm 0.37}$. As mentioned in Section 3.1, the *p*, *q* mapping for quenched galaxies is empirically $\Sigma_1 \sim M_*^{1.0} R_e^{-0.5}$. Our model then predicts

$M_{BH} \sim \Sigma_1^{1.76} \sim M_*^{1.76} R_e^{-0.88}$ on the quenched ridgeline, which agrees with the Hopkins exponents within their errors. We therefore propose that the BHFP for quenched galaxies, like the plane of van den Bosch (2016) for all galaxies, is really a stand-in for $\Sigma_1$ (or $\sigma_e$). Similar calculations (not shown) can explain the slopes of BHFPs fitted to $\sigma_e-R_e$ and $\sigma_e-M_*$ by Hopkins et al. (2007a) and by Aller & Richstone (2007).

Yet another expression offered for the quenched BHFP is $M_{BH} \sim (M_*\sigma_e^2)^{\alpha}$, where Hopkins et al. (2007a, 2007b) find $\alpha \simeq 0.7$ and Aller & Richstone (2007) find $\alpha \simeq 0.6$. The quantity $M_*\sigma_e^2$ is proportional to galaxy binding energy per particle, which these authors proposed as the property that regulates BH mass. This was motivated by the merger scenario of Hopkins et al. (2007b), in which gas that is brought to the center of the galaxy builds the BH until the BH energy unbinds the surrounding gas. In spirit, this is similar to the unbinding criterion in Equation (23), but we use gas-binding energy to the halo, $M_{vir} V_{vir}^2 \sim M_{vir}^{5/3}$, whereas Hopkins-Aller use gas-binding energy to the galaxy.

A problem noted in both works is that $\alpha$ ought to be 1.0 in a simple binding-energy argument, not 0.6-0.7 as seen. However, the binding-energy model again sheds light. From Section 5, halo gas-binding energy $E_{bind} \sim M_{vir} V_{vir}^2 \sim M_{vir}^{5/3}$. Setting $\sigma_e^2 \sim M_*/R_e$, noting that $R_e \sim M_*^{0.44}$ along the quenching boundary, and substituting $M_* \sim M_{vir}^{1.75}$ on the SMHM relation yield $M_{vir}^{5/3} \sim (M_*\sigma_e^2)^{0.61} \sim E_{bind}$. In other words, the unexpected exponent with the value of 0.6-0.7 in the model of Hopkins et al. (2007b) and Aller & Richstone (2007) is simply that needed to convert halo gas-binding energy into galaxy gas-binding energy!

In summary, fitting BH scaling relations to different samples has historically yielded slopes that differ by more than their formal errors. This is true whether single-variable or planar fits are used. The broken-power-law model in Figure 4 helps to explain this by positing different scaling laws (and $R_e$ scatter) for different evolutionary phases. Fitting any single law to the entire population gives an average. Better fits can be obtained by fitting each evolutionary stage separately or by including a new term in $\Delta\log\text{SSFR}$ to reflect a galaxy's location in the green valley (see Item 1 in Section 7.5). Finally, it appears that planar fits to two variables are really stand-ins for fits to $\Sigma_1$ (or $\sigma_e$), which are the fundamental quantities that describe $M_{BH}$.

## 7.8. Black Hole–Galaxy “Coevolution”: the Origin and Meaning of Black Hole Scaling Laws

Because of its highly transparent nature, the binding-energy model is a good tool for probing the nature of the so-called BH–galaxy “coevolution.” In our opinion, this word has been used too loosely by astronomers simply to signify that the properties of two objects (e.g., BHs and galaxies) are correlated. What we really want to know is the nature and arrow of causality: does A cause B, does B cause A, are A and B both caused separately by C, or do A and B reciprocally influence one other throughout time? The dictionary defines coevolution as “the reciprocal influence of one evolving entity on another evolving entity,” and we adopt that standard here.

It is then illuminating to compare the BH growth processes in the binding-energy model to those in other models, and detailed descriptions of BH growth mechanisms exist in the literature for the following hydro simulations: RAMSES/Horizon-AGN (Dubois et al. 2015; Habouzit et al. 2017;
Dekel et al. 2019), EAGLE (Bower et al. 2017; Davies et al. 2019; Oppenheimer et al. 2020), and FIRE (Anglés-Alcázar et al. 2017). We do not include IllustrisTNG because it uses two different BH feedback modes (thermal-quasar mode and kinetic-radio mode; Weinberger et al. 2018; Habouzit et al. 2019; Nelson et al. 2019), which would complicate the discussion.

Three different evolutionary phases are considered:

(1) **BH growth on the star-forming main sequence:** In the binding-energy model, all star-forming galaxies evolve along the BH growth law $M_{BH} = \text{const.} \times \Sigma_l^{1.76}$ (blue line in Figure 4), which is an evolutionary track. Growth is steady on average but may be episodic on short timescales. Because $M_*$ and $R_e$ determine $\Sigma_l$ through the p, q mapping, BHs are rigidly linked to their host galaxies during this phase. If halos also determine $R_e$ (through halo mass and concentration, for example), BHs are rigidly linked to their halos as well. In short, the properties of galaxies are determined by halos and the properties of BHs are in turn determined by galaxies, but there is no reciprocity at either level and so, strictly speaking, there is no coevolution.

In the above-cited hydro models, BH growth during the star-forming phase is severely limited by bursty stellar feedback, which periodically ejects gas from the center. Growth overall is episodic and inefficient. Instantaneous BH growth takes place by Bondi accretion (RAMSES, Horizon-AGN, EAGLE) or is torque limited (FIRE), but the time-averaged rate is more likely to be limited by the rate at which gas refills in the central regions. The central stellar density is low, the central potential is not very deep or well focused, and BH wandering is strong (Bellovary et al. 2019; Dekel et al. 2019), which may further reduce BH gas accretion. BHs are typically started with seed masses of around $10^5 M_\odot$, and growth beyond that in this phase is typically only 3× even though galaxy masses increase by up to 100×. Owing to the different physical processes governing central BHs and host galaxies, their properties during the star-forming phase are not closely linked, and we conclude that BHs and galaxies do not coevolve—they are not even closely correlated.

(2) **Entry to the green valley; transition to rapid BH growth:** In the binding-energy model, both of these phenomena occur together, when the BH energy output passes $4 \times E_{\text{bind}}$. The BH has been growing steadily and is just beginning to perceptibly reduce the ability of gas in the halo to fall onto the galaxy. Keying the location of this transition to a fixed fraction of halo binding energy closely matches the observed, sloping quenching boundaries (blue versus red lines in Figure 1). Galaxy, halo, and BH properties are perfectly correlated at the transition, as in the previous star-forming phase, but now BHs are backreacting on their halos and galaxies. This reciprocal influence satisfies the definition of coevolution, and so BHs, galaxies, and their halos are coevolving as galaxies enter the green valley.

In the cited hydro models, the transition to rapid BH growth occurs before entry into the green valley, when the stellar-feedback-driven wind fails. The stated reason is different in the different models (failure to reach escape velocity (Dekel et al. 2019), winds no longer buoyant (Bower et al. 2017; McAlpine et al. 2018), major mergers (McAlpine et al. 2018), and deepening central potential well (Anglés-Alcázar et al. 2017; Dekel et al. 2017)), but the net result is that gas starts to collect in the center of the galaxy, unleashing BH growth. This

transition tends to occur in a small range of bulge mass near $10^{9-10} M_\odot$, a halo mass near $10^{12} M_\odot$, or a fixed halo virial temperature $T_{\text{vir}}$ near $10^{6-7}$ K, all of which are plausibly related to the kinetic energy content of supernova winds. Because this is a fixed number (Dekel & Silk 1986), the theory tends to pick out a single stellar or halo mass where BH growth is launched, thus creating a momentary correspondence between the BH, galaxy, and halo but not true coevolution as there is no reciprocity at this point.

(3) **BH growth in the green valley; final transition to fully quenched:** In the binding-energy model, the BH growth trajectory versus $\Sigma_l$ in the green valley steepens. The log slope of the blue track in Figure 4 is 1.76, whereas the log slope of the green valley vectors is 1.5 dex/0.2 dex = 7.5 (for these numbers, see Figure 5; with minimum BH growth of 1.0 dex in the green valley, the ratio would be 5.0). The steep tracks needed by the model closely resemble the steep growth tracks in the hydro simulations, but the detailed physics shaping those tracks is not understood. Meanwhile, global gas content and star formation are falling due to AGN feedback. A reasonable inference, then, is that BH, galaxy, and halo are continuing to “talk” to one another, and this coevolution is plausibly what makes all green-valley vectors parallel and of the same length. Ultimately, galaxies arrive on the fully quenched ridgeline when $E_{\text{BH}}/E_{\text{bind}} \simeq 20$, which enforces absolute correspondence between halos and BHs again at that point. After that, halos continue to gain mass and drift away from both BHs and galaxies, and coevolution ends (except possibly for “maintenance mode”).

In hydro models, it is less certain what sets the rate of BH growth in the green valley or the moment of final quenching. The picture is perhaps clearest for EAGLE (Davies et al. 2019; Oppenheimer et al. 2020), in which fractional halo gas mass, $f_{\text{CGM}}$, appears to be a barometer of the degree of quenching. According to this parameter, entry into the green valley corresponds to $E_{\text{BH}}/E_{\text{bind}} \simeq 2$, and complete loss of halo gas corresponds to $E_{\text{BH}}/E_{\text{bind}} \simeq 20$ (Oppenheimer et al. 2020, Figure 2(b)). These numbers are not far from the binding-energy model, which gives 4 and 20, respectively. On that basis, it is plausible that coevolution proceeds in the green valley for EAGLE much as it does in the binding-energy model. A similar case can be made for IllustrisTNG, where $E_{\text{BH}}/E_{\text{bind}}$ again tends to be a fixed factor at full quenching (Terrazas et al. 2020).

Based on the above, we see that real coevolution does occur in the binding-energy model, but only in the green valley, where it plausibly is the agent that produces growth vectors of similar length and slope for all galaxies (see Figure 4). But the ultimate uniformity of quenched galaxies also depends on the uniformity of galaxies as they enter the green valley, which is impressed in our model by the initial scaling law $M_{\text{BH}} = \text{const.} \times \Sigma_l^{1.76}$ for star-forming galaxies. Thus, in the binding-energy model, the final BH regularities of galaxies on the quenched ridgeline arise from the combination of a reciprocal, coevolving phase in the green valley that builds on a unidirectional, causal phase for star-forming galaxies.

This insight motivates a new look, in our model, at the nature of the two $M_{\text{BH}} \sim \Sigma_l$ scaling laws, for star-forming and quenched galaxies separately. In the star-forming phase, $\Sigma_l$ determines $M_{\text{BH}}$ single-handedly according to $M_{\text{BH}} = \text{const.} \times \Sigma_l^{1.76}$, and the scaling relation is an evolutionary track that all galaxies must follow (blue line in Figure 4). Eventually,