\section{Introduction}\label{Sec: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 \citep[see review by][]{Somerville15}, 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$ \citep{Bell04, Faber07, Huang13, Muzzin13, Ilbert13, Tomczak14}, and the vast majority of that growth is in central galaxies, not satellites \citep{Tinker13}. 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.\footnote{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.}

Many studies have modeled galaxy quenching, from semi-analytic models 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 \emph{full panoply} of physical data exhibited by quenching galaxies vs\@.~ 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) 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~\ref{SubSec:Review} and  Figure~\ref{Fig:Mstar_Sig1_SMA1}) 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 black holes 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 \emph{and sloping} implies that at least one more parameter in addition to mass is involved in quenching. We show in this paper that \emph{taking this second parameter to be effective radius and tying black hole 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 black holes 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 black hole masses are computed from their halo properties using simple rules.  As galaxies enter the green valley, we find that black hole growth accelerates, which leads ultimately to full quenching.  The full history of black hole growth from star-forming to quenched is sketched in the ``cartoon'' shown in Figure~\ref{Fig:Cartoon}.

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 \citep[e.g.,][]{Ostriker77,Blumenthal84,Birnboim03,Dekel06,Cattaneo06,Keres09}. The other is ``AGN quenching'' by feedback from growing black holes (BHs), which either ejects gas from galaxies  (``ejective feedback'') \citep[e.g.,][]{Silk98, Springel06, Hopkins06} or heats or disturbs the halo gas in such a way as to prevent it from falling in (``preventive feedback'', ``strangulation'') \citep[e.g.,][]{Croton06,Bower06,Somerville08,Somerville12,Henriques15,Peng15}.\footnote{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 \citep{Somerville15}.}

Both halo quenching and AGN quenching have considerable observational support.  On the halo side, N-body simulations provide reliable models of halo growth \citep[e.g.,][]{Behroozi13}, 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 \citep[e.g.,][]{Arnaud10}. On the AGN side, nearly all quenched galaxies have massive spheroids \citep{Bell08, Bluck14}, and these spheroids harbor massive BHs \citep{Magorrian98, Ferrarese00, Gebhardt00}.  AGN activity peaks just as galaxies are entering quiescence \citep[e.g.,][]{Schawinski07,Schawinski10,Nandra07,Aird19}.  Tight correlations exist between BH mass and host spheroid properties, such as bulge mass ($\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_{\rm bulge}}$: \citealt{Haring04,Kormendy13}) and effective velocity dispersion ($\ensuremath{M_{\rm BH}}\text{-}\ensuremath{\sigma_{\rm e}}$: \citealt{Gebhardt00, Ferrarese00, Kormendy13,Ho14,vandenBosch16}), signaling some sort of close interplay between spheroid-building and BH growth.  Finally, matter falling onto a black hole 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 \citep[e.g.,][]{Silk98,Ostriker05,Oppenheimer18}, enough to either heat, blow away, or otherwise significantly disturb the gas and induce quenching.  For reviews of these topics see \citet{Cattaneo09} and \citet[hereafter KH13]{Kormendy13}.

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 $\ensuremath{M_*} \sim 10^{10}~\ensuremath{\rm 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 $\ensuremath{M_{\rm crit}} \sim 10^{11\text{-}12}~\ensuremath{\rm M_{\odot}}$ \citep[e.g.,][]{Birnboim03,Keres05,Dekel06}, and therefore at a corresponding stellar mass because the stellar-mass/halo-mass relation is tight \citep{Wechsler18}. Halo quenching therefore tends to predict a too-sharp boundary in stellar mass near $\ensuremath{M_*} \sim 3\times10^{10}~\ensuremath{\rm 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'' (\citealp{Somerville18}; also called ``quasar mode'' \citep{Croton06}, ``bright mode'' \citep{Somerville08}, or ``high mode'' \citep{Pillepich18}), AGN mass-accretion is rapid and is often triggered by major mergers \citep{Hopkins06} and/or disk-gas instabilities \citep{Cattaneo06}. Since 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'' (\citealp{Somerville18}; also called ``radio mode'' \citep{Croton06}, ``kinetic mode'' \citep{Weinberger17}, or ``low mode'' \citep{Pillepich18}), 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 \citep[e.g.,][]{Bower06,Croton06, Cattaneo06}, and so in practice, jet-like feedback also needs to wait until the halo mass exceeds $\ensuremath{M_{\rm crit}}$.  Prospects may be better in kinetic mode simulations that are more isotropic \citep[e.g., IllustrisTNG,][]{Pillepich18}, 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, \citet{Peng10} 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 \citep[e.g.,][]{Woo17}, 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 \emph{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 \citep{Fang18}. 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.


\subsection{Previous Work}
\label{SubSec:PreviousWork}

Many studies have remarked on the striking structural parameter correlations of both quenched and quenching galaxies. Let $\ensuremath{R_{\rm e}}$ be the optical effective radius containing half the light and $\ensuremath{\Sigma_{\rm e}}$ be the effective mass surface density given by $\ensuremath{\Sigma_{\rm e}} = 0.5\ensuremath{M_*}/\pi\ensuremath{R_{\rm e}}^2$.  \citet{Kauffmann03} showed that quenched SDSS galaxies have higher $\ensuremath{\Sigma_{\rm e}}$ than star-forming galaxies at the same mass, and \citet{Franx08} showed the same trend in distant quenched galaxies.  The corresponding offset ridgeline of quenched galaxies to smaller $\ensuremath{R_{\rm e}}$ at fixed $\ensuremath{M_*}$ was shown for SDSS galaxies by \citet{Shen03}, \citet{Omand14}, and \citet{Haines17} and for distant galaxies by  \citet[hereafter vdW09]{vanderWel09}, \citet[hereafter vdW14]{vanderWel14}, and \citet{Haines17}.

\citet{Kauffmann03,Kauffmann06,Kauffmann12} further noted that the boundary between quenched and star-forming galaxies corresponds closely to a line of constant effective surface density: $\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^{0.5}$.  Separately, \citet{Bell08}, \citet{Kauffmann12}, and \citet{Bluck14} showed that local quenched galaxies have high S\'{e}rsic~ indices characteristic of bulges, and the same was shown for distant galaxies by \citet{Bell12} and \citet{Wuyts12}. 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.

\citet{Cheung12} introduced the use of \emph{central} properties to predict quenching. After testing numerous quantities, they concluded that high stellar-mass surface density within a projected radius of 1~kpc ($\Sig1$) is the best predictor of quenching at $z \sim 0.7$. \citet[hereafter F13]{Fang13} found a similar correlation between $\Sig1$ and stellar mass for SDSS galaxies, with quenched and green valley galaxies occupying a narrow sloping ridgeline at the top of the $\Sig1\text{-}\ensuremath{M_*}$ distribution.  With an rms scatter of only 0.16~dex \citepalias{Fang13}, the $\Sig1\text{-}\ensuremath{M_*}$ ridgeline is among the tightest scaling relations known for quenched galaxies.  \citetalias{Fang13} further showed that specific star-formation rate (SSFR) declines sharply near this ridgeline and suggested that the residual $\Delta\Sig1$ could be a \emph{marker} for quenching, possibly signaling a threshold in BH feedback. \citet{Wake12} (also \citealp{Bluck16,Teimoorinia16}) found a similar result for SDSS galaxies using central velocity dispersion $\sigma_1$ within 1~kpc (which \citetalias{Fang13} showed correlates closely with $\Sig1$). Finally, \citet{Luo19} showed that high $\Delta\Sig1$ correlates with high bulge-to-total ratio ($B/T$), consistent with the known quenching-spheroid connection.

\citet{vanDokkum14F} demonstrated that the threshold relationship between $\Sig1$ and SSFR found by \citetalias{Fang13} obtains throughout cosmic time for high-mass galaxies.  Lower-mass galaxies were added by \citet{Tacchella15}, \citet{Barro17}, \citet{Whitaker17}, \citet{Lee18}, and \citet{Mosleh18}, who found that the  $\Sig1\text{-}\ensuremath{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 \emph{in unison} by $\sim0.3$~dex from $z = 3$ to now.\footnote{We note parenthetically that, while $\Sig1$ and $\ensuremath{\Sigma_{\rm e}}$ are both surface densities, their behaviors vs\@.~ mass and time are different \citep[hereafter B17]{Barro17}. The outer envelopes of quenched galaxies grow due to dry mergers after quenching, which increase $\ensuremath{R_{\rm e}}$ and decrease $\ensuremath{\Sigma_{\rm e}}$ while leaving $\Sig1$ roughly unchanged \citep[e.g.,][]{Daddi05,Trujillo06,vanDokkum08,Bezanson09,Naab09,vanDokkum10,Oser12,Szomoru12}.  The quenched-galaxy relation is therefore broader in $\ensuremath{\Sigma_{\rm e}}\text{-}\ensuremath{M_*}$ than in $\Sig1\text{-}\ensuremath{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 $\ensuremath{\Sigma_{\rm e}}\text{-}\ensuremath{M_*}$ as in $\Sig1\text{-}\ensuremath{M_*}$ \citepalias{Barro17}.  $\Sig1$ and $\ensuremath{\Sigma_{\rm e}}$ are also different for star-forming galaxies, and the difference makes $\Sig1$ more suitable for predicting black-hole mass.  We show below that $\Sig1 \sim  \ensuremath{M_*}/\ensuremath{R_{\rm e}}$ (Section~\ref{SubSec:InputPowerLaws}), whereas $\ensuremath{\Sigma_{\rm e}} \sim \ensuremath{M_*}/\ensuremath{R_{\rm e}}^2$.  Because BH mass scales approximately as $\ensuremath{\sigma_{\rm e}}^4$  (e.g., \citetalias{Kormendy13}, \citealt{vandenBosch16}) and $\ensuremath{\sigma_{\rm e}}^4 \sim (\ensuremath{M_*}/\ensuremath{R_{\rm e}})^2 \sim \Sig1^2$, a close power-law relation is predicted between $\ensuremath{M_{\rm BH}}$ and $\Sig1$, as borne out below. In contrast, no simple power-law relation is predicted between $\ensuremath{M_{\rm BH}}$ and $\ensuremath{\Sigma_{\rm e}}$.}

Various quenching models have tried to reproduce the structural differences between quenched and star-forming galaxies.  \citet[hereafter vdW09]{vanderWel09F} introduced a model in which galaxies evolve in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and quench when they reach a boundary set by their effective velocity dispersion, $\ensuremath{\sigma_{\rm e}}$.  Since quenched galaxies with higher $\ensuremath{\sigma_{\rm e}}$ today are older at fixed $\ensuremath{M_*}$ \citep{Bernardi05, Graves09}, \citetalias{vanderWel09} reasoned that the threshold $\ensuremath{\sigma_{\rm e}}$ for quenching must have been higher in the past. These facts motivated them to adopt a \emph{moving boundary} that shifts upwards with time in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ (and downwards in $\ensuremath{\sigma_{\rm e}}\text{-}\ensuremath{M_*}$). This  moving boundary was actually detected by \citet{Haines17} in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$.  However, neither of these works specified how galaxies approach the boundary before quenching.  This was added by \citet[hereafter vD15]{vanDokkum15}, who envisioned that galaxies \emph{evolve on parallel tracks} through $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ while star-forming and then quench when meeting a moving boundary in $\ensuremath{\sigma_{\rm e}}$. A similar picture for the motion of galaxies in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and the origin of quenched galaxies was presented by \citet{Cappellari13,Cappellari16}.

\citetalias{Barro17} and \citet{Tacchella15} considered the evolution of galaxies in $\Sig1\text{-}\ensuremath{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 $\Sig1$ for their mass, quench, and then pile up along the $\Sig1\text{-}\ensuremath{M_*}$ ridgeline.  Because the ridgeline in $\Sig1\text{-}\ensuremath{M_*}$ moves down with time, this quenching boundary also moves down.  Thus, like \citetalias{vanderWel09} and \citetalias{vanDokkum15}, \citetalias{Barro17} envisioned a moving boundary but modeled it as a power law in $\Sig1\text{-}\ensuremath{M_*}$, not as a line at constant $\ensuremath{\sigma_{\rm e}}$.  None of these papers offered any physical explanation for their boundaries -- they were simply set empirically to match the data.

\citet{Omand14} introduced yet another aspect of quenching by focusing on $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$. They observed that star-forming SDSS galaxies populate a shallow ridgeline in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ with rms scatter $\sigma_{\rm Re} \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 (\citealp{Omand14}; see also \citealp{Shen03}, \citetalias{vanderWel09}, \citetalias{vanderWel14}, \citealp{Lange16}). Contours of equal quenched fraction are sloping stripes running parallel to the quenched ridgeline, and quenched galaxies therefore \emph{co-exist with star-forming galaxies at the same $\ensuremath{M_*}$}. This fact had been noticed before \citep[e.g.,][]{Kauffmann03}, but Omand et al.~seized on it as evidence that quenching is controlled by at least one more parameter than mass and suggested that galaxies are a \emph{two-parameter family} whose quenching is governed by $\ensuremath{R_{\rm e}}$ as well as $\ensuremath{M_*}$.

Further evidence that both $\ensuremath{R_{\rm e}}$ and $\ensuremath{M_*}$ play a role came from the work of \citet{vandenBosch16} on local black holes.  Starting with the $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{\sigma_{\rm e}}$ relation, he noted that $\ensuremath{\sigma_{\rm e}}^2 \sim \ensuremath{M_*}/\ensuremath{R_{\rm e}}$ from the virial theorem and proposed that BHs populate a ``fundamental plane'' in  $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$.  In words, $\ensuremath{R_{\rm e}}$ can be thought of as cleaning up the scatter in $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ that appears when star-forming galaxies are added to quenched galaxies \citepalias{Kormendy13} in the sense that \emph{smaller-radii galaxies must have larger black holes}.  Similar arguments were made by \citet{Krajnovic18}. Both works stressed that $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ is a map of black hole mass, an insight that will prove basic to our model.

Finally, \citet[hereafter L16]{Lilly16} presented yet another quenching model. Galaxies evolve along parallel log-log tracks in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, as in \citetalias{vanderWel14} and \citetalias{vanDokkum15}, but quenching depends on $\ensuremath{M_*}$ probabilistically, increasing with $\ensuremath{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, \citetalias{Lilly16} showed that the contours of equal quenched fraction would actually be tilted despite the fact that quenching formally depends on $\ensuremath{M_*}$ only.  This model is discussed further in Section~\ref{SubSec:L16}.

\subsection{Collected Properties of Quenched vs\@.~ Star-forming Galaxies}
\label{SubSec:Review}

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~\ref{Tb_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 \emph{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 \citep[e.g.,][]{Ilbert10,Huang13,Muzzin13}\footnote{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.}.  This tail is generally due to satellite quenching, not AGN quenching \citep{Tinker13,Fossati17,Guo17} 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 $\ensuremath{M_*}/\ensuremath{\rm M_{\odot}} \simeq 10.7$ from $z = 3$ to now\footnote{\citet{Bundy06} and \citet{Haines17} 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., \citetalias{Barro17}, \citet{Straatman16}, \citet{Haines17}, \citet{Pan16}, and \citet{Fang12}.} and its substantial width, approaching 1 full dex in stellar mass. Line 3 of Table~\ref{Tb_1} also notes the continuing steady flow of galaxies though the high-mass channel down to recent times (e.g., \citealt{Straatman15}, F18).

\begin{table*}[!htbp]
\caption{\label{Tb_1}Observations To Be Explained$^{\rm a}$}
\renewcommand\arraystretch{1.5}
{\centering
\begin{tabularx}{\textwidth}{cXl}
\toprule
Item No.$^{\rm b}$ & Observation & References\\
\midrule
1 & The stellar mass distribution in the high-mass green valley peaks at approximately the same value  & \citet{Ilbert10}\\
& $\sim 10^{10.7}~\ensuremath{\rm M_{\odot}}$ from $z=3$ to now. & \citet{Huang13}\\
& & \citet{Barro17}\\
2 & The high-mass green valley is $\sim 1$~dex wide in stellar mass from $z = 3$ to now. & \citet{Ilbert10}\\
& & \citet{Huang13}\\
& & \citet{Barro17}\\
3 & There is a steady flow of high-mass central galaxies across the quenching boundary from $z = 3$ to now. & \citet{Bell04}\\
& &\citet{Faber07}\\
& &\citet{Tinker13}\\
& & \citet{Tomczak14}\\
4 & The ridgeline of star-forming galaxies in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ has constant slope $\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^{0.21}$ from $z = 3$ to now. & \citet{vanderWel14}\\
5 & The zero point of the star-forming ridgeline in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ has increased by $\sim 0.4$~dex from $z = 3$ to now. & \citet{vanderWel14}\\
$[6]$ & The ridgeline of quenched galaxies in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ has constant slope $\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^{0.75}$ from $z = 3$ to now. & \citet{vanderWel14}\\
$[7]$ & The zero point of the quenched ridgeline in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ has increased by $\sim 0.75$ dex from $z = 3$ to now. & \citet{vanderWel14}\\
8 & The ridgeline of star-forming galaxies in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ has constant slope $\Sig1 \sim \ensuremath{M_*}^{0.88}$ from  $z = 2.5$ to 0.5. & \citet{Barro17}\\
9 & The zero point of the star-forming ridgeline in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ has decreased by $\sim 0.30$~dex from $z = 2.5$ to 0.5. & \citet{Barro17}\\
$[10]$ & The ridgeline of quenched galaxies in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ has constant slope $\Sig1 \sim \ensuremath{M_*}^{0.66}$ from $z = 2.5$ to 0.5. & \citet{Barro17}\\
$[11]$ & The zero point of the quenched ridgeline in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ has decreased by $\sim 0.30$~dex from $z = 2.5$ to 0.5. & \citet{Barro17}\\
$[12]$ & The zero points of the quenched and AGN galaxy ridgelines in $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ appear to be decreasing with time & \citet{Bennert11}\\
& by similar amounts.  & \citet{Ding17}\\
13 & $\ensuremath{M_{\rm BH}}$ can be expressed as a power law of either $\ensuremath{M_*}/\ensuremath{R_{\rm e}}$ or $\ensuremath{\sigma_{\rm e}}$, but the slopes of quenched galaxies are  & \citet{Kormendy13}\\
&  shallower than the slopes of all galaxies taken together. & \citet{vandenBosch16}\\
14 & The ridgelines of quenched galaxies are narrow in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ but not so narrow in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$. & \citet{Fang13}\\
& & \citet{Barro17}\\
& & \citet{Kormendy13}\\
& & \citet{vanderWel14}\\
15 & The quenched ridgelines in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ show age stratification, with older galaxies  & \citet{Tacchella17}\\
& having higher values of $\Sig1$ and $\ensuremath{M_{\rm BH}}$ at fixed $\ensuremath{M_*}$. \\
& & \citet{Luo19}\\
$[16]$ & The zero points of quenched and star-forming galaxies have evolved together in lockstep in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ since $z = 3$. & \citet{Barro17}\\
$[17]$ & The zero points of quenched and star-forming galaxies have evolved differently in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ since $z = 3$. & \citet{vanderWel14}\\
18 & Star-forming galaxies mix with quenched galaxies on the quenched ridgeline in $\Sig1$ vs\@.~ $\ensuremath{M_*}$ today but lie    & \citet{Fang13}\\
& well below the quenched ridgeline in $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$. & \citet{Barro17}\\
& & \citet{Terrazas17}\\
\bottomrule
\end{tabularx}}

\footnotesize{$^{\rm a}$ Referring to central galaxies only.\\}
\footnotesize{$^{\rm b}$ 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~\ref{Fig:Mstar_Sig1_SMA1}), for which there is only one adjustable parameter, $f_{\rm 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~\ref{SubSec:Successes}.}
\end{table*}

Lines 4-7 in Table~\ref{Tb_1} describe the properties of quenched and star-forming ridgelines in the $\ensuremath{R_{\rm e}}\text{-}\ensuremath{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 $\ensuremath{R_{\rm e}}\text{-}\ensuremath{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 $\Sig1\text{-}\ensuremath{M_*}$.  Important facts are again the constancy of both ridgeline slopes over time, accompanied (in this case) by \emph{decreasing} zero points.

Line~12 of Table~\ref{Tb_1} says that the ridgeline zero point of quenched galaxies in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ also decreases with time. Since 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 $\Sig1\text{-}\ensuremath{M_*}$. The distribution in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ of star-forming AGN galaxies has been studied locally by \citet{Ho14} and \citet{Reines15} and in distant galaxies by \citet{Ding17,Ding20}.  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 $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ distribution for all galaxies maintains the same shape and moves down in bulk, just like that in $\Sig1\text{-}\ensuremath{M_*}$.  Line~13 points out that there seem to exist various scaling laws linking $\ensuremath{M_{\rm BH}}$ to both $\ensuremath{M_*}/\ensuremath{R_{\rm e}}$ and $\ensuremath{\sigma_{\rm e}}$ and that the slopes of these laws differ depending on whether all galaxies or only quenched ridgeline galaxies are included (\citetalias{Kormendy13};\citealp{vandenBosch16}).

Line~14 refers to the width of the quenched ridgelines, which are quite narrow in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ at fixed $\ensuremath{M_*}$ but not so narrow in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, 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 pile-ups 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 $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ is somewhat broader than the other two due to variable amounts of galaxy-wide dissipation during quenching \citep[e.g.,][]{Zolotov15,Tacchella16a} and to variable amounts of dry merging after quenching (\citealp[e.g.,][]{Naab09}; \citetalias{vanderWel09}), both of which tend to decouple global $\ensuremath{R_{\rm e}}$ from central $\Sig1$ and $\ensuremath{M_{\rm BH}}$.  The narrowness of the quenched ridgelines in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ therefore suggests that those relations more closely preserve the memory of quenching; significantly, both are properties that relate to the \emph{centers} of galaxies, not to their outer parts.  However, line~15 points out that recent data are beginning to reveal \emph{age stratification} even within both of these ridgelines, the quenched stellar populations above the ridgelines being older (\citealp{Tacchella17,Martin18,Luo19}).  The authors attribute this to ``progenitor bias'' \citep{vanDokkum96,Carollo13, Poggianti13} 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 vs\@.~ star-forming galaxies in the different spaces.  In $\Sig1\text{-}\ensuremath{M_*}$, ridgelines of both populations move down together in lock-step so that the relative positions of quenched vs\@.~ star-forming galaxies are preserved and the $\Sig1\text{-}\ensuremath{M_*}$ diagram \emph{looks the same at all times} except for this bulk downward shift (cf.~Figure~2 in \citetalias{Barro17}).  However, in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, the ridgeline of quenched galaxies moves up \emph{faster} than the ridgeline of star-forming galaxies and the two populations begin to overlap, as in \citetalias{vanderWel14}.

\begin{figure*}[htbp]
\centering
\includegraphics[scale=0.53]{Fig1_a.png}
\includegraphics[scale=0.53]{Fig1_b.png}
\caption{\label{Fig:Mstar_Sig1_SMA1} Data from all five CANDELS fields divided by redshift, patterned after similar figures in \citetalias{vanderWel09} and \citetalias{Barro17}.  The data and sample are described in Section~\ref{Sec:Data}; only face-on galaxies with $b/a \geq 0.5$ are used.  Satellites and centrals are included even though the model technically refers to centrals only.  Panel a) (top row):  Effective radius $\ensuremath{R_{\rm e}}$ along the semi-major axis (SMA) vs\@.~ stellar mass, colored by $\Sig1$ for star-forming galaxies.  Panel b) (second row):  The same as panel a) but for quenched galaxies. Panel c) (third row): $\Sig1$ vs\@.~ stellar mass colored by $\ensuremath{R_{\rm 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~\ref{Sec:PowerLawModel}. The blue lines are predicted by the halo gas binding-energy model in Section~\ref{Sec:ThermalModel}, which says that entry to the green valley occurs when total emitted BH energy ($\ensuremath{E_{\rm BH}}$) equals $\times4$ the binding energy of hot halo gas. The agreement between the red and blue lines is the major result of this paper.}
\end{figure*}

The first column in Table~\ref{Tb_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 $\ensuremath{R_{\rm e}}$ for star-forming galaxies before they reach the quenching boundary and that $\ensuremath{R_{\rm e}}$ is not modeled later, when dissipation and dry merging become important. Since lines~6, 7, and 17 all involve quenched galaxies, their values of $\ensuremath{R_{\rm e}}$ are not accurately known, which is why they are given brackets.  Other brackets are described in Section~\ref{SubSec:Successes}.

Figure~\ref{Fig:Mstar_Sig1_SMA1} illustrates many of the properties listed in Table~\ref{Tb_1}. This figure borrows heavily from similar figures in \citetalias{vanderWel14} and \citetalias{Barro17} using CANDELS data. The upper two panels compare star-forming and quenched galaxies in the space of $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$, and the lower two panels compare the same samples in the space of $\Sig1$ vs\@.~ $\ensuremath{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~\ref{Sec:PowerLawModel}, 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 $\Sig1$, we see that virtually all quenched galaxies lie on a ridgeline above the red quenching boundary and their numbers are growing with time, while the vast majority of star-forming galaxies lie below the same line.  Galaxies approach the quenching boundary from below-left in this plane, as shown in Figure~\ref{Fig:ThreePanel} below.  A smattering of star-forming galaxies also lies above the line, which is the so-called ``elbow effect'' pointed out by \citetalias{Fang13} and \citetalias{Barro17} in which star-forming and quenched galaxies co-exist at nearly the same values of $\Sig1$ (cf.~also \citealp{Whitaker17,Lee18}).  This is mentioned in line~18 of Table~\ref{Tb_1} and signifies a highly non-linear relation between $\Sig1$ and the decline of star-formation in the GV.  Line~18 also says that the downward trend in star-formation rate vs\@.~ $\ensuremath{M_{\rm BH}}$ is more gradual (straighter) \citep[cf.~Figure 3,][]{Terrazas17} so that star-forming galaxies \emph{never} lie on the quenched ridgeline in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ but always lie well below it \citep[Figure 2,][]{Terrazas16}.  These different behaviors of $\Sig1$ vs\@.~ $\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ near quenching are illustrated in Figures~\ref{Fig:ThreePanel} and \ref{Fig:Elbow} below, and matching them will prove to be a major clue to black-hole growth rates in the green valley.  Note that the quenched ridgelines, and hence also the red lines, maintain constant slope with time in $\Sig1\text{-}\ensuremath{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~\ref{Tb_1}.

The upper two panels in Figure~\ref{Fig:Mstar_Sig1_SMA1} show the same objects in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, where the pattern is reversed -- virtually all quenched galaxies now fall \emph{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~\ref{Fig:ThreePanel}.  Note that the quenched ridgelines, and hence also the red lines, again maintain constant slope but now move \emph{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 lock-step the way they are in $\Sig1\text{-}\ensuremath{M_*}$. Collectively these are properties 4-7, 14, and 17 in Table~\ref{Tb_1}.

Three more features of Figure~\ref{Fig:Mstar_Sig1_SMA1} are important.  First, one sees that quenched galaxies cover roughly a full dex in mass that is centered at about $10^{10.7}~\ensuremath{\rm M_{\odot}}$ at all redshifts.  These are lines~1 and 2 in Table~\ref{Tb_1}.  Second, there is a dearth of quenched galaxies near and below $10^{10}~\ensuremath{\rm 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~\ref{Fig:Mstar_Sig1_SMA1} 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}~\ensuremath{\rm 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~\ref{Tb_1} and illustrated in Figure~\ref{Fig:Mstar_Sig1_SMA1}.


\subsection{Outline for the Remainder of the Paper}
\label{SubSec:Outline}

The first step in this paper, in Section~\ref{Sec:PowerLawModel}, is an empirical parametrization 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 \emph{total energy quota} needed by the halo to stop cooling: $\ensuremath{E_{\rm BH}} = \ensuremath{E_{\rm quench}}$, and we derive an empirical power-law expression for $\ensuremath{E_{\rm quench}}$ as a function of $\ensuremath{M_{\rm vir}}$ and redshift from the observed (red) boundaries in Figure~\ref{Fig:Mstar_Sig1_SMA1}. Starting with this empirical formulation for $\ensuremath{E_{\rm 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~\ref{Sec:ThermalModel}, compares this empirical halo-energy quota to the gas binding energies of real halos.  This comparison is shown by the blue lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}, 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_{\rm 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 \citep[e.g.,][]{Silk98,Ostriker05,Booth10,Booth11,Bower17,Oppenheimer18,Davies19a,Davies19b,Oppenheimer19}.  Our major new result is that this simple model explains the observed quenching boundaries in Figure~\ref{Fig:Mstar_Sig1_SMA1} quite well, along with many other data in Table~\ref{Tb_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 \emph{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 $\Sig1$ during the star-forming phase, which is suggested by empirical data. This rule implies that black holes grow modestly but steadily (averaging over possible short-term high-accretion episodes) while galaxies form stars.\footnote{Black-hole growth as a power of $\Sig1$ during the star-forming phase is the default model for most of the paper.  However, the Discussion section considers a complementary picture in which black holes remain smaller than this throughout most of the star-forming phase and only begin to grow rapidly later (see Section~\ref{SubSec:BHGrowth}).}  The third feature replaces this early, modest power law with a steeper power law (vs\@.~ $\Sig1$) in the green valley (see Figure~\ref{Fig:Cartoon}).  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 \citet{Omand14,vandenBosch16,Krajnovic18}, is to model star-forming galaxies as a \emph{two-parameter family} labeled by $\ensuremath{R_{\rm e}}$ and $\ensuremath{M_*}$. This says that BH mass depends on galaxy radius as well as $\ensuremath{M_*}$ -- larger galaxies have smaller BHs because they have lower central $\Sig1$-- 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~\ref{Fig:Mstar_Sig1_SMA1} and \ref{Fig:ThreePanel}.

Section~\ref{Sec:Variants} explores four variants of the basic model. All fail, but each one teaches an important lesson, as discussed in Appendix~\ref{App:A4}.  Section~\ref{Sec:Discussion} is the Discussion.  Topics include the model's success in matching the data goals in Table~\ref{Tb_1}, the origin of scatter in effective radii, how quenching happens in the binding-energy model compared to other models, compaction, the connection between black holes and dark halos, black hole scaling laws, and what it means to say that black holes  and galaxies ``co-evolve''.  Finally, Section~\ref{Sec:Summary} reviews the basic findings and important take-away 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 \citep[cf.][]{Chauke19}.  Because of its use of $\Sig1$ to track black hole growth, the model cannot be used before $z \sim 3$, where many galaxies fit entirely within 1~kpc.  Likewise, no description 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; predict mass functions or other structure functions; or predict AGN frequencies. These topics are left for future work.

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

\section{Data}\label{Sec:Data}

The data used in this paper come from the multi-wavelength and ancillary datasets produced by the Cosmic Assembly Near-Infrared Deep Extragalactic Legacy Survey on the \emph{Hubble Space Telescope} \citep[CANDELS; ][]{Grogin11,Koekemoer11}. All five fields are used (GOODS-S, UDS, COSMOS, EGS, and GOODS-N) with $H \leq 24.5$~mag.  AB magnitudes, spectroscopic/photometric redshifts, stellar masses, star-formation rates, and $A_{\rm V}$ attenuation values are taken from the following catalogs: \citet{Guo13,Galametz13,Nayyeri17,Stefanon17,Barro19}. Methodologies for the stellar masses and $A_{\rm V}$ values are described in \citet{Santini15}. $\ensuremath{R_{\rm e}}$ values and axis ratios are measured using GALFIT by \citet{vanderWel12}. The selection and completeness of the sample are discussed in Appendix~\ref{App:A1}.

An important parameter in this paper is central stellar density, $\Sig1$, whose measurement
is explained in \citetalias{Barro17}.  Briefly, mass profiles were constructed from multiband \emph{HST} images, and raw values of \Sig1\ 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 \ensuremath{R_{\rm e}}\ and S\'{e}rsic~ index.   Details are given in Appendix C of \citetalias{Barro17}.
Two separate sets of values were computed using slightly different techniques by Barro and Chen, and they agree well. All data including \Sig1\ are available at the \emph{Rainbow} database \citep{Barro11}\footnote{http://arcoiris.ucsc.edu//Rainbow\_navigator\_public/}.  The values of \Sig1 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.

\section{Power-Laws in the Empirical Model}\label{Sec:PowerLawModel}

\subsection{Input Power Laws}\label{SubSec:InputPowerLaws}

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~\ref{App:A3}.

Figure~\ref{Fig:ThreePanel} 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 comprised of $\ensuremath{M_*}$, $\ensuremath{R_{\rm e}}$, $\Sig1$, and $\ensuremath{M_{\rm BH}}$.  The three panels in Figure~\ref{Fig:ThreePanel} show the three 2-D projections of this four-dimensional space onto $\ensuremath{M_*}$.  Eventually galaxies encounter a boundary in this 4-D space, as shown by the dashed lines, where they start to quench.  Crossing the boundary is equated with \emph{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 $\Sig1\text{-}\ensuremath{M_*}$ distribution of SDSS galaxies from \citet{Luo19} 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 black holes from \citet{Terrazas16}.

\begin{figure*}[htbp]
\centering
\includegraphics[scale=0.56]{Fig2.png}
\caption{\label{Fig:ThreePanel} 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 comprised of $\ensuremath{M_*}$, $\ensuremath{R_{\rm e}}$, $\Sig1$, and $\ensuremath{M_{\rm BH}}$. The panels show projections onto $\ensuremath{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 \ensuremath{M_*}\ is involved in quenching. Panel a: $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$.  The points are face-on CANDELS galaxies with $b/a > 0.5$ in the redshift range $z = 0.6\text{-}0.9$. Blue points are star-forming galaxies on the SFMS with $\Delta$logSSFR above $-0.45$~dex; red points are fully quenched galaxies with  $\Delta$logSSFR below $-1.0$~dex.  Green valley galaxies are not shown.  Panel b: $\Sig1$ vs\@.~ $\ensuremath{M_*}$.  The objects are the same as in panel a.  Panel c: $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$.  This diagram reproduces Figure~2 from \citet{Terrazas16}, which plots black-hole masses for 90 nearby galaxies.  Blue points are star-forming galaxies with the same selection criteria as in \citet{Terrazas16}; 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 \emph{move with time}.}
\end{figure*}

\begin{figure*}[htbp]
\centering
\includegraphics[scale=0.53]{Fig3.png}
\caption{\label{Fig:Mstar_Sig1_SMA2} $\ensuremath{R_{\rm e}}$ and $\Sig1$ vs\@.~ $\ensuremath{M_*}$ from all five CANDELS fields for star-forming galaxies only.  This figure differs from Figure~\ref{Fig:Mstar_Sig1_SMA1} in showing only star-forming galaxies and to a lower mass limit.  Galaxies satisfy $H \leq 24.5$ mag, axis ratio $b/a \geq 0.5$, and are classed as star-forming using the same criterion as in Figure~\ref{Fig:Mstar_Sig1_SMA1}.  Upper panel: Effective radius $\ensuremath{R_{\rm e}}$ along the semi-major axis (SMA) vs\@.~ stellar mass, colored by $\Sig1$.  Lower panel: $\Sig1$ vs\@.~ stellar mass colored by $\ensuremath{R_{\rm e}}$.  Loci of constant $\Sig1$ (above) and constant $\ensuremath{R_{\rm e}}$ (below) are visible as colored stripes.  These indicate the presence of a plane in $\ensuremath{R_{\rm e}}\text{-}\Sig1\text{-}\ensuremath{M_*}$ space for star-forming galaxies, and the two plots are different projections of it. $\Sig1 \sim \ensuremath{M_*}^q \ensuremath{R_{\rm e}}^{-p} \sim \ensuremath{M_*}^{1.1}\ensuremath{R_{\rm e}}^{-1}$  is the best mapping of this plane from one projection to another (Appendix~\ref{App:A1}).  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 (\citealp{Barro13,Barro14,Nelson14,Williams14}; \citetalias{vanDokkum15}).}
\end{figure*}

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~\ref{Sec:Variants} 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~\ref{Sec:Variants}).

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

$\bullet$ A \emph{structural mapping} relates one set of structural variables to another.  The main example is the mapping $\Sig1 \sim \ensuremath{M_*}^q\ensuremath{R_{\rm e}}^{-p}$, which maps the $\ensuremath{R_{\rm e}} \text{-} \ensuremath{M_*}$ plane onto the $\Sig1 \text{-} \ensuremath{M_*}$ plane.  The variables $p$ and $q$ depend on S\'{e}rsic~ $n$ \citepalias{Barro17}, and hence a given $p$,$q$ mapping applies only to galaxies with similar S\'{e}rsic~ indices.  The needed mapping is for star-forming galaxies, which generally have $n$ in the range 1-2 \citepalias{Barro17}.

$\bullet$ A \emph{quenched ridgeline} is the fit at a given redshift to the locus of galaxies that have quenched.  It is narrowest and best defined in $\Sig1\text{-}\ensuremath{M_*}$, where it has the slope 0.66 at all redshifts \citepalias{Barro17}.  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.

$\bullet$ A \emph{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 $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ plane (e.g., \citetalias{vanderWel14}), but it exists also in the $\Sig1\text{-}\ensuremath{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.

$\bullet$ A \emph{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 $\ensuremath{M_{\rm BH}} \sim \Sig1^v$ in star-forming galaxies and that total effective energy radiated during the course of BH mass growth varies as $\ensuremath{E_{\rm BH}} \sim \eta \ensuremath{M_{\rm BH}} c^2$.

$\bullet$ An \emph{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 $\ensuremath{M_*} \sim \ensuremath{M_{\rm vir}}^s$. Since 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.

$\bullet$ A \emph{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 $\Sig1\text{-}\ensuremath{M_*}$, $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$, which together form a single line in the master 4-D 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:

\begin{enumerate}
\item \label{item_BasicPowerLaw:1} \emph{SMHM relation for star-forming galaxies:} This evolutionary track is parameterized as
    \begin{equation}\label{Eq:MstarMvir}
    \ensuremath{M_*} = {\rm const.}\times\ensuremath{M_{\rm vir}}^s,
    \end{equation}
where $\ensuremath{M_*}$ is total stellar mass and $\ensuremath{M_{\rm vir}}$ is halo virial mass.  $\ensuremath{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~\ref{SubSec:EvolTracks}.  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 \citet[hereafter RP17]{Rodriguez17} from abundance matching (see Figure~\ref{Fig:HMSM}).  Other works give slightly different results, and Variant~2 in Section~\ref{Sec:Variants} 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 \citetalias{Rodriguez17}.  We note that this is very close to the value $s \sim 5/3$ predicted from the energetics of SNae feedback vs\@.~ halo potential well depth \citep{Dekel86,Dekel03}.  The SMHM law below the knee applies to star-forming galaxies only, and since no scatter is assumed, it is simply referred to as the SMHM relation.  Above the knee, the SMHM or its inverse, the HMSM, is used as appropriate.\footnote{The knee occurs at approximately \(\ensuremath{M_*} \sim 10^{10.7} ~\ensuremath{\rm M_{\odot}}\), which is not quite massive enough to encompass all star-forming galaxies.  This problem is addressed in Variant~1 in Section~\ref{Sec:Variants}.}  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.

\item \label{item_BasicPowerLaw:2} \emph{Median galaxy effective radius, $\ensuremath{\langle R_{\rm e}\rangle}$, vs\@.~ halo virial radius, $\ensuremath{R_{\rm vir}}$, for star-forming galaxies:} This evolutionary track assumes that the median $\ensuremath{\langle R_{\rm e}\rangle}$ of star-forming galaxies at a given stellar mass is the same fixed fraction of the halo virial radius, $\ensuremath{R_{\rm vir}}$, at all masses and redshifts ($\ensuremath{R_{\rm e}}$ is the observed half-light semi-axis major from \citet{vanderWel12}, called SMA in Figures~\ref{Fig:Mstar_Sig1_SMA1}-\ref{Fig:Mstar_Sig1_SMA2}).  Thus:
    \begin{equation}\label{Eq:ReRvir}
    \ensuremath{\langle R_{\rm e}\rangle} = {\rm const.}\times\ensuremath{R_{\rm vir}}.
    \end{equation}
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\% (\citealp{Kravtsov13,Shibuya15,Huang17}; \citetalias{Rodriguez17}; \citealp{Somerville18}).\footnote{All observed radii in this paper are projected 2-D semi-major-axis half-light radii.  These can differ significantly and systematically from 3-D half-mass stellar radii, as shown in \citet{Somerville18}.  Eq.~\ref{Eq:ReRvir} is obeyed more accurately by 3-D half-mass radii, but the differences are unimportant given the approximate nature of our treatment.}
\end{enumerate}

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

\begin{enumerate} [resume]
\item \label{item_BasicPowerLaw:3} \emph{$\Sig1$ vs\@.~ $\ensuremath{R_{\rm e}}$ and $\ensuremath{M_*}$ for star-forming galaxies:} This relation is an empirical structural mapping of the form
    \begin{equation}\label{Eq:Sig1MstarRe}
    \Sig1 = {\rm const.}\times\ensuremath{M_*}^q\ensuremath{R_{\rm e}}^{-p}.
    \end{equation}
    where $\Sig1$ is the projected stellar mass within 1~kpc of the center of the galaxy.  Figure~\ref{Fig:Mstar_Sig1_SMA2}a plots $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ for CANDELS star-forming galaxies color-coded by $\Sig1$, and Figure~\ref{Fig:Mstar_Sig1_SMA2}b plots $\Sig1$ vs\@.~ $\ensuremath{M_*}$ color-coded by $\ensuremath{R_{\rm e}}$.  The slanting colored stripes in both panels indicate the existence of a plane in $\Sig1\text{-}\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ space (a similar plane was shown for VIPERS galaxies at $z = 0.5$-1.0 by \citet{Haines17}). Fitted values for the exponents of this plane are $p = 1.00\pm0.01$, $q = 1.10\pm0.01$, which are essentially constant with redshift, as is the zero point.  Comparisons of data to predictions from this relation are shown in Figure~\ref{Fig:A1} in Appendix~\ref{App:A1}, where rms residuals in $\ensuremath{R_{\rm e}}$ and $\Sig1$ are both $\sim 0.2$~dex.  This $p$, $q$ mapping applies to star-forming galaxies only and is valid because their S\'{e}rsic~ indices are all in the range $n = 1\text{-}2$ \citepalias{Barro17}.  (The analogous mapping for quenched galaxies has $p = 0.50$, $q = 1.10$ and is less tight owing to the more variable S\'{e}rsic~ indices of quenched galaxies.)

\item \label{item_BasicPowerLaw:4} \emph{$\Sig1$ vs\@.~ central velocity dispersion, $\sigma_1$:} \citetalias{Fang13} empirically demonstrated the scaling law that
    \begin{equation}\label{Eq:Sig1sig1}
    \Sig1 = {\rm const.}\times\sigma_1^{2.0\pm0.2},
    \end{equation}
    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 $\ensuremath{M_*} \sim R \sigma^2/G$ applied to the mass inside a sphere if we note that
    \begin{equation*}
    \ensuremath{M_*} \sim \Sig1\times\pi R^2
    \end{equation*}
    and evaluate all quantities at $R = 1$~kpc.  It then follows that  $\Sig1 \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.  Eq.~\ref{Eq:Sig1sig1} is assumed to be valid for star-forming and quenched galaxies with the same constant at all redshifts.  The main use of Eq.~\ref{Eq:Sig1sig1} is not for actual computations but to help motivate what comes next.
\end{enumerate}

\begin{figure*}[htbp]
\centering
\includegraphics[scale=1.0]{Fig4.png}
\caption{\label{Fig:Cartoon} Cartoon showing scaling laws for $\ensuremath{M_{\rm BH}}$ vs\@.~ $\Sig1$.  In the default model, star-forming galaxies enter along the blue track at lower left, turn off into the green valley at some point along the steep green tracks, and come to rest at the top of the tracks along the quenched ridgeline (red dashed line).  The quenching boundaries shown as the red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1} represent entry into the green valley and correspond to the bends in the tracks.  The three green vectors correspond to galaxies with the same starting values of $\ensuremath{M_*}$ at a given epoch but different values of $\ensuremath{R_{\rm e}}/\ensuremath{\langle R_{\rm e}\rangle}$: small-radii galaxies turn off sooner and produce quenched galaxies with smaller $\ensuremath{M_*}$, $\ensuremath{M_{\rm BH}}$, $\ensuremath{R_{\rm e}}$, and $\Sig1$, etc.  The effect of the quenching boundary moving in time is illustrated by the green parallelogram, which shows the location of the green valley at an earlier time, when quenching values of $\Sig1$ and $\ensuremath{M_{\rm BH}}$ were larger; this parallelogram slides down along the scaling laws with time.  Except for this, no other aspects of the figure change with time.  Section~\ref{SubSec:BoundaryCrossing} cites evidence that the green-valley portions of the tracks (the green vectors) have the same slope and length independent of mass and time.  If so, the star-forming and quenched galaxy ridgelines will be parallel, and the observed slope of the upper relation can be used to determine the slope of the lower relation.  This is the assumption employed in Item~\ref{item_BasicPowerLaw:5} (this section).  The derivations of the slopes and zero points of the quenched and star-forming relations are fully described in Appendix~\ref{App:A3}.  Our fitted values of $p$, $q$ say that $\Sig1 \approx \ensuremath{M_*}/\ensuremath{R_{\rm e}}$ for star-forming galaxies, which is also roughly $\ensuremath{\sigma_{\rm e}}^2$ from \citet{vandenBosch16}, and thus also roughly $\sigma_1^2$ if $\sigma_1 \approx \ensuremath{\sigma_{\rm e}}$.  We have therefore used these approximate relations to label the X-axis as the log of all of these quantities (ignoring the fact that $\Sig1 \approx \ensuremath{M_*}/\ensuremath{R_{\rm e}}$ is not quite true for quenched galaxies).  The black dashed line shows schematically what happens if one attempts to fit the BH masses of all galaxies with a single power law, as in \citet{vandenBosch16}.  The slope of the black dashed line (in the same units) is $\sim 2.6$.  The value of $\Sig1$ for M31 is based on the enclosed mass at 1~kpc in the M31 mass model by \citet{Blanadiaz18}, and the value of $\Sig1$ for the Milky Way is based on a rotation speed of $175~\rm{km~sec^{-1}}$ at 1~kpc \citep{Courteau14,Blandhawthorn16}.  The dashed green lines below the blue track show alternative growth histories for BHs in star-forming galaxies, as discussed in Section~\ref{SubSec:BHGrowth}.}
\end{figure*}

The last step is to specify rules for BH growth.  The basic assumption is that BH mass scales as a power law of  $\Sig1$, but with potentially different slopes and zero points for star-forming, quenched, and green valley galaxies. The adopted laws are shown in Figure~\ref{Fig:Cartoon}.  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 $\Sig1$ for different-$\ensuremath{R_{\rm e}}$ galaxies at different times (see below), which creates multiple tracks through the GV 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.


\begin{figure*}[htbp]
\centering
\includegraphics[scale=0.53]{XXXX.png}
\caption{\label{Fig:Elbow} Cartoon illustrating the different evolution of $\Sig1$ and $\ensuremath{M_{\rm BH}}$ through the green valley.  Panel a: $\Delta \log$\Sig1\ is the residual of $\log \Sig1$ with respect to the quenched ridgeline in $\Sig1\text{-}\ensuremath{M_*}$, and $\Delta \log SSFR$ is the residual of $\log SSFR$ with respect to the star-forming main sequence.  $\Sig1$ shows a characteristic ``elbow'' shape whereby $SSFR$ falls sharply but $\Sig1$ grows only slightly \citepalias{Fang13, Barro17} through the green valley.  Panel b: The analogous plot for $\Delta \log \ensuremath{M_{\rm BH}}$ based on data from  \citet[][Figure 3]{Terrazas17}.  In contrast to $\Sig1$ in panel a, black hole mass grows considerably in the green valley.  The delta-quantities above the figures are used to plot the green-valley vectors in Figure~\ref{Fig:Cartoon}.}
\end{figure*}


The series of broken power-law tracks in Figure~\ref{Fig:Cartoon} is motivated by the following logic:  1) Black holes in local quenched galaxies obey a power-law relation $\ensuremath{M_{\rm BH}} \sim \ensuremath{\sigma_{\rm e}}^{u}$ \citepalias{Kormendy13}.  We assume $\sigma_1 \sim \ensuremath{\sigma_{\rm e}}$, in which case $\ensuremath{M_{\rm BH}} \sim \Sig1^{v}$, where $v = u/2.0$ from Item~\ref{item_BasicPowerLaw:4} above.  The value of $v$ is the slope of the upper, quenched ridgeline in Figure~\ref{Fig:Cartoon}.  2) A steep track in the GV is then indicated by the fact that star-forming and quenched galaxies can have similar $\Sig1$ values today at the same stellar mass but that $\ensuremath{M_{\rm BH}}$ differs strongly for these same two populations (line~18 of Table~\ref{Tb_1}).  This phenomenon is visible by comparing panels b and c of Figure~\ref{Fig:ThreePanel} and is illustrated more graphically in Figure~\ref{Fig:Elbow}.  Together these imply that $\Sig1$ largely stops increasing in the GV (hence the ``elbow") whereas $\ensuremath{M_{\rm BH}}$ keeps growing. 3) Based on a feature seen in SDSS galaxies \citep[Figure 5b,][]{Luo19}, the entry to the GV (i.e., the quenching boundary) is taken to be parallel to the quenched ridgeline in $\Sig1\text{-}\ensuremath{M_*}$ but is displaced downward from it by a fixed 0.2~dex at all stellar masses.  This motivates the increase of 0.2~dex that has been applied to $\Sig1$ across the GV for all galaxies in Figure~\ref{Fig:Cartoon} (green vectors; see also Figure~\ref{Fig:Elbow}a).  4) Likewise, the separation between the quenched and star-forming ridgeline in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ is observed to be $\sim 1.5$~dex at all stellar masses \citep[Figure 3]{Terrazas17}, which we take to be the amount of BH mass growth across the GV (see illustration in Figure~\ref{Fig:Elbow}b).  5) Assumptions 3) and 4) applied to all galaxies then imply that the lengths and slopes of the green GV tracks are independent of galaxy mass today and thus are all parallel.  6) Finally, point 5) means that the slope of the lower track for star-forming galaxies is the same as the slope of the upper scaling law for quenched galaxies, and hence can be determined from it.\footnote{Three further aspects of the cartoon model in Figure~\ref{Fig:Cartoon} merit comment: 1) The evolution in BH mass is steep vs\@.~ $\Sig1$, but this does not mean that it is necessarily steep in time.  The total time needed to cross the GV is not specified and could be long.  2) Similarly, since AGN luminosities depend on BH mass accretion \emph{rates}, they also cannot be calculated without additional assumptions about evolutionary times.  3) The fractional BH mass growth in the GV is taken to be the offset in BH mass between the star-forming and quenched populations in \citet[Figure 3]{Terrazas17}, which are separated by about $\sim$1.5 dex.  However, this snapshot of mass differences now does not equal the BH growth factor \emph{for an individual galaxy} if boundaries and ridgelines have been moving over time.  The sense is that 1.5~dex would be an overestimate of the growth factor, but the correction cannot be more than a few tenths of a dex (Section~\ref{SubSec:BoundaryCrossing}) and has not been applied in Figure~\ref{Fig:Cartoon}. It is discussed further in Section~\ref{Sec:Discussion}.}

The main justification for the above framework is the empirical observation that $\ensuremath{M_{\rm BH}} \sim \sigma_e^{4\text{-}5}$ in many types of galaxies (\citetalias{Kormendy13}; \citealt{Ho14}), which translates to a power law in $\Sig1$ from Eq.~\ref{Eq:Sig1sig1}.  Many papers have tried to derive the $\sigma^4$ law analytically \citep[e.g.,][]{Silk98,Ostriker00,Ciotti01,Adams01,Colgate03,King03,Cen07,Cen15,King15}, but there is as yet no generally accepted explanation.  As time passes, evidence is actually growing for two \emph{different} laws with different zero points, one for quenched galaxies \citepalias{Kormendy13} and one for star-forming galaxies \citep{Ho14,Reines15,Savorgnan16}, 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~\ref{Fig:Cartoon}.  More discussion of the adopted laws and their derivation is given in Section~\ref{Sec:Discussion}.

The preceding paragraphs have motivated the law that $\ensuremath{M_{\rm BH}} = \rm{const.} \times \Sig1^v$ for star-forming galaxies and that BHs grow more steeply relative to $\Sig1$ in the GV.  We now extrapolate those rules and apply them to star-forming and GV galaxies at all redshifts. The main argument in favor of this assumption is the similar appearance of the $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, $\Sig1\text{-}\ensuremath{M_*}$, and  $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{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~\ref{Fig:Cartoon} interprets the blue BH scaling law for star-forming galaxies \emph{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~\ref{SubSec:BHGrowth}.

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

\begin{enumerate} [resume]
\item \label{item_BasicPowerLaw:5} \emph{BH mass vs\@.~ $\Sig1$ for star-forming galaxies:}
    \begin{equation}\label{Eq:MbhSig1}
    \ensuremath{M_{\rm BH}} = {\rm const.}\times\Sig1^v,
    \end{equation}
    where $u$ is the slope of the quenched ridgeline, and $v = u/2.0$.  This is the blue line in Figure~\ref{Fig:Cartoon}, which is assumed to hold for all star-forming galaxies at all redshifts $z<3$.  Since $\ensuremath{M_{\rm BH}}$ is uniquely specified by $\Sig1$, this means that $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ maps onto $\Sig1\text{-}\ensuremath{M_*}$ and that $\ensuremath{M_{\rm BH}}$ is now also joined to the space of  $\Sig1\text{-}\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ for star-forming galaxies.
\end{enumerate}

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

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

To model the scatter, we multiply each starting radius by a scatter factor $R_{\rm e,s}/\langle R_{\rm e}\rangle$, where $\langle R_{\rm e}\rangle$ is the mean value of $\ensuremath{R_{\rm e}}$ for $10^9~\ensuremath{\rm M_{\odot}}$ galaxies at $z_{\rm s}$.  This scatter factor is assumed to be constant throughout the star-forming phase, and its logarithm is assumed to be Gaussianly distributed with $\sigma_{\ensuremath{R_{\rm e}}} = 0.2$~dex.   Every galaxy is now labeled by two numbers, the redshift $z_{\rm s}$ at which $\ensuremath{M_*} = 10^9~\ensuremath{\rm M_{\odot}}$ and its value $\ensuremath{\Delta} \log R_{\rm e,s}$ at that redshift.  The latter is also equivalent to the scatter in $\ensuremath{M_{\rm BH}}$ and $\Sig1$ since the previous equations imply:
\begin{equation}\label{Eq:Sig1Re}
\Delta \log \Sigma_{\rm 1,s} = -p\Delta \log R_{\rm e,s},
\end{equation}
and
\begin{equation}\label{Eq:Sig1Mbh}
\Delta \log M_{\rm BH,s} = v\Delta \log \Sigma_{\rm 1,s},
\end{equation}
and
\begin{equation}\label{Eq:ReMbh}
\Delta \log M_{\rm BH,s} = -vp\Delta \log R_{\rm e,s}.
\end{equation}
Note that galaxies with larger $\ensuremath{R_{\rm e}}$ at fixed $\ensuremath{M_*}$ have smaller BHs.

The next step is to evolve galaxies forward in time.  To do this, another crucial assumption is made, namely, that \emph{the star-formation rate at a given mass is independent of $\ensuremath{R_{\rm e}}$}.  This fact is now established from several studies on both the SDSS \citep{Omand14,Lin19} and CANDELS/3D-HST (\citealp{Whitaker17}; \citetalias{Fang18}; \citealp{Lin19}).  Galaxies with the same starting $\ensuremath{M_*}$ but different starting radii therefore grow in $\ensuremath{M_*}$ at the same rate, and each starting mass cohort stays together.  Since all relationships are power laws, it follows that $\ensuremath{\Delta} \log R_{\rm e,s}$, $\ensuremath{\Delta} \log \Sigma_{\rm 1,s}$, and $\ensuremath{\Delta} \log M_{\rm BH,s}$ are constant through time. Hence, galaxies at a given $\ensuremath{M_*}$ that start out with large radii always have smaller-than-average black holes, and vice versa.  If quenching requires reaching some threshold BH mass that does not vary too rapidly with $\ensuremath{M_*}$, large galaxies will tend to do this later and attain a larger stellar mass before quenching. The fact that starting radius confers an $\ensuremath{M_{\rm BH}}$ enhancement to some galaxies and a penalty to others is what creates the \emph{sloped quenching boundaries} in $\ensuremath{R_{\rm e}}$\text{-}$\ensuremath{M_*}$, $\Sig1$\text{-}$\ensuremath{M_*}$, and $\ensuremath{M_{\rm BH}}$\text{-}$\ensuremath{M_*}$ in Figure~\ref{Fig:ThreePanel}.

The final step is to specify the equation of the quenching boundary.  For this, our first-stage, empirical power-law model assumes that the \emph{total effective energy} emitted by a BH during its lifetime is equal to some \emph{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 \emph{halo quenching-energy quota}, and in the empirical model we parameterize it as a simple power of halo mass.  Thus
\begin{equation}
\ensuremath{E_{\rm BH}} = \ensuremath{E_{\rm quench}},
\end{equation}
or
\begin{equation}\label{Eq:MbhMvir}
\eta \ensuremath{M_{\rm BH}} c^2 = {\rm const.}\times\ensuremath{M_{\rm vir}}^t,
\end{equation}
where $\eta$ is an efficiency factor for converting BH accretion energy into useful feedback energy.  We reiterate that crossing over the boundary corresponds to \emph{entry} into the GV, and thus that BHs are just starting to affect their halos and galaxies at that point.  The use of $\ensuremath{E_{\rm 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 $\ensuremath{M_{\rm 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 $\Sig1$ 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 \emph{slopes} of the quenched ridgelines in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$, and the constant \emph{must vary with time} in order to match the moving boundaries in those same spaces.


\subsection{Quenching Boundary Power Laws}\label{SubSec:QuenchSlope}

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
\begin{equation*}
\ensuremath{M_{\rm BH}} \sim \ensuremath{M_{\rm vir}}^t.
\end{equation*}
by Eq.~\ref{Eq:MbhMvir}.  But
\begin{equation*}
\ensuremath{M_{\rm BH}} \sim \Sig1^v,
\end{equation*}
from Eq.~\ref{Eq:MbhSig1}, and
\begin{equation*}
\ensuremath{M_{\rm vir}} \sim \ensuremath{M_*}^{1/s}.
\end{equation*}
from Eq.~\ref{Eq:MstarMvir}.  Substituting, we have
\begin{equation}\label{Eq:Sig1MstarSlope}
\Sig1 \sim \ensuremath{M_*}^{t/sv},
\end{equation}
which gives the slope, $t/sv$, of the instantaneous boundary in $\Sig1\text{-}\ensuremath{M_*}$.  Now apply the mapping
\begin{equation*}
\Sig1 \sim \ensuremath{M_*}^q\ensuremath{R_{\rm e}}^{-p}
\end{equation*}
from Eq.~\ref{Eq:Sig1MstarRe} to obtain
\begin{equation}\label{Eq:ReffMstarSlope}
\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^{q/p - t/svp},
\end{equation}
which gives the slope, $q/p - t/svp$, of the instantaneous boundary in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$.
Finally, $\ensuremath{M_{\rm BH}} \sim \Sig1^v$ gives
\begin{equation}\label{Eq:MbhMstarSlope}
\ensuremath{M_{\rm BH}} \sim \ensuremath{M_*}^{t/s},
\end{equation}
which is the slope, $t/s$, of the instantaneous boundary in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{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:

\begin{enumerate}
\item \label{item_QuenchSlope:1} \emph{Slope of the quenching boundary in $\Sig1$ vs\@.~ $\ensuremath{M_*}$:} Even though the spaces $\Sig1\text{-}\ensuremath{M_*}$, $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ are all closely related (for star-forming galaxies), the quenched ridgeline is narrowest and best determined in $\Sig1\text{-}\ensuremath{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 $\Sig1\text{-}\ensuremath{M_*}$ distribution found by \citealp{Luo19}) and therefore has the same slope as the ridgeline.  \citetalias{Barro17} found the slope of this ridgeline to be nearly constant at 0.66 from $z = 3$ to now.  Eq.~\ref{Eq:Sig1MstarSlope} now gives: \footnote{\citetalias{Barro17} 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~\ref{Fig:Mstar_Sig1_SMA1}.  However, the ridgeline slope is insensitive to this definition since green valley galaxies and fully quenched galaxies occupy nearly the same regions in $\Sig1\text{-}\ensuremath{M_*}$ \citepalias{Fang13,Barro17}.}
    \begin{equation}\label{Eq:Sig1Mstar}
    \frac{\Delta \log \Sig1}{\Delta \log \ensuremath{M_*}} = \frac{t}{vs} = 0.66\pm0.02,
    \end{equation}
where the error is based on the scatter in four independent measurements at four different redshifts in \citetalias{Barro17}.  This is the first constraint.\footnote{For simplicity, the shape of the quenched ridgeline in $\Sig1\text{-}\ensuremath{M_*}$ is assumed here to be a pure power law, as found by \citetalias{Barro17}.  Recent work by \citet{Luo19} for SDSS galaxies indicates some downward curvature in $\Sig1\text{-}\ensuremath{M_*}$ at high masses, which may actually improve quenching predictions for high-mass galaxies, as discussed under Variant~1 in Section~\ref{SubSec:BrokenPLModel}.}

As discussed below, we do not use the quenched ridgeline slope in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ as an observational constraint because of the sensitivity of effective radii to dry merging and other post-quenching processes.  However, the slope of the \emph{quenching boundary} in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ will be needed in future.  Eq.~\ref{Eq:ReffMstarSlope} now gives:
    \begin{eqnarray}\label{Eq:ReffMstar}
    \frac{\Delta \log \ensuremath{R_{\rm e}}}{\Delta \log \ensuremath{M_*}} & = & \frac{q}{p}  - \frac{t}{vps} \nonumber\\
                                         & = & 1.10 - \frac{t}{vs}  = 0.44\pm0.02,
    \end{eqnarray}
    where the second line follows from inserting $t/vs = 0.66\pm0.02$ and $p = 1.0\pm0.01$, $q = 1.1\pm0.01$ from Appendix~\ref{App:A1}.

    Note that, although the quenching boundary and quenched ridgeline are assumed to be parallel in $\Sig1\text{-}\ensuremath{M_*}$, they are not parallel when mapped into $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ because quenched galaxies have different S\'{e}rsic~ indices and $\ensuremath{R_{\rm e}}$ due to dry merging.  This is a major reason why the quenched ridgeline in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ has the steep slope of 0.75 \citepalias{vanderWel14} whereas the quenching boundary slope in Eq.~\ref{Eq:ReffMstar} is only 0.44.

  \item \label{item_QuenchSlope:2} \emph{Slope of the quenching boundary in $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$:} As shown in Figure~\ref{Fig:ThreePanel}c here and more clearly in Figure~3 of \citet{Terrazas17}, the separation between fully quenched and star-forming galaxies in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ is a fixed $\sim 1.5$~dex today for galaxies above $10^{10}~\ensuremath{\rm 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 $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$, as it is in $\Sig1\text{-}\ensuremath{M_*}$, and the quenched slope in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ can therefore be used as another observational constraint on the unknown exponents $s$, $t$, and $v$. The local observed $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ relation consists of a narrow ridgeline of fully quenched galaxies with star-forming galaxies scattering below it (\citetalias{Kormendy13}; \citealp{vandenBosch16, Terrazas16, Krajnovic18}).  We equate the observed quenched ridgeline slope to that of local ellipticals and bulges, which is 1.16 from \citetalias{Kormendy13}.  Thus:
    \begin{equation}\label{Eq:MbhMstar}
    \frac{\Delta \log \ensuremath{M_{\rm BH}}}{\Delta \log \ensuremath{M_*}} = \frac{t}{s} = 1.16\pm0.08,
    \end{equation}
    where the error comes from \citetalias{Kormendy13}.  Since 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.

\item \label{item_QuenchSlope:3} Finally, the slope of the SMHM relation, $s$, is observed to be $1.75\pm0.17$ based on the abundance matching results of \citetalias{Rodriguez17} and \citet{Behroozi19}.  This is the third constraint.
\end{enumerate}

Constraints are now available from Items~\ref{item_QuenchSlope:1}, \ref{item_QuenchSlope:2}, and \ref{item_QuenchSlope:3}, and the solution is uniquely determined.  The combination $s = 1.75\pm0.17$, $t/s = 1.16\pm0.08$, and $t/vs = 0.66\pm0.02$ implies $t = 1.16 \times 1.75 = 2.03\pm0.24$, and $v = 1.16/0.66 = 1.76\pm0.13$, and from this last we derive $u = 2.0v = 3.52\pm0.26$, as was employed in Eq.~\ref{Eq:MbhSig1} 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~\ref{Fig:Cartoon} is actually correct and whether the same $\ensuremath{M_{\rm BH}}\text{-}\Sig1$ scaling relations obtain at all redshifts.  The general constancy of the cartoon relations is supported by the fixed power-law appearance of the $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ distributions back in time (Lines 12 and 16 in Table~\ref{Tb_1}), and by correct predictions for both evolving zero points when combined with the moving boundary hypothesis, as discussed in Section~\ref{SubSec:Successes}.  In short, the constant slopes and offsets of quenched and star-forming galaxies in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{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 \emph{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~\ref{Fig:Cartoon}, but the quenched slope may have evolved after quenching due to mergers.  For example, to match the scaling laws of quenched galaxies today requires $\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^2$ along the merging trajectory (\citealp[e.g.,][]{Naab09}; \citetalias{vanderWel09}), and thus $\ensuremath{\sigma_{\rm e}} \sim \ensuremath{M_*}^{-0.5}$, i.e., $\ensuremath{\sigma_{\rm e}}$ actually declines as galaxies merge!  Both of these vectors (in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ and $\ensuremath{\sigma_{\rm e}}$ vs\@.~ $\ensuremath{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 $\ensuremath{\sigma_{\rm e}}$ still more \citep{Boylan06, Lauer07}.  Some evidence of steepening in $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{\sigma_{\rm e}}$ at high masses has been noted (\citealt{Lauer07, Gultekin09, Krajnovic18}; \citetalias{Kormendy13}).  To summarize, effective quantities like $\ensuremath{R_{\rm e}}$ and $\ensuremath{\sigma_{\rm 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 $\Sig1\text{-}\ensuremath{M_*}$, since centers are more immune to merging than envelopes and the slope of $\Sig1\text{-}\ensuremath{M_*}$ is actually measured to be constant over time \citepalias{Barro17}.  In $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$, if pre-merger BHs merge along with their host galaxies and if BH mass is conserved, galaxies would tend to move \emph{along} this ridgeline since 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 $\ensuremath{M_{\rm BH}} \sim \sigma_1^{3.52\pm0.26}$. This value stems from $v = 1.76$, since $u$ is $2.0 \times v$.  But $v$ is simply the ratio of the slopes of the quenched ridgelines in $\Sig1\text{-}\ensuremath{M_*}$ and  $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$, both of which are well determined.  Alternative values for $u$ are $4.40\pm0.30$ for quenched ridgeline galaxies from \citetalias{Kormendy13} and $5.35\pm0.23$ for all galaxies together from \citet{vandenBosch16}.  At first glance, these seem incompatible with our value of 3.52 and also with each other. (These laws use $\ensuremath{\sigma_{\rm e}}$ while we use $\sigma_1$, but the difference is probably not large.)  However, as described in Section~\ref{SubSec:PredictionsTests} below, the Figure~\ref{Fig:Cartoon} 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: $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_{\rm vir}}^t$ (the halo energy quota needed to halt gas cooling), $\ensuremath{M_*} \sim \ensuremath{M_{\rm vir}}^s$ (the SMHM relation), and $\ensuremath{M_{\rm BH}} \sim \Sig1^v$ (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 $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$, which are observed, and $s$ is separately measured from abundance matching.  This makes three constraints, which together determine the three exponents.

\section{The Evolving Empirical Power-law Model}\label{Sec:EvolPLModel}

\subsection{Model for Galaxy Radii; Evolutionary Tracks}\label{SubSec:EvolTracks}

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  \emph{evolving} power-law model, so-called because the halo energy quota, $\ensuremath{E_{\rm quench}}$, is represented as a power law in $\ensuremath{M_{\rm vir}}$ with exponent $t$ set to 2.0, and a moving boundary is introduced to match the moving (red) quenching boundaries in Figure~\ref{Fig:Mstar_Sig1_SMA1}.  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 $\ensuremath{R_{\rm vir}}$ vs\@.~ $\ensuremath{M_{\rm vir}}$ to derive evolutionary tracks for \emph{central star-forming galaxies} in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, $\Sig1\text{-}\ensuremath{M_*}$, and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$.  Ingredients and assumptions in the model are listed in Table~\ref{Tb_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 $\ensuremath{M_*}$ and $\ensuremath{R_{\rm 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 \citet{Rodriguez15} and \citetalias{Rodriguez17}.  Figure~\ref{Fig:Trajectoryz} plots the resulting stellar and virial parameters for three representative different-mass halos and their associated galaxies vs\@.~ time.  The adopted constant value for the ratio $\ensuremath{\langle R_{\rm e}\rangle/R_{\rm vir}} = 0.02$ comes from the abundance matching results in \citetalias{Rodriguez17} for star-forming galaxies below $\ensuremath{M_*} = 10^{10.5}~\ensuremath{\rm M_{\odot}}$ (see Item~\ref{item_BasicPowerLaw:2} in Section~\ref{SubSec:InputPowerLaws} and Figure~\ref{Fig:HMSM} in Appendix~\ref{App:A4}).  See also the discussion of radii evolution in \citet{Somerville18}, which justifies a constant value.

\begin{table*}[!htbp]
\renewcommand\arraystretch{1.5}
{\centering
\caption{\label{Tb_2}Model Ingredients and Assumptions}
\begin{tabularx}{\textwidth}{cXl}
\toprule
No. & Item & References\\
\midrule
1 & Halo masses grow according to mean trajectories versus time derived from N-body simulations. & \citet{Rodriguez17}\\
2 & Central galaxies evolve along a power-law SMHM relation while star-forming, with constant slope $s = 1.75$, constant zeropoint, and no scatter. & \citet{Rodriguez17}\\
3 & The median effective radii $\ensuremath{\langle R_{\rm e}\rangle}$ of star-forming galaxies are $0.02 \times \ensuremath{R_{\rm vir}}$ at all masses and redshifts. & \citet{Huang17}\\
& & \citet{Somerville18}\\
4 & There is significant real scatter in both $\ensuremath{R_{\rm e}}$ and $\Sig1$ among star-forming galaxies at fixed $\ensuremath{M_*}$ and $\ensuremath{R_{\rm vir}}$. & \citet{Omand14}\\
& & \citet{vanderWel14}\\
& & \citet{Fang13}\\
& & \citet{Barro17}\\
5 & A galaxy's offset in radius $\ensuremath{\Delta}\log \ensuremath{R_{\rm e}}$ from the median $\ensuremath{\langle R_{\rm e}\rangle}$ is constant during the star-forming phase. SF galaxies evolve along parallel tracks in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$. & Assumed, Section~\ref{SubSec:EvolTracks}\\
6 & For star-forming galaxies, $\Sig1 = {\rm const.}~\ensuremath{M_*}^q \times \ensuremath{R_{\rm e}}^{-p}$, where $q = 1.1$ and $p = 1.0$. & Fitted to CANDELS data in Appendix~\ref{App:A1}. \\
7 & The star-formation rate of galaxies on the SF main sequence is independent of $\ensuremath{R_{\rm e}}$ and $\Sig1$ at fixed $\ensuremath{M_*}$. & \citet{Omand14}\\
& & \citet{Whitaker17}\\
& & \citet{Fang18}\\
& & \citet{Lin19}\\
8 & $\Delta \log \Sig1 = 0.2$~dex and $\Delta \log \ensuremath{M_{\rm BH}} = 1.5$~dex across the green valley, independent of stellar mass and   & \citet{Luo19}\\
& redshift. & \citet{Terrazas17}\\
& & Assumed in Section~\ref{SubSec:InputPowerLaws}\\
& & and discussed in Section~\ref{SubSec:BoundaryCrossing}\\
9 & Black-hole masses in star-forming galaxies grow as $\ensuremath{M_{\rm BH}} = {\rm const.}~\Sig1^{1.76}$ with the same zero point and slope at all redshifts and no scatter.$^{\rm a}$ & Assumed, Section~\ref{SubSec:InputPowerLaws}\\
10 & The effective energy radiated from the black hole is $\ensuremath{E_{\rm BH}} = \eta \ensuremath{M_{\rm BH}} c^2$, where $\eta = 0.01$. & Assumed, Section~\ref{SubSec:EvoTracksQuenchBound}\\
11 & 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 $\ensuremath{E_{\rm quench}} = k_{\rm e}(z) \left(\frac{\ensuremath{M_{\rm vir}}}{10^{12}~\ensuremath{\rm M_{\odot}}}\right)^t$ with no scatter. & Assumed, Section~\ref{SubSec:EvoTracksQuenchBound}\\
12 & 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 \ensuremath{E_{\rm bind}} = 4 \times 1/2~f_{bary}\ensuremath{M_{\rm vir}}\ensuremath{V_{\rm vir}}^2$ with no scatter. & Assumed, Section~\ref{Sec:ThermalModel} \\
13 & Mergers during the star-forming and GV phases are ignored, and galaxies do not rejuvenate in the GV. & Assumed throughout\\
\bottomrule
\end{tabularx}}

\footnotesize{$^{\rm a}$ This is the default assumption in both the empirical power-law and binding-energy models.  An alternative version in which black holes are smaller than this in star-forming galaxies and \emph{cross over} this scaling law as they enter the green valley is discussed in Section~\ref{SubSec:BHGrowth}.}
\end{table*}

\begin{figure}[htbp]
\centering
\includegraphics[scale=0.44]{Fig5.png}
\caption{\label{Fig:Trajectoryz} Evolutionary tracks for $\ensuremath{M_*}$, $\ensuremath{M_{\rm vir}}$, $\ensuremath{\langle R_{\rm e}\rangle}$, $\ensuremath{R_{\rm vir}}$, $\ensuremath{T_{\rm vir}}$, and $\ensuremath{V_{\rm vir}}$ for three sample star-forming central galaxies vs\@.~ redshift.  $\ensuremath{M_{\rm vir}}$ and $\ensuremath{R_{\rm vir}}$ are based on the evolutionary model for dark halos summarized in \citetalias{Rodriguez17}.  $\ensuremath{M_*}$ is derived from $\ensuremath{M_{\rm vir}}$ using the SMHM relation $\ensuremath{M_*} \sim \ensuremath{M_{\rm vir}}^s$ with $s = 1.75$ from Section~\ref{SubSec:InputPowerLaws}.  The constant value $\ensuremath{\langle R_{\rm e}\rangle}/\ensuremath{R_{\rm vir}} = 0.02$ is taken from \citetalias{Rodriguez17} (see Item~\ref{item_BasicPowerLaw:2} in Section~\ref{SubSec:InputPowerLaws} and Figure~\ref{Fig:HMSM}).  Tracks start at the starting redshift, $z_{\rm s}$, which is the redshift where stellar mass $\ensuremath{M_*} = 10^9~\ensuremath{\rm M_{\odot}}$.  The quantity $\ensuremath{V_{\rm vir}}$ is used in Section~\ref{Sec:ThermalModel}, which discusses the halo binding-energy interpretation of the halo energy quota.  The quantity $\ensuremath{T_{\rm vir}}$ is used to compute heating and cooling rates in Variants~3 and 4 in Appendix~\ref{SubSec:AppRatesorTotals}.}
\end{figure}

The conventions we are using for $\ensuremath{M_{\rm vir}}$ and $\ensuremath{R_{\rm vir}}$ are based on the definition of halo overdensity $\Delta_{\rm vir}$ by \citet{Bryan98} for an isothermal sphere mass profile.  The formula relating virial mass and radius is:
\begin{equation}\label{Eq:RvirMvir}
\ensuremath{M_{\rm vir}} = \frac{4\pi}{3}\Delta_{\rm vir}\rho_{\rm m}\ensuremath{R_{\rm vir}}^3,
\end{equation}
where $\rho_{\rm m}$ is the mean matter density $\Omega_{\rm m}\rho_{\rm c}$ \citep{Rodriguez16_2} and $\Delta_{\rm vir}$ is given by \citet{Bryan98}.  The quantity $\ensuremath{V_{\rm vir}}$ is defined such that $\ensuremath{M_{\rm vir}} \equiv \ensuremath{V_{\rm vir}}^2 \ensuremath{R_{\rm vir}}/G$.

The tracks in Figure~\ref{Fig:Trajectoryz} are used to evolve central star-forming galaxies in the $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ plane through time.  In order to check them and to develop handy power-law approximations, Figure~\ref{Fig:TrajectoryM} replots them vs\@.~ $\ensuremath{M_{\rm vir}}$ and $\ensuremath{M_*}$, where the relations in Eq.~\ref{Eq:Sig1MstarRe} and Eq.~\ref{Eq:MbhSig1} are used to derive $\langle\Sig1\rangle$ and $\langle\ensuremath{M_{\rm BH}}\rangle$ from $\ensuremath{\langle R_{\rm e}\rangle}$ and $\ensuremath{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 in the halo quenching condition  as in Eq.~\ref{Eq:EheatEbh} and is discussed in Section~\ref{Sec:Variants} and Appendix~\ref{App:A4}.

\begin{figure}[htbp]
\centering
\includegraphics[scale=0.44]{Fig6.png}
\caption{\label{Fig:TrajectoryM} The same evolutionary tracks from Figure~\ref{Fig:Trajectoryz} for central star-forming galaxies are replotted vs\@.~ mass, with $\langle\Sig1\rangle$ and $\langle\ensuremath{M_{\rm BH}}\rangle$ added.  The black solid lines are for the empirical power-law model with fixed SMHM zeropoint.  The black dashed lines are analogous tracks for the Variant~2 model, which has an evolving zeropoint in the SMHM relation (Section~\ref{SubSec:EvolvingSMHMZP}).  The corresponding red lines show approximate power-law slopes (with arbitrary zero points) that are valid for galaxies in the range $\ensuremath{M_{\rm vir}} = 10^{11.5\text{-}12.5}~\ensuremath{\rm M_{\odot}}$ and $z = 3$ to 0.5.  The first exponent refers to the empirical power-law model with fixed SMHM zero point; the second exponent refers to the Variant~2 model with the evolving SMHM zero point.  These slopes are useful for rough estimates; the real curves are used when computing model predictions.}
\end{figure}

The predicted tracks for $\ensuremath{\langle R_{\rm e}\rangle}$ in Figure~\ref{Fig:TrajectoryM} show the motions of single galaxies through time, whereas plots of observed $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ are snapshots of different galaxies at the same time.  Such snapshots are seen to populate a broad range in $\ensuremath{R_{\rm e}}$, have a shallower slope than the individual evolutionary tracks, and have a ridgeline zero point that evolves upwards in time \citepalias{vanderWel14}.  By definition, the peak of the ridgeline is the locus of galaxies with $\ensuremath{R_{\rm e}} = \ensuremath{\langle R_{\rm e}\rangle}$.  The properties of a ridgeline can therefore be calculated from the formula for $\ensuremath{\langle R_{\rm e}\rangle}$ versus time from \citetalias{Rodriguez17}.  From Eq.~\ref{Eq:MstarMvir}, Eq.~\ref{Eq:ReRvir}, and Eq.~\ref{Eq:RvirMvir}, we get:
\begin{eqnarray}\label{Eq:ReMstar}
\ensuremath{\langle R_{\rm e}\rangle} & \sim & \ensuremath{R_{\rm vir}} \nonumber \\
          & \sim & \ensuremath{M_{\rm vir}}^{1/3}\Delta_{\rm vir}^{-1/3}\rho_{\rm m}^{-1/3} \nonumber \\
          & \sim & \ensuremath{M_*}^{1/(3s)}\Delta_{\rm vir}^{-1/3}\rho_{\rm m}^{-1/3}.
\end{eqnarray}
Thus the predicted ridgeline slope $\Delta \log \ensuremath{\langle R_{\rm e}\rangle}/\Delta \log \ensuremath{M_*} = 1/(3\times1.75) = 0.19$, and the zero point evolves as $\Delta_{\rm vir}^{-1/3}\rho_{\rm 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 \citet{Bryan98} prescription for $\Delta_{\rm vir}$ is used.  Both trends agree well with observed data as noted by \citetalias{vanderWel14}, 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 \citetalias{vanDokkum15}, who reviewed observational and theoretical data on evolutionary tracks for star-forming galaxies in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$.  They found an average evolutionary track slope $\Delta \log \ensuremath{R_{\rm e}}/\Delta \log \ensuremath{M_*} = 0.3$, in reasonable agreement with our average track slope of 0.40 in Figure~\ref{Fig:TrajectoryM}.

A final consistency check can be made by comparing to slopes in $\Sig1\text{-}\ensuremath{M_*}$ using the $p,q$ mapping.  A given log slope $m$ in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ corresponds to slope $q-pm$  in $\Sig1\text{-}\ensuremath{M_*}$.  For $p = 1$, $q = 1.1$, the predicted ridgeline slope for star-forming galaxies is then $1.1 - 1.0\times0.19 = 0.91$.  This compares well to the fitted ridgeline slope in \citetalias{Barro17}, which averages 0.88 at all redshifts.  The corresponding evolutionary track slope is $1.1 - 1.0\times0.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}~\ensuremath{\rm M_{\odot}}$ at $z = 2.5$ in \citetalias{Barro17}, computing a stellar mass increase to $10^{10.2}~\ensuremath{\rm M_{\odot}}$ at $z = 0.75$ \citepalias{Rodriguez17}, and placing it in the middle of the star-forming cloud in $\Sig1\text{-}\ensuremath{M_*}$ of \citetalias{Barro17} at the new redshift and mass.  The resulting changes are 1.30~dex in $\ensuremath{M_*}$ and 0.90~dex in $\Sig1$, for a net evolutionary track slope of 0.69, in excellent agreement with the predicted slope (inferred from $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$) of 0.70.

Finally, star-forming galaxies are assumed to evolve along the power-law SMHM relation (Eq.~\ref{Eq:MstarMvir}) until they enter the green valley.

In summary, a simple halo-based model for $\ensuremath{M_*}$ and $\ensuremath{R_{\rm 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 $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and $\Sig1\text{-}\ensuremath{M_*}$.

\subsection{Starting Conditions for Each Galaxy}\label{SubSec:GalIniCon}

In order to ensure that galaxies cross the quenching boundary in the right place, star-forming galaxies must be properly positioned in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ to start with. Power-law evolutionary tracks are now in place, but each model galaxy needs to be provided with appropriate starting conditions for entry onto its track. We start by choosing a galaxy stellar mass today and use the observed halo-mass-to-stellar-mass (HMSM) relation at $z = 0$ from \citetalias{Rodriguez17} to calculate a corresponding halo mass.\footnote{The HMSM relation is the inverse of the SMHM relation, being the mean value of $\ensuremath{M_{\rm vir}}$ observed at fixed $\ensuremath{M_*}$ today.  It is taken from \citetalias{Rodriguez17} and is not a power law.}  This halo mass is then propagated backwards in time using the \citetalias{Rodriguez17} halo mass evolution formula together with the power-law SMHM relation, $\ensuremath{M_*} = k_{\rm m} \ensuremath{M_{\rm vir}}^s$, until $\ensuremath{M_*} = 10^9~\ensuremath{\rm M_{\odot}}$, which sets the starting redshift, $z_{\rm s}$.  This halo mass, stellar mass, and virial radius are then evolved forward using the halo evolution model of \citetalias{Rodriguez17} 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_{\rm m}$ is a constant.

It remains to scatter the starting effective radii.  A mean starting $\langle R_{\rm e,s} \rangle$ at $z_{\rm s}$ for galaxies with $\ensuremath{M_*} = 10^9~\ensuremath{\rm M_{\odot}}$ is computed from the $\ensuremath{R_{\rm vir}}$ formula in Eq.~\ref{Eq:ReRvir}, and galaxies are scattered above and below this value with a Gaussian distribution having $\sigma \log (\ensuremath{R_{\rm e}}/\ensuremath{\langle R_{\rm e}\rangle}) = 0.20$~dex, which is the observed rms scatter in the radii of star-forming galaxies (\citetalias{vanderWel14}; \citealp{Somerville18}).  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 $\ensuremath{R_{\rm e}}$ at all redshifts (\citealp{Omand14,Whitaker17}; \citetalias{Fang18}; \citealp{Lin19}) so that all galaxies in a given $\ensuremath{M_*}$ cohort stay together in mass regardless of their starting $\ensuremath{R_{\rm e}}$.


\subsection{Evolutionary Tracks in Relation to the Quenching Boundary}\label{SubSec:EvoTracksQuenchBound}

When the above evolutionary tracks were computed, a surprising fact emerged.  We had anticipated that galaxies would evolve through the $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and $\Sig1\text{-}\ensuremath{M_*}$ planes toward the sloping boundary, where specific star formation rates fall -0.45 dex below the 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 $\Sig1\text{-}\ensuremath{M_*}$ has slope 0.66 \citepalias{Barro17} whereas the estimated evolutionary track has slope 0.70, i.e., nearly the same. Likewise, the quenching boundary in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{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 \emph{moving  boundary}. A downward motion of $-0.3$~dex in the $\Sig1\text{-}\ensuremath{M_*}$ quenched ridgeline zero point was measured by \citetalias{Barro17} from $z = 2.5$ to 0.5, and the corresponding upward motion in $\ensuremath{R_{\rm e}}$ is $+0.3$~dex.  Our quenching boundaries are tied to these ridgelines.  This motion, which seems small at first sight, is actually the \emph{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 $\Sig1$ and from below in $\ensuremath{R_{\rm e}}$.

To make further progress, this motion needs to be inserted into the expression for the quenching boundary.  It was shown above in Eq.~\ref{Eq:MbhMvir} that equating the total energy emitted by the BH to a power of halo mass matches the \emph{slope} of the quenching boundary if $t = 2.0$.  We now take that power-law relation and add a moving zero point:
\begin{equation}\label{Eq:EheatEbh}
  \ensuremath{E_{\rm BH}} = \ensuremath{E_{\rm quench}}
\end{equation}
or
\begin{equation}\label{Eq:EbhEhalo}
  \eta \frac{\ensuremath{M_{\rm BH}}}{10^8~\ensuremath{\rm M_{\odot}}} c^2 = k_{\rm e}(z) \left(\frac{\ensuremath{M_{\rm vir}}}{10^{12}~\ensuremath{\rm M_{\odot}}}\right)^t,
\end{equation}
where $k_{\rm e}(z)$ varies with time and we have taken $\eta = 0.01$ \citep[cf.][]{Oppenheimer19}.  The left side is the effective emitted BH energy, $\ensuremath{E_{\rm BH}}$, and the right side is the halo energy quota for quenching, $\ensuremath{E_{\rm quench}}$.  From $\ensuremath{M_{\rm BH}} \sim \Sig1^v$ (Eq.~\ref{Eq:Sig1Mbh}), we find that the change $\Delta \log k_{\rm e}(z)$ with time is equal to $v$ times the change in the zero point of the quiescent $\Sig1$ ridgeline in \citetalias{Barro17} at fixed $\ensuremath{M_*}$. A fit of this zero point vs\@.~ redshift, called $B^Q$ in \citetalias{Barro17}, is shown in Figure~\ref{Fig:A2} in Appendix~\ref{App:A2}.  It yields $\log B^Q = 9.44 + 0.66 \log (H(z)/H(0))$, and thus $k_{\rm e}(z) = k_{\rm e}(0) \times h(z)^{0.66v}$.

The constant $k_{\rm 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 $\Sig1\text{-}\ensuremath{M_*}$ based on the data on SDSS green valley galaxies of \citet{Luo19}, and we apply that same shift at all redshifts.  The value of $k_{\rm e}(0)$ is then found to be $10^{17.66}$ km$^2$ sec$^{-2}$.  Finally, the moving boundary in  $\Sig1\text{-}\ensuremath{M_*}$ is mapped into $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ using the $p,q$ mapping.  This prescription for $k_{\rm e}(z)$ (and $t$) perfectly matches the red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}.

We pause here to make two points.  The first is that the moving zero point in Eq.~\ref{Eq:EbhEhalo} is needed to match the empirically determined quenching boundaries (red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}), which move with time.  However, a more physical model for quenching based on real halo properties is coming up in Section~\ref{Sec:ThermalModel}, which also naturally produces a moving 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 Eq.~\ref{Eq:EheatEbh} 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~\ref{Sec:ThermalModel}.  On the other hand, one could omit the labels and view Eq.~\ref{Eq:EbhEhalo} simply as an empirical relation between $\ensuremath{M_{\rm BH}}$ and a power of $\ensuremath{M_{\rm 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 $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ are illustrated in Figure~\ref{Fig:TrajectoryRM}.  These tracks use the actual trajectories of halos, not the  approximate power-laws.  The black line is for an average galaxy with radius $\ensuremath{R_{\rm e}} = \ensuremath{\langle R_{\rm e}\rangle}$ and $\ensuremath{M_*} = 10^{10.6}~\ensuremath{\rm M_{\odot}}$ today.  The blue tracks above and below this line are galaxies of the same starting mass but scattered up and down in $\ensuremath{R_{\rm e}}$ by $\pm0.75\sigma_{\rm Re}$ and $\pm1.5\sigma_{\rm Re}$, where $\sigma_{\rm Re} = 0.2$~dex (two extra values are added on the low side).  Tracks with lower values of $\ensuremath{R_{\rm e}}$ have higher BH masses and vice versa.  The observed distribution of galaxies at $z = 0.5\text{-}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 \emph{zero point} that causes galaxies to cross it.

\begin{figure}[htbp]
\centering
\includegraphics[scale=0.6]{Fig7.png}
\caption{\label{Fig:TrajectoryRM} Sample evolutionary paths for a selection of galaxies of a given stellar mass in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ in the empirical power-law model.  The black line is for an average galaxy with radius $\ensuremath{R_{\rm e}} = \ensuremath{\langle R_{\rm e}\rangle}$ and $\ensuremath{M_*} = 10^{10.6}~\ensuremath{\rm 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 $\ensuremath{R_{\rm e}}$ by $\pm0.75~\sigma_{\rm Re}$ and $\pm 1.5 \sigma_{\rm Re}$ (two extra values are added on the low side).  Tracks with lower values of $\ensuremath{R_{\rm e}}$ have higher BH masses.  The observed distribution of galaxies is shown for comparison at $z = 0.5\text{-}1.0$, colored by $\Sig1$.  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 $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ after entering the green valley is not covered by the model and is not shown.}
\end{figure}

To summarize, a galaxy is identified as quenching perceptibly (i.e., entering the green valley) when the quenching condition in Eq.~\ref{Eq:EbhEhalo} is satisfied, and the redshift when 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 $\Sig1\text{-}\ensuremath{M_*}$ at that redshift, which is taken to be 0.2~dex below the observed quenched ridgeline.

All the power-law expressions used in the model and their slopes and zero points are collected and explained for reference in Appendix~\ref{App:A3}.


\subsection{Results for the Empirical Power-Law Model}\label{SubSec:EvolModelResults}

\begin{figure*}[htbp]
\centering
\includegraphics[scale=0.66]{Fig8.png}
\caption{\label{Fig:BHEnergy} The BH energy output ($\ensuremath{E_{\rm BH}}$, rainbow bands) and the halo quenching-energy quota ($\ensuremath{E_{\rm 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.  $\ensuremath{E_{\rm BH}}$ is the total effective energy, $\eta\ensuremath{M_{\rm BH}} c^2$, emitted by the black hole to that epoch, and $\ensuremath{E_{\rm quench}}$ is the total energy quota that needs to be absorbed by the halo as galaxies enter the green valley (Eq.~\ref{Eq:EbhEhalo}).  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: since galaxies above the lines are in the green valley, their real values of $\ensuremath{M_{\rm BH}}$ and $\ensuremath{E_{\rm BH}}$ would be different than shown.  Star-forming galaxies do not actually populate these locations; see text.)  The width of the rainbow bands corresponds to scattering the initial values of $R_{\rm e,s}$ by $ \pm 1.5\sigma_{\rm Re} = 0.30$~dex, which scatters $\ensuremath{M_{\rm BH}}$ by $\pm 0.52$~dex from Eq.~\ref{Eq:ReMbh}.   Galaxies with smaller $\ensuremath{R_{\rm e}}$'s (and bigger BHs) at a given $\ensuremath{M_*}$ quench sooner.  Notice the acute ``angle of attack'' between $\ensuremath{E_{\rm BH}}$ and $\ensuremath{E_{\rm 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 $\ensuremath{E_{\rm quench}}$ to match the observed quenching boundary.  The blue solid lines show the halo quenching-energy quota, 4$\times$\ensuremath{E_{\rm bind}}, from the binding-energy model in Section~\ref{Sec:ThermalModel}.  It is seen that $\ensuremath{E_{\rm quench}}$ follows $4 \times \ensuremath{E_{\rm bind}}$ very closely, showing that the empirical power-law energy quota we have used to quench a halo (Eq.~\ref{Eq:EbhEhalo}) 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~\ref{Fig:Mstar_Sig1_SMA1}.}
\end{figure*}

\begin{figure*}[hbtp]
\centering
\includegraphics[scale=0.58]{Fig9.png}
\caption{\label{Fig:SimuCartoon}  Model star-forming galaxies from the empirical power-law model in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$ (panel a) and in $\Sig1$ vs\@.~ $\ensuremath{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~\ref{Fig:Mstar_Sig1_SMA1} (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 Eq.~\ref{Eq:EbhEhalo} to evolve appropriately.  The model galaxy points have been scattered by a total of $\pm 1.5\sigma_{\rm Re}$ in $R_{\rm e,s}$, and the shaded regions outline the corresponding $\pm 1.5\sigma_{\rm Re}$ areas that are covered by star-forming galaxies in \citetalias{vanderWel14}.  The similarity of the two distributions indicates that our method for choosing starting values of $R_{\rm e,s}$ (Section~\ref{SubSec:GalIniCon}) 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 $\Sig1$ (or up in $\ensuremath{R_{\rm 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 \ref{Sec:Variants}.}
\end{figure*}

Results for the empirical  power-law model are shown in Figures~\ref{Fig:BHEnergy} and \ref{Fig:SimuCartoon}.  Figure~\ref{Fig:BHEnergy} plots the BH heating energy, $\ensuremath{E_{\rm BH}}$, and the halo heating-energy quota, $\ensuremath{E_{\rm quench}}$, vs\@.~ redshift for three different halo-galaxy combinations.  The curves are labeled with ending stellar and halo mass today.  The red dashed lines show $\ensuremath{E_{\rm quench}}$ from  Eq.~\ref{Eq:EbhEhalo}, which depends only on halo mass and is therefore shown as a single line for each $\ensuremath{M_*}$.  The rainbow-colored bands show $\ensuremath{E_{\rm BH}}$; their spread reflects the range of BH masses coming from a $\pm 1.5\sigma$ scatter in $\ensuremath{R_{\rm e}}$ (smaller $\ensuremath{R_{\rm e}}$ and larger BHs are at the upper edges).  A given galaxy quenches when its $\ensuremath{E_{\rm BH}}$ line crosses above $\ensuremath{E_{\rm quench}}$.  We stress that $\ensuremath{E_{\rm 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 $\ensuremath{E_{\rm BH}}$ trajectories have been plotted ignoring the fact that galaxies perceptibly quench when they cross over the $\ensuremath{E_{\rm quench}}$ lines and their subsequent $\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}$ evolution switches to the green valley vectors in Figure~\ref{Fig:Cartoon}.  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 $\ensuremath{E_{\rm 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 declining value of $k_{\rm e}(z)$ in the $\ensuremath{E_{\rm quench}}$ zero point (Eq.~\ref{Eq:EbhEhalo}).  Without this bend, $\ensuremath{E_{\rm quench}}$ would continue to rise with $\ensuremath{M_{\rm vir}}$, and galaxies would quench very late or not at all. Note that small-$\ensuremath{R_{\rm e}}$ galaxies with bigger BHs quench sooner (upper envelopes).  Finally one sees that a critical parameter is the ``angle of attack'' that a given $\ensuremath{E_{\rm BH}}$ trajectory makes with $\ensuremath{E_{\rm quench}}$. We mentioned this angle before in the $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ diagrams; the concept here is similar -- if the angle of attack is large, all galaxies of a given starting $\ensuremath{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 $\ensuremath{R_{\rm e}}$, and the extent of the quenching locus in $\ensuremath{M_*}$ will be large. Matching the observed large width of the green valley requires acute angles of attack, as discussed further below.\footnote{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.}

Figure~\ref{Fig:SimuCartoon} shows corresponding snapshots of galaxies in the $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and $\Sig1\text{-}\ensuremath{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 $\pm0.75\sigma_{\rm Re}$ and $\pm 1.5\sigma_{\rm Re}$. They are superimposed on a shaded region denoting the observed $\pm 1.5 \sigma_{\rm Re}$ range of observed star-forming galaxies at that redshift from \citetalias{vanderWel14}.  The model galaxy points populate the observed range of galaxies fairly well, which validates our scheme for starting galaxies at $10^9~\ensuremath{\rm M_{\odot}}$ with $\ensuremath{R_{\rm e}} = 0.02 \times \ensuremath{R_{\rm vir}}$ and growing $\ensuremath{R_{\rm e}}$ in proportion to $\ensuremath{R_{\rm 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.  The straight red lines underneath them are the quenched ridgelines of \citetalias{Barro17} displaced downward by 0.2~dex in $\Sig1$ to represent the adopted quenching boundary (entry into the green valley). They are the same as the red lines shown in Figure~\ref{Fig:Mstar_Sig1_SMA1}.  The good match between the red lines and the green points verifies that the zero point and redshift evolution of the quenching expression in Eq.~\ref{Eq:EbhEhalo} 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~\ref{Sec:PowerLawModel} and are $t = 2.03$, $v = 1.76$, $s = 1.75$, $p = 1.0$, and $q = 1.1$.

\section{Interpretation of $\ensuremath{E_{\rm quench}}$ in Terms of Halo-Gas Binding Energy}\label{Sec:ThermalModel}

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 $\ensuremath{E_{\rm quench}} = k_{\rm e}(z) \left(\frac{\ensuremath{M_{\rm vir}}}{10^{12}~\ensuremath{\rm M_{\odot}}}\right)^t$ in Eq.~\ref{Eq:EbhEhalo}.\footnote{In adopting this criterion, we have neglected both radiative losses in the AGN feedback and possible energy transfer to the halo by supernovae \citep[e.g.][]{Dekel86} and by gravitational heating \citep[e.g.,][]{Dekel08,Johansson09}.  Recent hydro simulations give some sense of the importance of the first two factors.   In both EAGLE \citep{Bower17} and TNG50 \citep{Nelson19}, successful AGN winds tend to conserve energy, which justifies the use of cumulative BH energy in Eqs.~\ref{Eq:EbhEth} and \ref{Eq:EbhEbind}.  In contrast, supernovae winds are prone to large radiative losses \citep[e.g.][]{Somerville15}, and adding them to the current framework would require a major energy-loss correction factor. In TNG50, the relative roles of supernovae vs\@.~ BHs are hard to separate near \ensuremath{M_*} = $10^{10.5}~\ensuremath{\rm M_{\odot}}$, but AGN winds clearly dominate above this \citep{Nelson19}.}  The formulation of $\ensuremath{E_{\rm quench}}$ as a power law of $\ensuremath{M_{\rm 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 \emph{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 $\ensuremath{E_{\rm BH}}$ is of the same order of magnitude as the halo-gas binding energy and thus that actual BHs emit enough power to influence their halos \citep{Silk98,Ostriker05,Bower17,Oppenheimer18, Davies19a,Oppenheimer19,Nelson19}.  If in addition we can show that halo-gas binding energy varies similarly to $\ensuremath{E_{\rm quench}}$ as a function of $\ensuremath{M_{\rm 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:
 \begin{equation}\label{Eq:GasBindingEnergy}
 \ensuremath{E_{\rm bind}} = \frac{1}{2} f_{\rm bary} \ensuremath{M_{\rm vir}} V_{\rm vir}^2,
 \end{equation}
 where $f_{\rm bary}$ is the fraction of $\ensuremath{M_{\rm vir}}$ in the form of hot gas, set equal to 0.1 throughout. Inserting this into Eq.~\ref{Eq:EbhEhalo} we now have a new quenching condition:
 \begin{equation}\label{Eq:EbhEth}
  \ensuremath{E_{\rm BH}} =  f_{\rm B} \times \ensuremath{E_{\rm bind}}
\end{equation}
or
\begin{equation}\label{Eq:EbhEbind}
  \eta \frac{\ensuremath{M_{\rm BH}}}{10^8~\ensuremath{\rm M_{\odot}}} c^2 = f_{\rm B} \times \frac{1}{2} f_{\rm bary} \ensuremath{M_{\rm vir}} V_{\rm vir}^2,
\end{equation}
where $f_{\rm B}$ is the multiple of halo-gas binding energy that needs to be supplied by the black hole at the quenching boundary.

The quantity $ f_{\rm B} \times \ensuremath{E_{\rm bind}}$ is plotted as the solid blue lines for the three typical halo masses in Figure~\ref{Fig:BHEnergy}, and the equivalent quenching boundaries for this new model are the blue lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}.  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_{\rm 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 $\ensuremath{E_{\rm quench}}$ to $\ensuremath{E_{\rm bind}}$.  At fixed $z$, $\ensuremath{E_{\rm quench}}$ varies as $\ensuremath{M_{\rm vir}}^2$, while $\ensuremath{E_{\rm bind}}$ varies as $\ensuremath{M_{\rm vir}}\ensuremath{V_{\rm vir}}^2 \sim \ensuremath{M_{\rm vir}}^{1.67}$, which is similar but slightly shallower.  This slight difference accounts for the shallower slope of the blue lines in Figure~\ref{Fig:Mstar_Sig1_SMA1} and also explains why the blue line lies above the red dashed line for small galaxies but below it for massive galaxies in Figure~\ref{Fig:BHEnergy}.\footnote{We note that the slope of the blue lines could be increased in Figure~\ref{Fig:Mstar_Sig1_SMA1} 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.}

For a given galaxy over time, it may also be shown using the above equations plus $\ensuremath{M_{\rm vir}} \sim \ensuremath{V_{\rm vir}}^3/h(z)^{1.11}$ that $\ensuremath{E_{\rm bind}}/\ensuremath{E_{\rm quench}}$ evolves as $\ensuremath{M_{\rm vir}}^{-0.33}h(z)^{0.74-0.66v} \sim \ensuremath{M_{\rm vir}}^{-0.33}h(z)^{-0.42}$.  But $\ensuremath{M_{\rm vir}}$\ grows approximately as $h(z)^{-1.50}$ over this redshift range, and so $\ensuremath{E_{\rm bind}}/\ensuremath{E_{\rm quench}} \sim h(z)^{0.08}$.  In other words, $\ensuremath{M_{\rm 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~\ref{Fig:BHEnergy}.  An additional insight comes from noting that $\ensuremath{E_{\rm quench}}$ varies as $k_{\rm e}(z)\ensuremath{M_{\rm vir}}^{2.0}$ while $\ensuremath{E_{\rm bind}}$ varies as $\ensuremath{M_{\rm vir}}^{2.0}\ensuremath{R_{\rm vir}}^{-1}$.  Thus, the evolving term $k_{\rm e}(z)$ in the empirical power-law model is needed to mimic the evolution of $\ensuremath{R_{\rm vir}}^{-1}$ over time.  The question of the zero-point evolution of the \emph{boundary} over time is different from the evolution of an individual galaxy and is considered in Section~\ref{SubSec:Successes}.

This is a good place to make a point about the uniformity of quenching.  If $f_{\rm B}$ varies from halo to halo, say, with rms scatter $\delta_{\rm {fB}}$ dex at a given epoch, this will induce the same amount of scatter in $\ensuremath{M_{\rm BH}}$ and a scatter of $\delta_{\rm {fB}}/v$ in $\Sig1$.  The $\Sig1$ ridgeline has rms scatter equal to 0.16~dex today \citepalias{Fang13}, which implies that $\delta_{\rm {fB}} = 0.28$ dex.  However, this is an upper limit since age stratification has been detected \citep{Tacchella17,Luo19}, so the real scatter in $\delta_{\rm {fB}}$ must be less.  The point is that the halo binding-energy criterion in Eq.~\ref{Eq:EbhEbind}, represented by $f_{\rm B}$, must be \emph{highly uniform from halo to halo} in order to create the tight ridgelines in $\Sig1\text{-}\ensuremath{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 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 vs\@.~ $\ensuremath{M_{\rm 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 \citet{Booth10,Booth11} and has been seen and analyzed in the EAGLE simulations \citep{Davies19a,Oppenheimer19} and in IllustrisTNG (\citealp{Terrazas19}).  Comparisons to these and other models are discussed in Section~\ref{Sec:Discussion}.

\section{Variants of the Empirical Power-Law Model}\label{Sec:Variants}

This section describes four variations on the empirical power-law model.  Summaries are presented here, and details are given in Appendix~\ref{App:A4}.  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~\ref{Sec:Discussion}.

\subsection{Variant~1: A Broken Power-Law for the SMHM Relation} \label{SubSec:BrokenPLModel}

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 $\ensuremath{M_*} > 10^{12}~\ensuremath{\rm M_{\odot}}$ at late redshifts in Figure~\ref{Fig:SimuCartoon}.  This tongue is due to the fact that we have used the curved HMSM relation from \citetalias{Rodriguez17} at $z = 0$ to find the starting redshifts, $z_{\rm s}$, but have calculated stellar masses going forward using the power-law SMHM relation of Eq.~\ref{Eq:MstarMvir}, $\ensuremath{M_*} = k_{\rm m} \ensuremath{M_{\rm vir}}^s$.  Since the latter keeps climbing with $\ensuremath{M_{\rm vir}}$, massive galaxies that start near the knee of the real SMHM relation today can wind up with higher stellar mass than their assumed initial mass.  Galaxies with $\ensuremath{M_*} > 10^{10.7}~\ensuremath{\rm M_{\odot}}$ and large $\ensuremath{R_{\rm 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}~\ensuremath{\rm M_{\odot}}$ but switches to a flatter power-law for SMHM above that. Results are shown in Figures~\ref{Fig:HMSM} and \ref{MEnergy_V1} in Appendix~\ref{SubSec:AppBrokenPLModel}.  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 \citep{Luo19}.  The quenched ridgeline in $\Sig1\text{-}\ensuremath{M_*}$ in this new sample is seen to bend over to lower $\Sig1$ 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}~\ensuremath{\rm M_{\odot}}$, and it is seen from Figure~\ref{Fig:SimuCartoon} 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 vs\@.~ 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 $\ensuremath{M_*}$ is difficult to predict.


\subsection{Variant~2: Move the Evolving Zero Point from $\ensuremath{E_{\rm quench}}$ to the SMHM Relation}\label{SubSec:EvolvingSMHMZP}

We have noted that the quenching boundary evolves in $\Sig1\text{-}\ensuremath{M_*}$, and the empirical  power-law model matches this by evolving the zero point, $k_{\rm e}(z)$, of $\ensuremath{E_{\rm quench}}$ in Eq.~\ref{Eq:EbhEhalo}.  An alternative approach is to move the evolving zero point to the SMHM relation in Eq.~\ref{Eq:MstarMvir} and remove it from $\ensuremath{E_{\rm quench}}$.  The expression that we have used for the SMHM relation in this case is derived in Appendix~\ref{SubSec:AppEvolvingSMHMZP}.  It makes the zero point of $\ensuremath{M_*}$ vs\@.~ $\ensuremath{M_{\rm vir}}$ larger by about 0.5~dex at $z = 0$ compared to $z = 2$.  This grows BHs more quickly, which allows them to overcome halos without resort to $k_{\rm e}(z)$.  Changes in the SMHM relation back in time that bear on this are described in Appendix~\ref{SubSec:AppEvolvingSMHMZP}.  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 vs\@.~ galaxy mass is different, as shown by the dashed lines in Figure~\ref{Fig:TrajectoryM}.  BHs grow more rapidly at late times, so quenching starts too late, i.e., not until around $z = 1$ (see Figure~\ref{MEnergy_V2} in Appendix~\ref{SubSec:AppEvolvingSMHMZP}). Trying to move the onset of quenching earlier by increasing $\eta$ in Eq.~\ref{Eq:EbhEhalo} produces over-quenching at late times.  Also, the more rapid growth in $\ensuremath{E_{\rm quench}}$ with time no longer parallels the slower growth in halo-gas energy content, $\ensuremath{E_{\rm bind}}$, as shown by the differing red and blue lines in Figure~\ref{MEnergy_V2} in Appendix~\ref{SubSec:AppEvolvingSMHMZP}, and so the appealing physical argument equating $\ensuremath{E_{\rm 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 $\ensuremath{E_{\rm 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.


\subsection{Variants~3 and 4: Mass Rates or Energy Totals?}\label{SubSec:RatesorTotals}

A major difference between the empirical  power-law model and other halo-heating models is that the latter often use heating and cooling \emph{rates} as opposed to cumulative energy totals as in Eq.~\ref{Eq:EbhEhalo}.  Use of rates is especially common in semi-analytic models \citep[e.g.,][]{Bower06,Croton06,Somerville08,Henriques15}.  Variants~3 and 4 explore whether Eq.~\ref{Eq:EbhEhalo} can be recast in terms of rates. In both variants, the halo energy quota, $\ensuremath{E_{\rm quench}}$, on the right side of Eq.~\ref{Eq:EbhEhalo} is replaced by the halo mass-cooling rate, $\ensuremath{\dot{m}_{\rm cool}}$, from \citet{Croton06} \citep[also][]{White91,Cole00,Somerville08,Guo11,Henriques15}.

To replace total BH energy on the left side of Eq.~\ref{Eq:EbhEhalo} 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 $\ensuremath{M_{\rm BH}}$ values from the empirical power-law model.  In Variant~4, the BH accretion expression introduced by \citet{Croton06} for radio mode is used, but as this has no dependence on \ensuremath{R_{\rm e}}, this model is only a one-parameter family.  Equations are given in Appendix~\ref{SubSec:AppRatesorTotals}, and the accretion rates vs\@.~ time are compared in Figure~\ref{Fig:Trajectory_BH}.  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~\ref{Fig:MRate_V3}-\ref{Fig:SimuSize_V4}.  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 2-D structural diagrams of galaxies, which automatically excludes one-parameter models that depend only on mass, like Variant 4.

\section{Discussion}\label{Sec:Discussion}

The Discussion has nine sections.  The first reviews and explains the successes that were  ascribed to the model in Table~\ref{Tb_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 latter open the way to thinking about the BH-halo connection, what is meant by  BH-galaxy ``co-evolution'', and the  implications of having \emph{three separate yet matching} quenching boundaries, in $\Sig1\text{-}\ensuremath{M_*}$,  $\Sig1\text{-}\ensuremath{M_{\rm BH}}$, and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$.  We close with a comparison of our picture to the \citetalias{Lilly16} quenching model.

\subsection{Scorecard for the Binding-Energy Model}\label{SubSec:Successes}

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

Column~1 in Table~\ref{Tb_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_{\rm B}$ in Eq.~\ref{Eq:EbhEbind}, which has been set equal to 4 to make the zero point of the blue lines match that of the red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1} on average.\footnote{The adopted value of $f_{\rm B}$ depends on the assumed value of the black hole efficiency factor $\eta$, which we have set to 0.01, and we are actually adjusting the ratio of these two quantities.}  Aside from the general insistence on power-laws, none of the observations in Table~\ref{Tb_1} were used to determine the five basic ingredients above, and hence all of the successes in Table~\ref{Tb_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 $\ensuremath{M_*}$ and $\ensuremath{\langle R_{\rm e}\rangle}$ derived from halo $\ensuremath{M_{\rm vir}}$ and $\ensuremath{R_{\rm vir}}$.  As pointed out by \citetalias{vanderWel14} and noted in Section~\ref{SubSec:EvolTracks}, this halo-based approach matches rather well the ridgeline slope, zero-point evolution, and evolutionary track slopes of star-forming galaxies in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and, through the $p$,$q$ mapping, also in $\Sig1\text{-}\ensuremath{M_*}$.  These are lines~4, 5, 8, and 9 in Table~\ref{Tb_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 shallow in $\Sig1\text{-}\ensuremath{M_*}$ and a bit too steep in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{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~\ref{Fig:Mstar_Sig1_SMA1}).  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 $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ relation for quenched galaxies is offset to higher values at early times  \citep{Peng06,Treu07,Woo08,Gu09,Decarli10,Bennert10,Bennert11,Ding17,Ding20}.  This offset has often been interpreted to mean that BHs form \emph{before} the stars of their host galaxies \citep[e.g.,][]{Peng06,Treu07,Bennert11,Ding17}, 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 black holes form \emph{after} most of the stars. The model has a different explanation for the offset, namely, that the quenching boundary in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ -- and therefore its associated ridgeline -- was higher in the past.  The predicted motion can be calculated from the observed motion of the $\Sig1\text{-}\ensuremath{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 $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$  \citep{Bennert11,Ding17,Ding20}.  The point is that offsets of $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ with time are not necessarily due to growth-time differences between stars and black holes but may instead reflect motions of the quenching boundary.\footnote{The fact of a moving zero point of $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_*}$ does not contradict our  fundamental assumption that the zero point of  $\ensuremath{M_{\rm BH}} = \rm{const.}  \Sig1^{1.76}$  remains the same since $\ensuremath{M_*}$ and $\Sig1$ are not the same quantity.} Tnhe same point was made by \citet{Caplar18}, who explored a model similar to ours in which black holes grow rapidly just as galaxies enter the green valley.

Moving on, line 16 in Table~\ref{Tb_1} says that the observed zero point of the quenched ridgeline in $\Sig1\text{-}\ensuremath{M_*}$ evolves downward as $h(z)^{0.66}$ (Appendix~\ref{App:A2}), while the observed zero point of the star-forming ridgeline in the same plane evolves downward by the same amount (\citetalias{Barro17}).  This is the key equality that maintains star-forming galaxies \emph{in the same relationship} to quenched galaxies in $\Sig1\text{-}\ensuremath{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 \emph{over the same range of masses at all redshifts}. This is seen in the predicted model quenching loci of Figure~\ref{Fig:SimuCartoon}, where the smallest galaxies always enter quenching near $10^{9.5}~\ensuremath{\rm M_{\odot}}$ and are fully quenched by $10^{10}~\ensuremath{\rm M_{\odot}}$, which is the observed lower edge of the high-mass quenching channel at all redshifts \citepalias{Barro17}.  Finally, this same constant vertical offset between the quenching boundary and the star-forming cloud ensures a \emph{smooth, steady flow of galaxies across the boundary at all times} (line~3 in Table~\ref{Tb_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 $\Sig1\text{-}\ensuremath{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 $\Sig1\text{-}\ensuremath{M_*}$ goes as $\ensuremath{R_{\rm e}}^{-p} \sim \ensuremath{R_{\rm e}}^{-1} \sim \ensuremath{R_{\rm vir}}^{-1} \sim h(z)^{0.74}$ (this last is from Eq.~\ref{Eq:ReMstar} in Section~\ref{SubSec:EvolTracks}).  The zero-point motion of the quenching boundary (and thus the quenched ridgeline) varies as $\ensuremath{E_{\rm bind}}^{1/v} \sim (\ensuremath{M_{\rm vir}}^2/\ensuremath{R_{\rm vir}})^{1/v}$, which at fixed $\ensuremath{M_*}$ goes as $\ensuremath{R_{\rm 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~\ref{Fig:Mstar_Sig1_SMA1}, leaving a residual in $\Sig1$ 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, black hole masses (at fixed $\ensuremath{M_*}$) grow smaller at late times because galaxies have larger $\ensuremath{R_{\rm e}}$.  However, halo binding energies (at fixed \ensuremath{M_{\rm vir}}) also become smaller with time because halos have larger \ensuremath{R_{\rm vir}}.  Together, these two factors cause the downward boundary motion.  This point is not new; it was made by \citet{Booth11}, who predicted that the quenching zero point of models in which BH energy output is balanced against halo-gas binding energy would decline with time. Our new insight is that the parallel evolution of black holes 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 black holes 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~\ref{SubSec:Quenching}).

A moving quenching boundary has another consequence.  If the evolution of individual galaxies \emph{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'' \citep[e.g.,][]{Carollo13, Poggianti13}.

\citetalias{Barro17} listed three factors that help quenched galaxies stay close to the downward-evolving ridgeline in $\Sig1\text{-}\ensuremath{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 $\Sig1\text{-}\ensuremath{M_*}$ will be only moderate.

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

Age stratification in the quenched ridgeline has been measured in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ by \citet{Martin18} and in $\Sig1\text{-}\ensuremath{M_*}$ by \citet{Tacchella17} and \citet{Luo19}.  Early-quenched galaxies fail to fully keep up with the moving boundary in both spaces, but the mismatch in $\Sig1\text{-}\ensuremath{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 $\ensuremath{M_*}$.  As described in Section~\ref{SubSec:EvolModelResults}, this is explained by the highly acute ``angles of attack'' of the evolutionary tracks relative to the quenching boundary in $\Sig1$ vs\@.~ $\ensuremath{M_*}$, $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$, and $\ensuremath{E_{\rm BH}}$ vs\@.~ $\ensuremath{E_{\rm quench}}$ (and also $\ensuremath{E_{\rm bind}}$), which enable star-forming galaxies with the same $\ensuremath{M_*}$ but only slightly different values of $\ensuremath{R_{\rm e}}$ and $\Sig1$ to evolve to very different final quenching masses.  These acute angles are a natural feature of the model, and line~2 is accordingly marked as a success.

Line~14 in Table~\ref{Tb_1} notes that the ridgelines of quenched galaxies are narrow in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ but broader in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, especially at high $z$.  The basic narrowness in $\Sig1$ and $\ensuremath{M_{\rm BH}}$ is due to the fact that galaxies stop evolving and pile up along these loci.  The broader width of $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ reflects the greater impacts on $\ensuremath{R_{\rm e}}$ of compaction (see Section~\ref{SubSec:Compaction}) and dry minor merging \citep[e.g.,][]{Daddi05,Trujillo06,vanDokkum08,Naab09,vanDokkum10,Szomoru12,Oser12}.  The evolution of the latter goes as $\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^2$, which is steeper than the quenched ridgeline itself ($\ensuremath{R_{\rm e}} \sim \ensuremath{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~\ref{Tb_1} says that $\ensuremath{M_{\rm BH}}$ varies roughly as some power $\ensuremath{M_*}^w \ensuremath{R_{\rm e}}^{-w}$ (or as $\ensuremath{\sigma_{\rm e}}^{2w}$) for quenched and star-forming galaxies together but as a smaller power than $w$ for quenched galaxies alone (compare slopes from \citet{vandenBosch16} and \citetalias{Kormendy13}).  This difference turns out to be the same problem as the low value of $u = 3.52$ in Section~\ref{SubSec:InputPowerLaws}.  Both issues are discussed in Section~\ref{SubSec:PredictionsTests}, 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 $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ than in $\Sig1\text{-}\ensuremath{M_*}$, and both lines~13 and 18 are marked as successes.

The above discussion covers all the items without brackets in Table~\ref{Tb_1}.  Lines~6 and 7 (with brackets) involve the slope and zero point of the quenched ridgeline in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$.  As these depend on how $\ensuremath{R_{\rm 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 $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ vs\@.~  $\Sig1\text{-}\ensuremath{M_*}$ (line~17) is attributable to the larger impact of post-quenching dry merging on $\ensuremath{R_{\rm e}}$ than on $\Sig1$, and so partial credit is given to that line as well.

This completes the review of all the items in Table~\ref{Tb_1}.  Yet another test would use X-ray luminosities to test whether star-forming galaxies with higher $\Sig1$ (and lower $\ensuremath{R_{\rm e}}$) at fixed $\ensuremath{M_*}$ have larger black holes. The findings of \citet{Ni19} suggest this, but their star-forming sample also includes green valley galaxies and thus is not quite pure enough.  A re-analysis of a pure main-sequence sample is underway (Kocevski et al., in prep.).  X-ray data for green valley galaxies are discussed in Section~\ref{SubSec:BoundaryCrossing}.

\subsection{The Origin of Scatter in Effective Radii}\label{SubSec:RadiusScatter}

A basic ingredient of the binding-energy model is that star-forming galaxies at a fixed stellar mass have different radii.  Radii differences correlate with central $\Sig1$ (Figure~\ref{Fig:Mstar_Sig1_SMA1}), and these induce differences in BH mass and BH feedback.  Where do radii differences come from?  Since the model assumes that all star-forming galaxies of the same stellar mass have the same dark-halo mass, $\ensuremath{R_{\rm e}}$ differences must come from some other parameter.

The connection between dark halos and galaxy radii was reviewed by \citet{Somerville18} and \citet{Wechsler18}.  \citet{Mo98} introduced the formula for the exponential disk scale length $r_{\rm d} = 1/\sqrt{2}\lambda\ensuremath{R_{\rm vir}}$, where $\lambda$ is the dimensionless halo spin parameter.  This suggested that $\lambda$-variations might cause scatter in $\ensuremath{R_{\rm e}}$.  However, the observed rms scatter in $\ensuremath{R_{\rm e}}$ is only 0.2~dex \citepalias{vanderWel14} compared to 0.25~dex for $\lambda$ in the models, which leaves no room for any other sources of scatter or measurement errors \citep{Desmond15}.  Hydro simulations also show that cooling gas gains or loses angular momentum unpredictably while falling onto galaxies  \citep{Pichon11,Stewart13,Danovich15,Stevens17,Jiang18}, blurring any connection between original and final angular momenta.  The net result is that galaxy radius is found to be almost uncorrelated with halo $\lambda$ in EAGLE \citep{Desmond17}, and  \citet{Jiang18} find at most poor correlations for VELA \citep{Ceverino14} and for NIHAO \citep{Wang15}.

An interesting alternative to $\lambda$ is halo concentration, $c$.  \citet{Jiang18} find that the radii of VELA and NIHAO galaxies obey $\ensuremath{R_{\rm e}} = 0.02 \ensuremath{R_{\rm vir}}(c/10)^{-0.7}$, with galaxies being smaller in higher-concentration halos.  This is plausible  since 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 black holes owing to their higher stellar densities.  For the purposes of this paper, the origin of radius scatter needs no explanation -- the $\Delta \log \ensuremath{R_{\rm 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 \ensuremath{R_{\rm e}}$ stay constant during the star-forming phase so that galaxies evolve on parallel tracks in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$.  However, a picture in which the two leading structural parameters of star-forming galaxies, $\ensuremath{M_*}$ and $\ensuremath{R_{\rm e}}$, arise directly from two leading parameters of dark halos, \ensuremath{M_{\rm 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 black holes, larger  values of $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$, and are more quenched today \citep{Davies19a,Davies19b,Oppenheimer19}.  The effect in EAGLE is large -- halos in the top quartile of $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$ were quenched several billion years ago while halos in the bottom quartile are still making stars today.  Furthermore, \citet{Davies19a} showed that it is the \emph{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 $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$.

These properties support our assumption in the binding-energy model that quenching is determined by the ratio $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$, not by instantaneous heating or cooling rates.  EAGLE quenching becomes perceptible when $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}} \sim 2$ and is complete when $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}} \sim 20$ whereas we predict the range 4-to-20 for this ratio (see Section~\ref{SubSec:BoundaryCrossing} below).  This is a promising level of agreement given the simplicity of our model.\footnote{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 \citep[e.g.,][]{Tumlinson17}, and so they must refill.  EAGLE and IllustrisTNG clusters in fact do refill \citep{Oppenheimer19,Davies19b}.  A different form of BH feedback must take over to to maintain quenching at the centers of gas-rich groups and clusters. This later evolutionary stage is not covered by the model.}

Guided by EAGLE and IllustrisTNG, it is tempting to create a revised binding-energy model by adding halo concentration to $\ensuremath{M_{\rm vir}}$ as a secondary parameter.  This addition of a second halo parameter in addition to $M_{\rm halo}$ would be an example of \emph{assembly bias} \citep{Wechsler18}.  However, doing so would raise important questions.  For one thing, N-body simulations show that mass accretion rates onto halos at fixed time vary with halo concentration \citep{Lin19}, which, if unopposed, might disturb the observed absence of any trend in star-formation rate vs\@.~ galaxy size at fixed $\ensuremath{M_*}$ (\citealp{Omand14,Whitaker17}; \citetalias{Fang18}; \citealp{Lin19}).  Indeed, star-formation rates \emph{do} vary with size in EAGLE \citep{Oppenheimer19}, in disagreement with the observations.  A dodge around this could posit weaker galactic winds in small galaxies on account of their higher escape velocities \citep{Dutton10,Lin19}.  But even if \emph{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 \emph{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 \emph{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.

\subsection{Quenching in the Binding-Energy Model Compared to Other Models}\label{SubSec:Quenching}

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:

\emph{Stellar feedback}. This is the most basic feedback mechanism, present in all models.  Energy/momentum from young massive stars and SNae drive a wind that expels gas from a galaxy during the star-forming phase.

\emph{Quasar mode} \citep{Croton06,  Somerville08, Henriques15}, also called ``bright mode'' or ``high mode'' \citep{Somerville08,Pillepich18} or ``radiative mode'' \citep{Somerville15}. This is the  principal mode of BH growth in many models, triggered by mergers \citep{Hopkins06}, disk instabilities \citep{Cattaneo06}, or other  unspecified secular processes \citep{Martin18}. Key features of quasar mode (in most models) are a high black-hole accretion rate, formation of a thin accretion disk, and high radiative efficiency from the accretion disk (see review by \citealt*{King15}).  It is generally quasar-mode winds that quench star formation in semi-analytic models and in many hydro models.

\emph{Radio mode} \citep{Croton06,  Somerville08, Henriques15}, also called ``kinetic mode'' or ``low mode'' \citep{Weinberger17,Pillepich18}, or ``jet mode'' \citep{Somerville15}.  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 black-hole accretion rate, formation of a thick accretion disk, and low radiative efficiency (see review by \citealt*{Yuan14}).  Radio mode may quench, as in IllustrisTNG \citep{Weinberger17,Davies19b,Terrazas19,Nelson19}, or it may maintain a low star-formation rate in gas-rich halos after other mechanisms have quenched.

\emph{Halo quenching} \citep{Rees77, Blumenthal84,Silk86,Birnboim03,Keres05,Dekel06,Dekel19}. As halos grow in mass, their gas becomes hotter and more dilute, which retards cooling.  A critical halo mass has been proposed, $\ensuremath{M_{\rm crit}} \sim 10^{11.5\text{-}12}~\ensuremath{\rm M_{\odot}}$ \citep{Dekel06,Dekel19}, 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'' \citep{Keres05,Dekel06,Keres09}, which are less able to penetrate to the central galaxy.

\emph{Late infall dry-up}. The specific accretion rate of gas onto halos declines after $z \sim 1$ in all cosmologies, regardless of $\Omega_{\rm m}$ or $\Omega_{\Lambda}$ \citep{Dekel14}.  Though not a feedback process per se, it lowers the amount of gas available for star formation \citep{Bouche10, Dave11, Dave12, Krumholz12, Dekel13, Lilly13, Dekel14, Forbes14, Feldmann15, Mitra15, Rodriguez16_1} and retards halo cooling by keeping gas densities low at late times.  Along with the transition from cold streams to hot streams \citep{Keres05,Dekel06}, this makes late halos easier to quench and would help to lower the zero point in the quenching boundary relation (red lines) in $\Sig1\text{-}\ensuremath{M_*}$, as seen in Figure~\ref{Fig:Mstar_Sig1_SMA1}.

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 \emph{rate} of BH gas accretion vs\@.~ the \emph{rate} of halo-gas cooling in some models \citep[e.g.,][]{Croton06,Bower06,Somerville08,Henriques15} whereas the binding-energy model integrates \emph{total} BH energy over time and compares it to the \emph{total} binding energy of the halo gas.  The impact of the difference was explored in Variants~3 and 4 in Section~\ref{Sec:Variants}, where it was shown that rate formulations do not feed gas to black holes at the right rate.  As noted above, two recent hydro simulations also seem to exhibit the total-energy criterion (EAGLE: \citet{Davies19a,Oppenheimer19} and IllustrisTNG: \citet{Davies19b,Terrazas19}).

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 \citet{Somerville08} SAM \citep{Pandya17} and in TNG50 \citep{Nelson19}. In contrast, entry into the green valley in the binding-energy model is the well-defined point when total AGN output surpasses $f_{\rm B} \times$halo-gas binding energy, where $f_{\rm B} \simeq 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\ensuremath{M_{*,\rm GV}}\rangle$, over time. A plausible value for this mass, cited in many works, is the value of $\ensuremath{M_*}$ at the critical halo cooling mass $\ensuremath{M_{\rm crit}}$ \citep[e.g.,][]{Birnboim03,Keres05,Dekel06,Dekel19}.  A second characteristic mass has also emerged in recent hydro models from balancing supernovae wind energy against galaxy or halo potential well depth \citep{Dubois15,Angles17,Bower17,Habouzit17,Davies19a,Dekel19,Oppenheimer19,Henriques19}. SNae wind failure allows central gas to accumulate \citep[e.g.,][]{Dekel86}, which fuels rapid BH growth and causes quenching. \citet{Dekel19} even argue that these two masses converge to create a single ``golden mass'' at early times.

In contrast, the constancy of $\langle\ensuremath{M_{*,\rm GV}}\rangle$ in our model has nothing to do with $\ensuremath{M_{\rm crit}}$ or SNae winds but comes instead from the fact that halos and galaxies evolve so as to keep the quenching boundary \emph{located in the same place relative to the blue cloud} over time.  This means that quenching galaxies always exit the blue cloud with the \emph{same characteristic mass}.  Granted, it is fortunate that $\ensuremath{M_{\rm crit}}$ for hot halos is near $10^{11.5\text{-}12.0}~\ensuremath{\rm M_{\odot}}$ and not larger, since having a hot halo may facilitate BH quenching.  Central SNae 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 SNae winds to keep galaxy centers gas free, but by the \emph{cumulative effect of black hole feedback on the halo gas}.  The characteristic masses set by halo cooling and by SNae 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 \emph{is} a characteristic halo mass in the binding-energy model, and the answer is yes.  This is the halo mass, $M_{\rm vir,cross}$ where the growing energy from the BH in a typical galaxy (with \ensuremath{R_{\rm e}} = \ensuremath{\langle R_{\rm e}\rangle}) during the star-forming phase crosses the binding energy of halo gas.  Scaling laws in the Appendix may be combined with the quenching condition in Eq.~\ref{Eq:EbhEbind} to yield:
\begin{eqnarray}\label{Eq:MvirCross}
\left(\frac{M_{\rm vir,cross}}{10^{12}~\ensuremath{\rm M_{\odot}}}\right)^{1.60v-1.67} & = & 9.5  k_{\rm{Re}}^{v}  k_{\rm{BH}}^{-1}   k_{\rm{\Sig1}}^{-v}   k_{\rm{Rvir}}^{v-1}   k_{\rm m}^{-1.1v} \frac{f_{\rm B}}{4} \times \nonumber\\
 & & \frac{f_{\rm hot}}{0.1} \times \left(\frac{\eta}{0.01}\right)^{-1} \times  h(z)^{0.74(v-1)},
\end{eqnarray}
where the constants $k$ are the zero points of the five relevant scaling laws from Appendix~\ref{App:A3} and we have inserted our best values for $p$, $q$, $s$, and $t$.   $M_{\rm {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 black-hole building are fully understood. Two more free parameters at the present time are $f_{\rm B}$ and $\eta$. These also should be expressible (someday) using fundamental constants.

For reference, we insert $k$-values from the Appendix, use $v = 1.76$, and substitute $\ensuremath{M_*}$ for $\ensuremath{M_{\rm vir}}$ to get the crossing mass in terms of $\ensuremath{M_*}$, which should equal to $\ensuremath{M_*}$ in the middle of the green valley.  We obtain:
\begin{eqnarray}\label{Eq:MstarCross}
\frac{M_{\rm *,cross}}{10^{10}~\ensuremath{\rm M_{\odot}}} & = & 4.0 \left( \frac{f_{\rm B}}{4} \right)^{1.54} \times \left( \frac{f_{\rm hot}}{0.1}  \right)^{1.54} \times \nonumber\\
  & & \left( \frac{\eta}{0.01} \right)^{-1.54} \times \left( \frac{h(z)}{2.3} \right)^{-0.86},
\end{eqnarray}
where $h(z)$ has been normalized by its value at $z = 1.5$, where $k_{\rm {Rvir}}$ is calculated.  $M_{\rm *,cross}$ therefore equals $10^{10.6}~\ensuremath{\rm M_{\odot}}$ at $z = 1.5$, which is a good match to the green valley stellar mass.  However, note that $M_{\rm *,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~\ref{Fig:Mstar_Sig1_SMA1} that was pointed out in Section~\ref{SubSec:Successes}.  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 SNae 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.  Since 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 $\Sig1\text{-}\ensuremath{M_*}$ in Figure~\ref{Fig:Mstar_Sig1_SMA1}.


\subsection{Compaction in the Binding-Energy Model}\label{SubSec:Compaction}

Implicit in the binding-energy picture is that galaxies, halos, and black holes 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 \emph{wrong} side of the boundary at high redshift but then disappear at low redshift.  This is especially true using $\ensuremath{R_{\rm e}}$ as a measure of density as opposed to $\Sig1$ (cf. red points below the boundary in the top row of Figure~\ref{Fig:Mstar_Sig1_SMA1}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$  (\citealp{Barro13,Barro14,Nelson14,Williams14,vanDokkum15}; \citetalias{Barro17}).  They are thought to be the star-forming precursors of the similarly dense, quenched ``red nugget'' galaxies \citep{Daddi05, Trujillo06, vanDokkum08, Damjanov11}, which are visible in Figure~\ref{Fig:Mstar_Sig1_SMA1} at the same redshifts. VELA simulations \citep{Ceverino14} 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 \citep{Dekel14,Zolotov15,Tacchella16b}. This scenario says that star-forming galaxies compactified to form blue nuggets, which then faded to become red nuggets. Simulations by \citet{Choi18} and by the IllustrisTNG team \citep{Habouzit19} produce similar objects.

Actually, VELA galaxies never fully quench owing to their lack of central black holes.  But X-ray data show a large population of X-ray sources near the quenching boundary at $z \sim 2$ \citep{Kocevski17}, right where maximum compaction occurs \citep{Dekel19}. If these accreting black holes were added to VELA, their feedback might achieve quenching.  In contrast to VELA, BHs \emph{are} present in IllustrisTNG and in the \citet{Choi18} models, and in these simulations compact blue nuggets quench.

A variety of factors can trigger compaction events in central galaxies, including mergers \citep{Barnes91, Hernquist91, Hopkins06}, violent-disk-instability-driven inflows \citep{Dekel14}, counter-rotating streams \citep{Danovich15}, a tri-axial halo core \citep{Tomassetti16}, and the return of recycled low-angular-momentum gas \citep{Elmegreen14}.  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 \citep{Tacchella16a}.

We take this as the basic picture for compaction in the binding-energy model, to which we add black holes.  As blue nuggets compactify, their central densities rise rapidly and their BHs grow according to the green valley vectors in Figure~\ref{Fig:Cartoon}. Compaction thus sweeps galaxies \emph{rapidly} over the quenching boundary at high redshift, at which point their growing BHs begin to impede star formation according to Eq.~\ref{Eq:EbhEbind}.

Adding black holes in this way solves four problems.  First, the inclusion of black holes makes galaxies quench.  Second, quenching according to the binding-energy equation produces \emph{sloping} quenching boundaries, which are seen at high redshift as well as at low redshift (Figure~\ref{Fig:Mstar_Sig1_SMA1}).  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 \emph{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., \citealt{Yesuf14,Barro17,Woo15,Wu18}; \citealt*{Woo19}).  Fast-track galaxies, sometimes called ``post-starburst galaxies'' \citep{Dressler83, Balogh99, Dressler99}, exhibit strong Balmer absorption lines and small radii, which indicate rapid quenching and strong gaseous dissipation \citep{Whitaker12,Yesuf14,Yano16,Almaini17,Wu18}.  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 black holes can grow and quench. Slow-track galaxies in contrast slide over the boundary gradually, and their star-formation rates never get strongly out of synch with their structural properties and BH masses.  Since the fast-track scenario requires lots of gas, the fraction of fast-track to slow-track galaxies should decrease with time, as observed (\citealp{Whitaker12,Yano16}; \citetalias{Barro17}).  In conclusion, although fast-track galaxies are not explicitly included in the binding-energy model, modest extensions to it seem able to accommodate them.


\subsection{Evolution Through the Green Valley}\label{SubSec:BoundaryCrossing}

The manner in which galaxies evolve through the green valley is important in the model because the slope of \emph{fully-quenched galaxies} is assumed for BH growth in \emph{star-forming galaxies} (Section~\ref{SubSec:InputPowerLaws}), 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 (cf.~Figure~\ref{Fig:Cartoon}).

Brief evidence to support the steepness and uniformity of GV tracks was presented in Section~\ref{SubSec:InputPowerLaws}, which we expand on here.  The first fact is the sharp ``elbow'' seen in $\Delta \log SSFR$ vs\@.~ $\Delta \log \Sig1$ by \citetalias{Fang13} and again in Figure~6b of \citetalias{Barro17}, in which $SSFR$ declines steeply through the GV while $\Sig1$ hardly changes (see cartoon version of this in Figure~\ref{Fig:Elbow}a here).  It is this pile-up at near-constant $\Sig1$ that creates the narrow ridgeline of quenched galaxies in $\Sig1\text{-}\ensuremath{M_*}$ (cf.~Figure~\ref{Fig:Mstar_Sig1_SMA1}).  The constant morphology of the elbow (\citetalias{Fang13,Barro17}; \citealp{Lee18,Mosleh18}) sets the growth of $\log \Sig1$ through the GV at roughly 0.2~dex at all masses and redshifts $z \leq 3$.  Consistent with this, the GV entry locus also lies a constant 0.2~dex below the quenched ridgeline for SDSS galaxies at all masses today \citep{Luo19}.  In summary, there seems to be a universal $\Sig1$ ``elbow diagram'' such that the growth of $\Sig1$ through the GV is small, only about $\sim0.2$~dex, while $SSFR$ falls by much more than this  (\citetalias{Fang13,Barro17}; \citealt{Luo19}).

The evolution of black-hole mass through the green valley is different.  \citet[][Figure 3]{Terrazas17} have shown that BH mass \emph{continues to grow} through the GV even as star formation declines.  That is, $\Delta \log SSFR$ vs\@.~ $\Delta \log \ensuremath{M_{\rm BH}}$ does not show a sharp elbow but rather a smooth, monotonic decline, with $\ensuremath{M_{\rm BH}}$ increasing across the green valley by about 1.5~dex (cf.~Figure~\ref{Fig:Elbow}b here).  The same difference is also seen in  Figure~\ref{Fig:ThreePanel}c, where star-forming and quenched galaxies are well separated, in contrast to $\Sig1$ in Figure~\ref{Fig:ThreePanel}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_{\rm X}$ rises by a factor of 2\text{-}3 compared to the star-forming main sequence \citep[Figure 10c]{Aird19}. 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 GV based on X-ray fluxes has also been deduced by \citet{Delvecchio19}, and a similar picture in which black holes accrete most of their mass just before quenching has been modeled by \citet{Caplar18}.

To summarize, in the binding-energy model, $\Sig1$ virtually stops growing as galaxies enter the green valley whereas $\ensuremath{M_{\rm BH}}$ grows considerably.  This is the behavior captured by the steep green tracks in Figure~\ref{Fig:Cartoon}.  The different tracks from left to right represent differences in $\ensuremath{R_{\rm e}}$ during the star-forming phase: galaxies with smaller initial $\ensuremath{R_{\rm e}}$ enter the green valley sooner at lower values of $\Sig1$ and $\ensuremath{M_{\rm 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~\ref{Fig:Cartoon} resembles BH growth in the following hydro models: RAMSES/Horizon-AGN \citep{Dubois15,Habouzit17,Dekel19}, EAGLE \citep{Bower17}, FIRE \citep{Angles17}, IllustrisTNG \citep{Habouzit19}, and MARVEL \citep{Bellovary19}.  BH growth in these models is retarded during star formation due to SNae feedback, wandering motions of BHs in shallow central potential wells, and low central gas densities. Retardation ends when SNae are no longer able to evacuate the central regions (for various reasons) and when growing central stellar surface density captures the black hole and concentrates the gas around it.  These conditions will also occur in the binding-energy model, since declining star formation at entry to the GV means less SNae feedback, and higher $\Sig1$ 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 black holes near the end of star formation in hydro models.

The late growth of black holes implies that BH masses lie \emph{below} the blue line in Figure~\ref{Fig:Cartoon} 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 \emph{boundary}  that is crossed by galaxies that are moving upward in Figure~\ref{Fig:Cartoon} into the green valley.  This would  profoundly alter the meaning of the $\ensuremath{M_{\rm BH}} \sim \Sig1^v$ scaling law, as discussed further in Section~\ref{SubSec:BHGrowth}.

In the above, we have interpreted both loci in Figure~\ref{Fig:Elbow} as evolutionary tracks that galaxies traverse while evolving from the blue cloud through the green valley and onto the red sequence.  Panel b, showing black-hole growth, is a re-rendering of data from \citet{Terrazas17}.   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 (\citealt{Bell04, Faber07, Huang13, Muzzin13, Ilbert13, Tomczak14,Straatman15}; \citetalias{Fang18}). 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 $\sim$0.5 dex in \ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}\ from $z = 3$ to now in  Section~\ref{SubSec:Successes}. Even if this full correction were applied, we would still have a net BH growth factor across the GV of 1.0~ dex, or $\times10$.  We take this as our estimate of the  \emph{minimum} amount of BH growth in the GV in what follows.

We can now use this growth factor to estimate how much energy BHs emit during their stay in the GV.  Galaxies enter the green valley when their effective black hole energy equals $f_{\rm B}$ $\times$ their halo-gas binding energy, where $f_{\rm B}$ is set to 4 to match the zero points of the blue lines to the red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}.  We assume for the sake of argument that $\ensuremath{M_*}$ grows like $\Sig1$ through the GV (0.2~dex) and that $\ensuremath{M_{\rm vir}}$ grows similarly, galaxies being near the knee of the SMHM relation at this point, where the log slope is near unity  (Figure~\ref{Fig:HMSM}).  Then $\ensuremath{E_{\rm bind}} \sim \ensuremath{M_{\rm vir}}^{5/3}$ grows by $5/3 \times 0.2 = 0.3$~dex (Eq.~\ref{Eq:GasBindingEnergy}), and the ratio $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$ grows by at least $1.0-0.3 = 0.7$~dex, which is a factor of 5.  This together with $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}  = 4$ on entry to the green valley then implies that $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}  \gtrsim 20$ when galaxies come to a halt on the quenched side.  The factor $\times$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 $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm 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 \citet{Terrazas17} sample of black holes is rather small, highly biased to massive galaxies (as shown in Figure~\ref{Fig:ThreePanel}c), and reflects only BHs today.  We start by listing seven different pieces of evidence that do favor uniform and steep black hole growth in the green valley over extended cosmic times:

1) Rodgr\'iguez-Puebla et al.~(in prep.) have assembled a sample of reliable BH masses that is considerably larger than the sample in \citet{Terrazas17}.  Using $\Delta\log SSFR$ to locate a galaxy's position in the green valley, as predicted by Figure~\ref{Fig:Cartoon}, they derive a single, unified scaling law that fits all BH masses.  The coefficient of  $\Delta\log SSFR$ in the new law is consistent with the slope and length of the GV vectors in Figure~\ref{Fig:Cartoon}.

2) The zero points of the $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{\sigma_{\rm e}}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ relations for star-forming galaxies are roughly 1~dex below those of quenched galaxies (\citetalias{Kormendy13}; \citealp{Ho14,Reines15,Savorgnan16,Terrazas16,Davis18}).  These offsets make sense if it is assumed that black holes grow steeply in the green valley as galaxies pass between the two relations.

3) A long-standing puzzle has been the very different black hole masses of the Milky Way ($4\times10^6~\ensuremath{\rm M_{\odot}}$) and M31 ($1.5\times10^8~\ensuremath{\rm M_{\odot}}$) \citepalias{Kormendy13} despite the small difference in their stellar masses of only 0.2-0.3~dex \citep{Sick14,Licquia15}.  No standard BH scaling law has been able to fit both galaxies.  The two objects are plotted in Figure~\ref{Fig:Cartoon}.  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 \citetalias{Kormendy13} 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 $\ensuremath{M_*}$, AGN frequency rises strongly in the green valley.  This is true for SDSS galaxies \citep[e.g.,][]{Schawinski07,Schawinski10,Guo19} and also for CANDELS galaxies \citep[e.g.,][]{Nandra07,Gu18,Aird19}.

5) \citetalias{Rodriguez17} used abundance matching to measure evolutionary tracks for $\Sig1$ and $\ensuremath{M_{\rm BH}}$ back in time.  Rodr\'iguez et al.~(in prep.) now find, using similar methods, that black holes in quenched galaxies follow evolutionary tracks that look much like the broken power-laws in Figure~\ref{Fig:Cartoon}.

6) In Section~\ref{SubSec:PredictionsTests}, we show that steep BH growth through the green valley (like that in Figure~\ref{Fig:Cartoon}) 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 \citet{vandenBosch16} when fitting all galaxies.

7) Lastly, there is the theoretical justification from the hydro simulations cited above, which predict that black hole growth is unleashed when bulges reach a threshold size, resulting in steep growth trajectories similar to the green vectors in Figure~\ref{Fig:Cartoon}.

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 \citet{Aird19}, which show only a modest increase in X-ray luminosity as galaxies enter the GV ($\times$2-3).  If black holes grow by at least $\times$10 in the GV, they accrete 10\% of their mass on the main sequence and 90\% of their mass in the GV (assuming no further growth after quenching).  If main sequence and GV lifetimes are similar (which if anything probably overestimates time in the GV), average X-ray fluxes should be roughly an order of magnitude higher in the GV, not the factor of $2\text{-}3$ observed.  In short, AGNs should not only be prevalent in the GV, they should also be \emph{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 \citet{Aird19} has not been corrected for Compton-thick sources with $N_{\rm H} > 10^{24}$~cm$^{-2}$ \citep{Aird18}. Obscured AGN are reviewed by \citet{Hickox18}, who conclude that their number could be large.  The second option is that BH mass accretion in the GV is radiatively inefficient, i.e., that black holes can accrete much mass without shining in X-rays.  This is thought to occur in ADAFs \citep{Yuan14}, but their mass accretion rates are low whereas our green valley sources would be accreting fast.  The point is that a \emph{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 GV 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 black hole mass and halo mass that exists at the end of quenching, where $\ensuremath{E_{\rm BH}} \sim 20\times\ensuremath{E_{\rm 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.


\subsection{The Black-Hole/Dark Halo Connection}\label{SubSec:BHsandHalos}

The nature of the link between black holes and halos is controversial.  A tight correlation between $\ensuremath{M_{\rm BH}}$ and $\ensuremath{M_{\rm vir}}$ was claimed by \citet{Ferrarese02}, \citet{Baes03}, and \citet{Bandara09} but was disputed by \citetalias{Kormendy13} 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~\ref{Fig:Cartoon} clearly imply different correlations for different kinds of galaxies.  The equation $\Sig1 \sim \ensuremath{M_*}^{0.88}$ is the observed ridgeline slope of star-forming galaxies at fixed epoch in $\Sig1\text{-}\ensuremath{M_*}$ \citepalias{Barro17}.  Converting this to BH mass using $\ensuremath{M_{\rm BH}} \sim \Sig1^{1.76}$ yields $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_*}^{1.55}$.  This is steeper than the ridgeline for quenched galaxies ($\ensuremath{M_{\rm BH}} \sim \ensuremath{M_*}^{1.16}$), lies below it, and should show more scatter because star-forming galaxies have different  $\ensuremath{R_{\rm e}}/\ensuremath{\langle R_{\rm e}\rangle}$, and therefore different BH masses at fixed $\ensuremath{M_*}$.  These predictions agree tolerably well with observations (\citetalias{Kormendy13}; \citealp{Reines15,Savorgnan16,Terrazas16,Davis18,Delvecchio19}).  Finally, using $\ensuremath{M_*} \sim \ensuremath{M_{\rm vir}}^{1.75}$ for star-forming galaxies predicts $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_{\rm vir}}^{2.71}$, but with large scatter due to $\ensuremath{R_{\rm 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 $\ensuremath{M_*} = 10^{10.7}~\ensuremath{\rm M_{\odot}}$ (Figure~\ref{Fig:HMSM}).  Since most quenched ridgeline galaxies are more massive than this \citepalias{Kormendy13}, we start with the quenched ridgeline slope $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_*}^{1.16}$ and use the SMHM relation $\ensuremath{M_*} \sim \ensuremath{M_{\rm vir}}^{0.35}$ above the knee to obtain $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_{\rm 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 $\ensuremath{M_{\rm vir}} \sim 10^{12.5}~\ensuremath{\rm M_{\odot}}$.  The net result is that BH mass is predicted to rise strongly in star-forming galaxies with halo mass below $\ensuremath{M_{\rm vir}} = 10^{12.5}~\ensuremath{\rm 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, \citet{Ferrarese02} found a single tight relation for all galaxies with slope 1.65-1.82 by estimating $\ensuremath{M_{\rm vir}}$ from the outer rotation curves of spirals in the range $\log \ensuremath{M_{\rm vir}} = 11.5\text{-}12.5$ and the outer kinematic profiles of ellipticals \citep{Kronawitter00} in the range $\log \ensuremath{M_{\rm vir}} = 12.5\text{-}13.5$.  \citet{Baes03} found a single slope of 1.27 from a larger sample of spirals but the same ellipticals from \citet{Kronawitter00}.  Finally, \citet{Bandara09} found a slope of 1.55 using strong lensing to estimate halo masses of mainly ellipticals in the range $\log \ensuremath{M_{\rm vir}} = 13.0\text{-}13.8$.  None of these results agrees with the predictions, the slopes being too uniform and equal to some compromise value between the two predicted extremes of 0.40 to 2.71.

The steep slopes found for the massive end in these papers clearly disagree with elementary estimates.  Large spirals like M31 have $\ensuremath{M_{\rm BH}}$'s a few times $10^8~\ensuremath{\rm M_{\odot}}$ and dark halos a few times $10^{12}~\ensuremath{\rm M_{\odot}}$, while large clusters like Coma have central $\ensuremath{M_{\rm BH}}$'s a few times $10^9~\ensuremath{\rm M_{\odot}}$ but dark halos a few times $10^{15}~\ensuremath{\rm M_{\odot}}$.  The differential log slope is about 1/3, 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 \citetalias{Kormendy13} and \citet{Delvecchio19}.  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 $\ensuremath{M_{\rm vir}}$ are too small, and the calculated slopes of $\ensuremath{M_{\rm BH}}$ vs\@.~ $\ensuremath{M_{\rm vir}}$ are too steep.

To summarize, it is not true that there is \emph{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 \citetalias{Fang13}, who guessed that a close coordination between halo and BH is needed in order to produce a narrow quenched $\Sig1\text{-}\ensuremath{M_*}$ ridgeline but wondered how the halo ``knows'' what the BH is doing. The answer is that the halo \emph{makes} the black hole during the star-forming phase and then is quenched by it when $\ensuremath{E_{\rm BH}} \sim 20\times \ensuremath{E_{\rm bind}}$ (Section~\ref{SubSec:BoundaryCrossing}), closing the loop.  Black-hole mass differences due to $\Delta \log \ensuremath{R_{\rm e}}$ cause extra scatter about the $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_{\rm vir}}$ relation for star-forming and green valley galaxies but do not fundamentally change this basic picture.

\subsection{Black Hole Scaling Laws and Fundamental Planes}\label{SubSec:PredictionsTests}

This section returns to the problem mentioned in Section~\ref{SubSec:QuenchSlope} that the value of $v$ in the model seems too small.  Specifically, our value $v = 1.76$ implies that $\ensuremath{M_{\rm BH}} \sim \ensuremath{\sigma_{\rm e}}^{3.52\pm0.26}$ (if $\sigma_1 \sim \ensuremath{\sigma_{\rm e}}$) whereas the observed relation (for classical bulges and ellipticals in \citetalias{Kormendy13}) is $\ensuremath{M_{\rm BH}} \sim \ensuremath{\sigma_{\rm e}}^{4.4\pm0.3}$.  This suggests that $v$ is too small by the factor $4.4/3.52 \simeq 1.25 \pm 0.12$.  \citet{vandenBosch16F} fitted all galaxies together and found an even steeper law $\ensuremath{M_{\rm BH}} \sim \ensuremath{\sigma_{\rm e}}^{5.35\pm0.23}$, for a difference with us of $5.35/3.52 = 1.5\pm 0.13$.  We don't have measurements of $\ensuremath{\sigma_{\rm e}}$ in the CANDELS database, but we do have measurements of $\ensuremath{M_*}$ and $\ensuremath{R_{\rm e}}$, which were also fitted by \citet{vandenBosch16}.  He found $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_*}^{2.91\pm0.14}\ensuremath{R_{\rm e}}^{-2.77\pm0.22}$ from fitting all galaxies, which is close to $\ensuremath{M_{\rm BH}} \sim \Sig1^{2.7\pm0.17}$ using the $p$,$q$ mapping.  This implies that $v$ is too small by the factor $2.7/1.76 \simeq 1.5\pm0.14$, i.e., the same as $\ensuremath{\sigma_{\rm e}}$. (Note: this last also implies that $\ensuremath{M_*}/\ensuremath{R_{\rm e}} \sim \ensuremath{\sigma_{\rm e}}^2.0$.)  The goal of this section is to understand these slope differences.

We have taken advantage of these rough correspondences between $\Sig1$, $\sigma_1$, $\ensuremath{\sigma_{\rm e}}$, and $\ensuremath{M_*}/\ensuremath{R_{\rm e}}$ to label the X-axis in Figure~\ref{Fig:Cartoon} 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~\ref{Fig:Cartoon}), whereas \citet{vandenBosch16} 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~\ref{Fig:Cartoon}.  A rough estimate is shown by the black dashed line in Figure~\ref{Fig:Cartoon}, which is schematically fit to  all galaxies.  The slope of the black line is 2.6, very close to the \citet{vandenBosch16} value of 2.7 above.  The construction is crude but illustrates clearly how major BH mass growth through the GV, accompanied by little or no growth in $\Sig1$, can produce a big increase in the average slope of $\ensuremath{M_{\rm BH}}$ vs\@.~ $\Sig1$ (or $\sigma$, or $\ensuremath{M_*}/\ensuremath{R_{\rm e}}$) when all galaxies are fitted together.

The difference between our slope and that of \citetalias{Kormendy13} 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 \citetalias{Kormendy13} may actually be in the green valley, and at the high-mass end, \citet{Lauer07} have argued that $\ensuremath{M_{\rm BH}}\text{-}\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~\ref{SubSec:QuenchSlope} deliberately chose the  $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ relation, not $\ensuremath{M_{\rm BH}}\text{-}\sigma$, for quenched galaxies, precisely because it is less sensitive to merging.)

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

Separately, \citet{Hopkins07a} fitted a plane to quenched ridgeline galaxies only and found $\ensuremath{M_{\rm BH}} \sim \ensuremath{M_*}^{1.78 \pm 0.40} \ensuremath{R_{\rm e}}^{-1.05 \pm 0.37}$.  As mentioned in Section~\ref{SubSec:InputPowerLaws}, the  $p$, $q$ mapping for quenched galaxies is empirically $\Sig1 \sim \ensuremath{M_*}^{1.0} \ensuremath{R_{\rm e}}^{-0.5}$.  Our model then predicts $\ensuremath{M_{\rm BH}} \sim \Sig1^{1.76} \sim \ensuremath{M_*}^{1.76} \ensuremath{R_{\rm 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 \citet{vandenBosch16} for all galaxies, is really a stand-in for $\Sig1$ (or $\ensuremath{\sigma_{\rm e}}$).  Similar calculations (not shown) can explain the slopes of BHFPs fitted to $\ensuremath{\sigma_{\rm e}}\text{-}\ensuremath{R_{\rm e}}$ and $\ensuremath{\sigma_{\rm e}}\text{-}\ensuremath{M_*}$ by \citet{Hopkins07a} and by \citet{Aller07}.

Yet another expression offered for the quenched BHFP is $\ensuremath{M_{\rm BH}} \sim (\ensuremath{M_*} \ensuremath{\sigma_{\rm e}}^2)^{\alpha}$, where \citet{Hopkins07a,Hopkins07b} find $\alpha \simeq 0.7$ and \citet{Aller07} find $\alpha \simeq 0.6$. The quantity $\ensuremath{M_*} \ensuremath{\sigma_{\rm e}}^2$ is proportional to \emph{galaxy} binding energy per particle, which these authors proposed as the property that regulates BH mass.  This was motivated by the merger scenario of \citet{Hopkins07b}, in which gas that is brought to the center of the galaxy builds the black hole until the BH energy unbinds the surrounding gas.  In spirit, this is similar to the unbinding criterion in Eq.~\ref{Eq:EbhEbind}, but we use gas-binding energy to the halo, $\ensuremath{M_{\rm vir}}\ensuremath{V_{\rm vir}}^2 \sim \ensuremath{M_{\rm 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~\ref{Sec:ThermalModel} halo gas-binding energy $\ensuremath{E_{\rm bind}} \sim  \ensuremath{M_{\rm vir}}\ensuremath{V_{\rm vir}}^2 \sim \ensuremath{M_{\rm vir}}^{5/3}$.  Setting $\ensuremath{\sigma_{\rm e}}^2 \sim \ensuremath{M_*}/\ensuremath{R_{\rm e}}$, noting that $\ensuremath{R_{\rm e}} \sim \ensuremath{M_*}^{0.44}$ along the quenching boundary, and substituting $\ensuremath{M_*} \sim \ensuremath{M_{\rm vir}}^{1.75}$ on the SMHM relation yields $\ensuremath{M_{\rm vir}}^{5/3} \sim (\ensuremath{M_*} \ensuremath{\sigma_{\rm e}}^2)^{0.61} \sim \ensuremath{E_{\rm bind}}$.  In other words, the unexpected exponent with the value 0.6-0.7 in the model of \citet{Hopkins07b} and \citet{Aller07} 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~\ref{Fig:Cartoon} helps to explain this by positing different scaling laws (and $\ensuremath{R_{\rm 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 SSFR$ to reflect a galaxy's location in the green valley (see Item 1 in  Section~\ref{SubSec:BoundaryCrossing}).  Finally, it appears that planar fits to two variables are really stand-ins for fits to $\Sig1$ (or $\ensuremath{\sigma_{\rm e}}$), which are the fundamental quantities that describe $\ensuremath{M_{\rm BH}}$.



\subsection{Black Hole-Galaxy ``Co-Evolution''; the Origin and Meaning of BH Scaling Laws}\label{SubSec:BHGrowth}

Because of its highly transparent nature, the binding-energy model is a good tool for  probing the nature of so-called BH-galaxy ``co-evolution''.  In our opinion, this word has been used too loosely by astronomers simply to signify that the properties of two objects (e.g., black holes 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 co-evolution as ``the \emph{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 \citep{Dubois15, Habouzit17, Dekel19}, EAGLE \citep{Bower17,Davies19a,Oppenheimer19}, and FIRE \citep{Angles17}.  We do not include IllustrisTNG because it uses two different BH feedback modes (thermal-quasar mode and kinetic-radio mode;  \citet{Weinberger18,Habouzit19,Nelson19}), which would complicate the discussion.

Three different evolutionary phases are considered:

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

In the above-cited \emph{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 \citep{Dekel19,Bellovary19}, which may further reduce BH gas accretion.  BHs are typically started with seed masses of around $10^5~\ensuremath{\rm M_{\odot}}$, and growth beyond that in this phase is typically only $\times3$ even though galaxy masses increase by up to $\times100$.  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 black holes and galaxies do {not} co-evolve -- they are not even closely correlated.

2) \emph{Entry to the green valley; transition to rapid BH growth:} In the \emph{binding-energy model}, both of these phenomena occur together, when the BH energy output passes $4 \times \ensuremath{E_{\rm 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 vs. red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}). Galaxy, halo, and black hole properties are perfectly correlated at the transition, as in the previous star-forming phase, but now black holes are back-reacting on their halos and galaxies. This reciprocal influence satisfies the definition of co-evolution, and so black holes, galaxies, and their halos \emph{are co-evolving} as galaxies enter the green valley.

In the cited \emph{hydro models}, the transition to rapid BH growth occurs \emph{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 \citep{Dekel19}, winds no longer buoyant \citep{Bower17,McAlpine18}, major mergers \citep{McAlpine18}, deepening central potential well \citep{Angles17,Dekel17}, 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\text{-}10}~\ensuremath{\rm M_{\odot}}$, a halo mass near $10^{12}~\ensuremath{\rm M_{\odot}}$,  or a fixed halo virial temperature $\ensuremath{T_{\rm vir}}$ near $10^{6{\text-}7}$~K, all of which are plausibly related to the kinetic energy content of supernovae winds.  Since this is a fixed number \citep{Dekel86}, the theory tends to pick out a single stellar or halo mass where BH growth is launched, thus creating a momentary correspondence between the black hole, galaxy, and halo but not true co-evolution since there is no reciprocity a this point.

3) \emph{BH growth in the green valley; final transition to fully quenched:} In the \emph{binding-energy model}, the BH growth trajectory vs\@.~ $\Sig1$ in the green valley steepens.  The log slope of the blue track in Figure~\ref{Fig:Cartoon} is 1.76, whereas the log slope of the green GV vectors is 1.5 dex/0.2 dex = 7.5 (for these numbers, see Figure~\ref{Fig:Elbow}; 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 black hole, galaxy, and halo are continuing to ``talk'' to one another, and \emph{this co-evolution is plausibly what makes all green valley vectors parallel and of the same length}.  Ultimately, galaxies arrive on the fully quenched ridgeline when $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}} \simeq 20$, which enforces absolute correspondence between halos and black holes again at that point.  After that, halos continue to gain mass and drift away from both black holes and galaxies, and co-evolution ends (except possibly for ``maintenance mode'').

In \emph{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 \citep{Davies19a,Oppenheimer19}, in which fractional halo gas mass, $f_{\rm CGM}$, appears to be a barometer of the degree of quenching.  According to this parameter, entry into the green valley corresponds to $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}} \simeq  2$, and complete loss of halo gas corresponds to $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}} \simeq  20$ \citep[][Figure 2b]{Oppenheimer19}.  These numbers are not far from the binding-energy model, which gives 4 and 20 respectively.  On that basis, it is plausible that co-evolution 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 $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$ again tends to be a fixed factor at full quenching \citep{Terrazas19}.

Based on the above, we see that real co-evolution 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 (cf.~Figure~\ref{Fig:Cartoon}).  But the ultimate uniformity of quenched galaxies also depends on the uniformity of galaxies \emph{as they enter the green valley}, which is impressed in our model by the initial scaling law $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$ for star-forming galaxies.  Thus, in the binding-energy model, the final black-hole regularities of galaxies on the quenched ridgeline arise from the combination of a reciprocal, co-evolving phase in the green valley that builds on a uni-directional, causal phase for star-forming galaxies.

This insight motivates a new look, in our model, at the nature of the two $\ensuremath{M_{\rm BH}} \sim \Sig1$ scaling laws, for star-forming and quenched galaxies separately.  In the star-forming phase, $\Sig1$ determines $\ensuremath{M_{\rm BH}}$ single-handedly according to $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$, and the scaling relation is an \emph{evolutionary track} that all galaxies must follow (blue line in Figure~\ref{Fig:Cartoon}). Eventually, galaxies arrive on the quenched ridgeline, and a second scaling law for $\ensuremath{M_{\rm BH}}$ vs\@.~ $\Sig1$ emerges having the same slope as the blue line but higher zero point (red line in Figure~\ref{Fig:Cartoon}).  However, despite the mathematical similarity of the second law to the first, its interpretation is entirely different -- galaxies do not evolve along this second line, and it is not an evolutionary track.   Rather, it is the \emph{boundary} where $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}} \simeq 20$, and $\Sig1$ along this line \emph{labels the paths that galaxies took through the green valley to get there}.

In summary, some scaling laws may represent unique evolutionary tracks but others may be boundaries linking points on \emph{different} tracks where a third quantity takes on a common value.  The scaling law for quenched galaxies in our model is of the second sort.  It is not a track like the blue line, and galaxies do not evolve along it;  rather, it connects galaxies where $\ensuremath{E_{\rm BH}}/\ensuremath{E_{\rm bind}}$ has the common value $\simeq20$.

With that as background, we now explore briefly whether the blue line might itself be a boundary rather than an evolutionary track.  Up to now, we have regarded the blue line as a track that all star-forming galaxies must evolve along while building their black holes.  But actually, BHs need not follow this track for their entire lives -- they only need to obey the requisite relations \emph{as they enter the green valley}.  The blue line would then express the relationship between $\ensuremath{M_{\rm BH}}$ and $\Sig1$ that needs to obtain only there.

This change of perspective allows us to address a discrepancy that we have thus far glossed over, which is a disagreement between the binding-energy model and state-of-the-art theory.  As stressed in Section~\ref{SubSec:InputPowerLaws}, a generally accepted \emph{analytical} explanation of the BH scaling laws is thus far lacking, which means that the best theory comes from hydro models.  However, these models disagree with the binding-energy model in that their black holes in star-forming galaxies lie ``dormant'' and evolve rapidly only later, approaching the blue line \emph{from below} (see the dashed green lines in Figure~\ref{Fig:Cartoon}). Galaxies no longer evolve along the blue line -- instead, it becomes a \emph{boundary scaling law} linking together different entry points to the green valley along different evolutionary tracks.

This new interpretation has several implications.  First, it means that there is no relation like $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$ for star-forming galaxies or, if there is, it is much broader and fuzzier than what we've assumed.  The zero point of the blue line was chosen to fit nearby star-forming galaxies with \emph{detected} black holes.  The new interpretation says that many star-forming galaxies with the same $\ensuremath{M_*}$ and $\Sig1$ have even smaller black holes.  In other words, the locus of measured black holes in star-forming galaxies, as in \citet{Ho14} and \citet{Reines15}, is only the \emph{upper envelope} of a much larger population with smaller black holes below this locus.  This conclusion agrees with estimates that the measured BH relation is biased high for star-forming galaxies, having missed a population of small black holes \citep{Shankar16}.  It would also explain the large scatter in the $\ensuremath{M_{\rm BH}} \sim \ensuremath{\sigma_{\rm e}}^4$ law that is seen among star-forming  galaxies  (\citetalias{Kormendy13}, \citealt{Ho14}) and the anomalously small BHs found in some of them \citep[e.g., M33,][]{Gebhardt00}.

With the abandonment of the blue line as an evolutionary track, a new theory for the nature of the boundary to the green valley is needed.  One possibility is the model of AGN winds by \citet{King03} and \citet{King15}, which says that the structure of quasar outflows fundamentally changes at a certain black hole mass.  According to this theory, winds change from being momentum-driven to energy-driven when
\begin{equation}\label{Eq:KingPounds1}
\frac{M_{\rm BH}}{10^8~\rm{M_{\odot}}}
\simeq 3.2 \times \left(\frac{f_{\rm g}}{0.16}\right) \sigma_{200}^4.
\end{equation}
Here 0.16 is the cosmic baryon-to-dark matter ratio and $f_{\rm g}$ is the total gas mass fraction in the central galaxy potential well (which is assumed to be isothermal with velocity dispersion $\sigma_{200}$).  This change makes winds much more efficient at driving gas out of the galaxy and halo.  Combining this with Eq.~\ref{Eq:Sig1_sig1} from Appendix~\ref{App:A3} gives
\begin{equation}\label{Eq:KingPounds2}
\frac{M_{\rm BH}}{10^8~\rm{M_{\odot}}}
\simeq 0.09\times \left(\frac{f_g}{0.16}\right) \left(\frac{\Sig1}{10^9~\ensuremath{\rm M_{\odot}}~{\rm kpc^{-2}}}\right)^{2.0}.
\end{equation}
This can be compared to the equation of the blue line in Figure~\ref{Fig:Cartoon}, given by Eq.~\ref{Eq:Mbh_Sig1} in Appendix~\ref{App:A3}, which is:
\begin{equation}\label{Eq:KingPounds3}
\frac{M_{\rm BH}}{10^8~\rm{M_{\odot}}} = 0.025\times\left(\frac{\Sig1}{10^9~\ensuremath{\rm M_{\odot}}~{\rm kpc^{-2}}}\right)^{1.76}.
\end{equation}
If the small difference between the two exponents can be ignored, we see that the King-Pounds expression will lie on or near the blue line provided
\begin{equation}\label{Eq:KingPounds4}
\frac{f_g}{0.16} = \frac{0.025}{0.09} \sim 0.3,
\end{equation}
or
\begin{equation}\label{Eq:KingPounds5}
{f_g}\sim 0.04,
\end{equation}
which is a reasonable value for the gas mass fraction at the center of a galaxy that is building a black hole.

The King-Pounds theory therefore predicts a wind-enhancement boundary in the vicinity of the blue line.\footnote{This point was also made by \citet{King03} and \citet{King15}.}  Hydro models approaching this boundary from below would experience a large increase in wind efficiency at this point such that feedback energy from the BH is suddenly transferred more efficiently to the halo gas.  With higher feedback, quenching becomes perceptible, and galaxies enter the green valley along the steep green tracks in Figure~\ref{Fig:Cartoon}.

We therefore have an \emph{analytic} theory that predicts a boundary in the right place at the entrance to the green valley, together with a plausible way of feeding that boundary (from below) as in the hydro simulations.  However, a major ingredient is still missing: the theory predicts a boundary in $\ensuremath{M_{\rm BH}}\text{-}\Sig1$ only, whereas matching quenching boundaries are also seen in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$  and $\Sig1\text{-}\ensuremath{M_*}$ (cf.~Figure~\ref{Fig:ThreePanel}bc).\footnote{The boundary in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ in Figure~\ref{Fig:ThreePanel}a is not fundamental, being a mapping of the $\Sig1\text{-}\ensuremath{M_*}$ boundary into $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ via the $p$ ,$q$ relation.}  In other words, to match all the data, a galaxy needs to cross over the King-Pounds boundary in $\ensuremath{M_{\rm BH}}\text{-}\Sig1$ \emph{with the right $\ensuremath{M_*}$}.

In summary, the existence of a quenching boundary involving \emph{three} separate variables ($\ensuremath{M_{\rm BH}}$, $\Sig1$, and $\ensuremath{M_*}$) is a major constraint whose significance has not been appreciated until now.  The original binding-energy model met this constraint by making the $\ensuremath{M_{\rm BH}}\text{-}\Sig1$ relation a tight evolutionary track obeyed by \emph{all} star-forming galaxies, which removed one degree of freedom.  If the $\ensuremath{M_{\rm BH}}\text{-}\Sig1$ relation is instead a boundary that can be crossed in different places by galaxies with different $\ensuremath{M_*}$, a separate process is needed to align $\Sig1$ and $\ensuremath{M_{\rm BH}}$ with the proper value of $\ensuremath{M_*}$ at the crossing point.  Meeting this second constraint is a new test for hydro models.


\subsection{Comparison to the Model of Lilly et al. (2016)}\label{SubSec:L16}

We close this Discussion by returning to the model of Lilly et al.~(2016; L16), which was mentioned in the Introduction. This model is distinguished by the assumption that the quenching probability is a function of $\ensuremath{M_*}$ only, not of, say, $\ensuremath{M_*}$ and $\ensuremath{R_{\rm e}}$. Despite this, the model still manages to reproduce the sloping contours of equal quenched fraction in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ by exploiting two factors: 1) galaxies that quench earlier have smaller $\ensuremath{R_{\rm e}}$ and thus lower the quenched ridgeline zero point at fixed $\ensuremath{M_*}$ (progenitor bias), and 2) $\ensuremath{R_{\rm e}}$ values are light-weighted and shrink upon fading in this model owing to lower-$M/L$ stars in their outskirts (also mentioned by \citetalias{Fang13}). These two effects together increase the fraction of quenched galaxies at smaller $\ensuremath{R_{\rm e}}$ at fixed $\ensuremath{M_*}$.

Despite this success, the model fails to explain other data:

1) The difference in $\ensuremath{R_{\rm e}}$ between quenched vs\@.~ star-forming galaxies in the model is smaller at high redshift than at low redshift (see \citetalias{Lilly16}, Figure 3), whereas the reverse is true in real data: quenched galaxies, particularly lower-mass ones, are \emph{much} smaller than star-forming galaxies at $z \sim 2$ than at later times (cf.~Figure~\ref{Fig:Mstar_Sig1_SMA1}). This is naturally accommodated in the binding-energy model because $\ensuremath{R_{\rm e}}$ is free to shrink more at early times due to compaction (Section~\ref{SubSec:Compaction}), whereas the \citetalias{Lilly16} model does not have this extra degree of freedom.

2) The \citetalias{Lilly16} model predicts a $\times10$ radial increase in SSFR from 1 to 10~kpc within star-forming galaxies of mass $10^{10.75}~\ensuremath{\rm M_{\odot}}$ at all redshifts, which is generated by the assumption of strong inside-out star formation. A gradient of this magnitude is needed in order to produce the required amount of $\ensuremath{R_{\rm e}}$ shrinkage upon fading (see above). Current data on distant galaxies do not support gradients of this magnitude within 10~kpc -- the SSFR profiles of most galaxies are flat out to this radius \citep{Wang17, Tacchella18,Woo19}.

3) The \citetalias{Lilly16} model does not include BHs and therefore does not explain the $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ scaling relation of quenched galaxies. Quenching is assumed to occur with uniform probability per unit stellar mass growth, but no mechanism is offered to achieve this and the role of BHs in quenching is unaddressed.

4) Most important, the \citetalias{Lilly16} model has not yet attempted to match the $\Sig1\text{-}\ensuremath{M_*}$ diagram.  $\Sig1$ is likely to remain rather constant through and after quenching in the \citetalias{Lilly16} picture. But $\Sig1$ is set by the $p$, $q$ mapping for star-forming galaxies, and hence by $\ensuremath{R_{\rm e}}$ and $\ensuremath{M_*}$.  If quenching is equally probable for all galaxies at a given mass regardless of $\ensuremath{R_{\rm e}}$, the $\Sig1$ distributions of green valley and star-forming galaxies at a given $\ensuremath{M_*}$ should be the same.  Instead, the $\Sig1$ distributions of real green-valley galaxies are \emph{shifted to much higher values and are much narrower} than those of star-forming galaxies at all redshifts (\citetalias{Fang13,vanDokkum15,Barro17}; \citealp{Woo17,Lee18,Luo19}).  The key point is that GV galaxies are observed \emph{close to the beginning of quenching}, so that the evolutionary time difference between them and their star-forming progenitor galaxies should be small.  Even so, progenitor bias might induce a small offset in \emph{average} $\Sig1$, but it cannot explain the large difference in \emph{widths}.  Moreover, as a mass density, $\Sig1$ is not affected by fading. Thus, $\Sig1$ sidesteps both effects that the \citetalias{Lilly16} model uses to convert quenching contours that are intrinsically vertical at fixed mass into sloping contours in $\ensuremath{R_{\rm e}}$ vs\@.~ $\ensuremath{M_*}$.  The sloping behavior of $\Sig1$ vs\@.~ $\ensuremath{M_*}$ signals that a second physical variable is involved in quenching beyond $\ensuremath{M_*}$, which we have suggested is $\ensuremath{R_{\rm e}}/\ensuremath{\langle R_{\rm e}\rangle}$.

Parenthetically, we note that \emph{all} quenching models that rely largely on stellar (or halo) mass meet the same objection, of which the obvious examples are the pure ``halo-quenching'' models \citep[e.g.,][]{Ostriker77,Blumenthal84,Cattaneo06,Dekel06,Zu15,Zu18}.  Being one-dimensional theories in which the only variable is mass, they inherently fail to reproduce the observed diversity of properties of galaxies at a given stellar mass.

\section{Summary}\label{Sec:Summary}

A schematic model is presented for the growth of BHs in \emph{central} star-forming galaxies and their role in causing central galaxies to quench. The simplicity of the model is useful for appreciating how various factors in galaxy evolution -- such as structural evolution, black-hole growth rules, and quenching criteria -- interact and how results depend on these inputs.  The impetus for the new model is the belief that the observed scaling laws of galaxies contain important clues to quenching that can be elucidated from a proper model with the right degrees of freedom.  The paper starts with a list of 18 separate observations for quenched and star-forming galaxies that a successful model should explain.  The final binding-energy model presented here explains 11 of these items fairly well and sheds useful light on the remaining seven.

The main feature that enables this success is a new degree of freedom in which larger star-forming galaxies at fixed mass have smaller BHs owing to their lower central densities, $\Sig1$.  A second key assumption is that galaxies enter the green valley when the total effective energy emitted by their central black holes equals some fixed multiple, $f_{\rm B}$, of the binding energy of the gas in their dark halos.  Together, these assumptions imply that large-$\ensuremath{R_{\rm e}}$ star-forming galaxies must evolve to higher stellar masses before AGN feedback terminates the gaseous infall from their halos, and they accordingly have higher values of $\ensuremath{M_{\rm BH}}$, $\Sig1$, and $\ensuremath{R_{\rm e}}$ once quenched.  In other words, the simple addition of an $\ensuremath{R_{\rm e}}$ dependence to BH mass explains the otherwise mysterious \emph{sloped} quenching boundaries in $\Sig1\text{-}\ensuremath{M_*}$, $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$, and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ whereby more massive quenched galaxies have larger central densities, larger effective radii, and larger BHs. For star-forming galaxies, these three spaces are different projections of the same 4-D parent space, and mappings are presented to transform among them: $\Sig1 \sim \ensuremath{R_{\rm e}}^{-p}\ensuremath{M_*}^q$ is the transformation between $\Sig1$ and $\ensuremath{R_{\rm e}}$, where $p = 1$ and $q = 1.1$, while $\ensuremath{M_{\rm BH}} \sim \Sig1^v$ is the transformation between $\Sig1$ and $\ensuremath{M_{\rm BH}}$, where $v = 1.76$.

As long as galaxies are star-forming, they are assumed to evolve along log-linear parallel tracks in this 4-D space.  Each track is labeled with two numbers: the starting redshift, $z_{\rm s}$, when the galaxy passes $\ensuremath{M_*} = 10^9~\ensuremath{\rm M_{\odot}}$, and the factor by which its radius $\ensuremath{R_{\rm e}}$ differs from the mean $\ensuremath{\langle R_{\rm e}\rangle}$ at that mass and time.  These two numbers make central star-forming galaxies a two-parameter family, which can also be labeled by $\ensuremath{M_*}$ (or equivalently, $\ensuremath{M_{\rm vir}}$) and $\ensuremath{R_{\rm e}}/\ensuremath{\langle R_{\rm e}\rangle}$ today.  The physical origin of the scatter in $\ensuremath{R_{\rm e}}$ is not known, but halo concentration and/or halo formation time are promising candidates. If further tests support this, then $\ensuremath{M_{\rm vir}}$ and $c$ (or formation time) would be the two principal halo parameters that predict the evolutionary paths of central galaxies.

All star-forming galaxies in a given slice of $\ensuremath{M_*}$ are assumed to have the same halo mass, $\ensuremath{M_{\rm vir}}$, at all times, and the mean effective $\ensuremath{\langle R_{\rm e}\rangle}$ in that slice is always $0.02 \times \ensuremath{R_{\rm vir}}$. The star-forming branch of the stellar-mass-halo-mass relation (SMHM) is represented by a power-law $\ensuremath{M_*} = {\rm const.} \times \ensuremath{M_{\rm vir}}^s$ with zero scatter, where $s = 1.75$ and the zero point is constant with time. Both slope and zero point are taken from abundance matching \citepalias{Rodriguez17}.  Another key assumption is that all galaxies in a single mass slice make stars at the same rate regardless of $\ensuremath{R_{\rm e}}/\ensuremath{\langle R_{\rm e}\rangle}$, in agreement with observational data.   These assumptions, together with the physics of evolving halos, yield mean evolutionary track slopes and zero points for star-forming galaxies in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ and $\Sig1\text{-}\ensuremath{M_*}$ that generally agree with observations.

Galaxies are assumed to enter the green valley when they pass over a boundary in $\Sig1\text{-}\ensuremath{M_*}$ (and matching boundaries in $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ or  $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$). The empirically observed boundary, shown as the red lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}, is based on the well-measured quenched ridgeline in $\Sig1\text{-}\ensuremath{M_*}$, which is offset downwards by 0.2~dex in $\Sig1$ to represent entry into the green valley.  In the evolving empirical power-law model in Section~\ref{Sec:EvolPLModel}, this boundary is defined as the point where the total energy deposited in the halo by the black hole, $\eta \ensuremath{M_{\rm BH}} c^2$, equals some empirical energy quota needed to impact the halo and stop gas infall. To fit the observed quenching boundary slopes in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ forces the energy quota to have the form $\ensuremath{E_{\rm quench}} = {\rm const.}\times h(z)^{0.66v} \ensuremath{M_{\rm vir}}^t$, where $h(z)$ is the normalized Hubble constant $H(z)/H(0)$ and $t = 2.03$. This prescription fits the red lines precisely, and the log-linear slope of the boundary in $\Sig1\text{-}\ensuremath{M_*}$ is $t/vs = 0.66$, in $\ensuremath{R_{\rm e}}-\ensuremath{M_*}$  it is $q-pt/vs = 0.44$, and in $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ it is $t/s = 1.16$.

The above empirical expression for $\ensuremath{E_{\rm quench}}$ is then compared to $4\times$ the total binding-energy content of the hot halo gas, \ensuremath{E_{\rm bind}}, and it is shown that their behaviors vs\@.~ mass and time are similar (red lines vs\@.~ blue lines in Figure~\ref{Fig:Mstar_Sig1_SMA1}).  Hence, nearly the same quenching boundary is predicted using this more physically-based expression for halo gas-binding energy.

Four variations on this basic picture are considered in Section~\ref{Sec:Variants}, including two that express the black hole-halo balance in terms of an energy-\emph{rate} balance instead of total-energy equality. Such rate prescriptions have been employed in SAM models in the past but fail badly (as implemented here) to match the observed width of the green valley, the sloping power-law nature of the quenching boundaries, or the steady flow of galaxies from the blue cloud to the red sequence since $z \sim 3$.

The formula for BH mass growth, $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$ with constant zeropoint, applies only to galaxies on the star-forming main sequence (i.e., down to $-0.45$~dex below the ridgeline). Evidence is given that the BH growth relation steepens as galaxies enter the green valley, indicating an acceleration in BH growth relative to $\Sig1$ at that point (Figure~\ref{Fig:Cartoon}).  Increased BH accretion in the green valley is also seen in X-ray luminosities and in hydro models.  Total BH mass growth while galaxies are in the green valley is at least $\times$10.  The star-forming slope, $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$, is reacquired when galaxies come to rest on the quenched ridgeline, but now with a different zero point that is $\sim$10 times larger than the star-forming zero point.  This scenario for BH growth is summarized in Figure~\ref{Fig:Cartoon}.

Rapid black hole growth in the green valley implies that AGN activity should be high there.  AGNs are indeed frequent in the green valley, but integrated X-ray luminosity is not as high as expected.  Perhaps AGNs are obscured, or accretion may be radiatively inefficient.  This conflict with the X-ray data is the biggest unanswered challenge for the binding-energy model.

A surprise from the model is that the derived galaxy evolutionary tracks in $\Sig1\text{-}\ensuremath{M_*}$ and $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ are nearly parallel to the quenching boundary and most galaxies would not cross it if it were stationary.  In practice, it is the decline in halo gas-binding energy (at fixed $\ensuremath{M_{\rm vir}}$) that brings the quenching boundary down and across the tracks.  The angle between the tracks and boundary remains acute, however, and the small value of this angle is responsible for the 1$^+$-dex width of the green valley.

An important requirement at the heart of any successful model is that the quenching boundary and the star-forming galaxy cloud must maintain a constant relationship to one another over time in $\Sig1\text{-}\ensuremath{M_*}$, i.e., their zero points must evolve together in lockstep \citepalias{Barro17}. This constancy is observed and is responsible for the unchanged appearance of the $\Sig1\text{-}\ensuremath{M_*}$ diagram over time (except for the downward evolution in both zero points) and is also why the green valley quenching channel has nearly the same mass at all redshifts.  An exact match between the boundary and the star-forming cloud is imposed by fiat in the empirically derived quenching boundary, $\ensuremath{E_{\rm quench}}$, in Eq.~\ref{Eq:EbhEhalo}.  The more physically motivated quenching boundary based on halo gas-binding energy, $\ensuremath{E_{\rm bind}}$ in Eq.~\ref{Eq:EbhEbind}, also moves down but not quite fast enough, and the average green-valley mass would vary slightly on that account from $-0.2$~dex at $z = 3$ to $+0.2$ dex at $z=0.5$.

For most of the paper, we assume that star-forming galaxies build their black holes by evolving \emph{along} the relation $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$ (blue line in Figure~\ref{Fig:Cartoon}).  If this is true, then BH masses in star-forming galaxies can be read off directly from the $\Sig1\text{-}\ensuremath{M_*}$ diagram using $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$, or from the $\ensuremath{R_{\rm e}}\text{-}\ensuremath{M_*}$ diagram using $\ensuremath{M_{\rm BH}} = \rm{const.}\times\ensuremath{R_{\rm e}}^{-1.76}\ensuremath{M_*}^{1.94}$.

However, inspired by hydro models, Section~\ref{SubSec:BHGrowth} explores an alternate picture in which black holes in star-forming galaxies spend most of their lives \emph{below} the blue line and grow rapidly only later. The wind theory of \citet{King03} and \citet{King15} predicts that AGN winds should suddenly become more efficient as galaxies cross a relation very near the blue line, setting the stage for entry into the green valley.  The blue line would then no longer be an evolutionary track but rather the boundary where the King-Pounds relation is satisfied.  This picture  is still incomplete, however, because galaxies must cross the boundary with the right combination of $\Sig1$, $\ensuremath{M_{\rm BH}}$, \emph{and $\ensuremath{M_*}$}, and the tight connection to $\ensuremath{M_*}$ is still theoretically unexplained.

The major takeaway messages from the binding-energy model are:

1) Central galaxies perceptibly quench when the total effective energy emitted by their black holes equals some constant multiple, $f_{\rm B}$, of the binding energy of gas in their dark halos.  Galaxies enter the green valley when $f_{\rm B}$ is near 4 (Section~\ref{Sec:ThermalModel}) and are fully quenched when $f_{\rm B}$ is near 20 (Section~\ref{SubSec:BoundaryCrossing}).

2) Effective radius is an important second parameter in shaping the life histories of galaxies and the masses of their black holes.  Promising candidates for the origin of radius scatter (at fixed stellar mass) are halo concentration and/or halo formation time (Section~\ref{SubSec:RadiusScatter}).

3) In the simple version of the theory, star-forming galaxies build their black holes \emph{along} the relation $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$ (blue line in Figure~\ref{Fig:Cartoon}).  The positive dependence on central stellar density means that black holes are smaller in galaxies with larger $\ensuremath{R_{\rm e}}$ at fixed mass (Section~\ref{SubSec:InputPowerLaws}).  Points (1) and (2) then explain the sloping boundaries between quenched and star-forming galaxies in key structural diagrams (Section~\ref{SubSec:EvoTracksQuenchBound}). In this interpretation, the blue line is an evolutionary track.

4) An alternate version of the theory says that black holes in star-forming galaxies spend most of their time \emph{below} the relation $\ensuremath{M_{\rm BH}} = \rm{const.}\times\Sig1^{1.76}$ and grow only later when SNae winds fail, as in hydro models.  The King-Pounds wind theory then predicts a jump in wind feedback efficiency when black holes reach the blue line, which would propel host galaxies into the green valley.  In this interpretation, the blue line is a boundary not an evolutionary track, and many small black holes should exist below the blue line in star-forming galaxies that have not yet been detected (Section~\ref{SubSec:BHGrowth}).

5) The quenching boundary moves with time in structural diagrams because the gas binding energies of halos at fixed $\ensuremath{M_{\rm vir}}$ are lower at late times (Section~\ref{SubSec:Successes}). The evolutionary tracks of star-forming galaxies are intrinsically parallel to the boundary, and it is mainly the sideways motion of the boundary that causes them to quench (Sections~\ref{SubSec:EvoTracksQuenchBound} and \ref{Sec:ThermalModel}).

6) The downward motion of the quenching boundary in $\Sig1\text{-}\ensuremath{M_*}$ is matched by a similar motion of the star-forming cloud. This holds the quenching and star-forming populations in the same relationship to one another, which keeps the green valley mass nearly constant. In physical terms, black hole masses are getting smaller with time owing to the growth in galaxy radii (at fixed $\ensuremath{M_*}$) while dark halos are getting easier to quench with time owing to the decline in halo-gas binding energy (at fixed $\ensuremath{M_{\rm vir}}$) (Section~\ref{SubSec:Successes}).  This approximate equality is why quenching has looked so similar since $z \sim 3$ (Section~\ref{SubSec:Successes}).

7) Black-hole scaling laws and their relationships to dark halos are well defined but are different for star-forming galaxies, the green valley, and quenched galaxies (Section~\ref{SubSec:PredictionsTests}). At least 90\% of black hole mass is gained in the green valley (Sections~\ref{SubSec:InputPowerLaws} and  \ref{SubSec:BoundaryCrossing}).

8) The quenching process is highly uniform from halo to halo.  The rms raw scatter in the factor $f_{\rm B}$ is at most 0.28~dex based on the observed scatter in the $\Sig1\text{-}\ensuremath{M_*}$ quenched ridgeline and will probably be considerably less than this once progenitor bias is fully allowed for (Section~\ref{Sec:ThermalModel}).

9) Despite the successes of the binding-energy model, the physical basis of the model remains incomplete.  The foundation of the model is the way that black holes grow in star-forming galaxies, and the default version assumes that a tight evolutionary track, $\ensuremath{M_{\rm BH}} = \rm{const.} \times \Sig1^{1.76}$, is followed by all  star-forming galaxies.  This assumption has empirical support but disagrees with recent hydro models, which predict more scatter and generally smaller black holes in star-forming galaxies.  An alternative theory predicts that AGN wind efficiency should increase strongly with $\Sig1$ in star-forming galaxies, but this creates a boundary in $\ensuremath{M_{\rm BH}}\text{-}\Sig1$, not an evolutionary track.  It explains the boundary in $\ensuremath{M_{\rm BH}}\text{-}\Sig1$ but not the associated boundaries in $\Sig1\text{-}\ensuremath{M_*}$ or $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ since $\ensuremath{M_*}$ is not linked to these other two quantities.  The existence of \emph{three separate, matched quenching boundaries} -- in $\Sig1\text{-}\ensuremath{M_{\rm BH}}$, $\Sig1\text{-}\ensuremath{M_*}$, and $\ensuremath{M_{\rm BH}}\text{-}\ensuremath{M_*}$ -- is  currently the biggest challenge to theoretical quenching models.


\begin{acknowledgments}
The major impetus for the model in this paper was the structural and star-forming data on galaxies from the CANDELS program on the \textit{Hubble Space Telescope}.  We thank the many workers who helped to create this unique database.  The main CANDELS observations were supported under program number HST-GO-12060, provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract NAS5-26555.  C. Z. acknowledges support from NFSC grants 11403016,11433003,11673018, Innovation Program 2019-01-07-00-02-E00032 of SMEC. C. Z. thanks the support by Science and Technology Commission of Shanghai Municipality to the Key Lab for Astrophysics, International Joint Lab and Project No. 18590780100.  Members of the CANDELS team at UCSC acknowledge support from NASA \textit{HST} grant GO-12060.10-A and from NSF grants AST-0808133 and AST-1615730.  We thank the anonymous referee for raising important questions that stimulated several significant and beneficial improvements.
\end{acknowledgments}

