∎

11institutetext:  Jiulin Wang 22institutetext: Mengmeng Song 33institutetext:
Yong Xia 44institutetext: LMIB of the Ministry of Education, School of
Mathematical Sciences, Beihang University, Beijing, 100191, P. R. China
44email<EMAIL_ADDRESS>(J. L. Wang<EMAIL_ADDRESS>(M. M.
Song<EMAIL_ADDRESS>(Y. Xia, corresponding author) )

# Trust-region and $p$-regularized subproblems: local nonglobal minimum is the
second smallest objective function value among all first-order stationary
points ††thanks: This research was supported by the Beijing Natural Science
Foundation, grant Z180005 and by the National Natural Science Foundation of
China under grants 11822103, 11571029, and 11771056.

Jiulin Wang  Mengmeng Song  Yong Xia

(Received: date / Accepted: date)

###### Abstract

The local nonglobal minimizer of trust-region subproblem, if it exists, is
shown to have the second smallest objective function value among all KKT
points. This new property is extended to $p$-regularized subproblem. As a
corollary, we show for the first time that finding the local nonglobal
minimizer of Nesterov-Polyak subproblem corresponds to a generalized
eigenvalue problem.

###### Keywords:

Trust-region subproblem $p$-regularized subproblem Nesterov-Polyak subproblem
Local nonglobal minimizer Monotonicity property

###### MSC:

90C20, 90C26, 90C30

††journal: JOTA

##  1 Introduction

The trust-region subproblem (TRS) is to minimize a (possibly nonconvex)
quadratic function over a Euclidean ball. It is an essential problem in each
iteration of the trust region method for solving nonlinear programming
problems Y15 . (TRS) has other applications in convex quadratic integer
programming B15 , bounded linear regression, best rank-$1$ tensor
approximation, and image deconvolution P20 .

Nonconvex (TRS) has the property of zero Lagrangian-duality gap, see the
recent survey on hidden convexity X20 . There is a necessary and sufficient
optimality condition for the global minimizer, established in the early 1980s
G81 ; S82 ; M83 . In 1994, Martínez M94  proved that (TRS) has at most one
local nonglobal minimizer. In 2020, Wang and Xia WX20  established the
necessary and sufficient optimality condition for the local nonglobal
minimizer.

Besides the second-order optimality conditions for local and global
minimizers, first-order optimality (KKT) conditions are deeply investigated
for (TRS), see for example LPR98 . Gander G80  showed that the objective
function values of KKT points of linear least squares problem with a quadratic
equality constraint are monotonic in terms of their Lagrangian multipliers. In
order to locate local solutions of the Celis-Dennis-Tapia (CDT) subproblem C85
which minimizes a quadratic function over the intersection of two ellipsoids,
Chen and Yuan CY00  studied the monotonicity property for KKT points of (CDT)
in a partial order relation of corresponding Lagrangian multipliers.

The $p$-regularized subproblem ($p$-RS) is an unconstrained optimization
problem of minimizing a (possibly nonconvex) quadratic function with an
additional $p$-th power regularization term of the norm of the variables GRT10
; Hsia17 . In literature, if the regularization term is cubic, ($p$-RS) is
known as Nesterov-Polyak subproblem N06 . For more references, we refer to W07
; C11 ; CGT11 ; L20  and references therein. If the regularization term is
quartic, ($p$-RS) reduces to the double-well potential optimization, which has
particular applications in solid mechanics and quantum mechanics F17 ; X17 .

In this paper, based on the monotonicity properties for first-order stationary
points, we prove that the local nonglobal minimizer of (TRS), if it exists,
has the second smallest objective function value among all KKT points. We then
extend the new property to ($p$-RS) and show that the local nonglobal
minimizer of (TRS) has the smallest objective function value among all
critical points. As a corollary, we show for the first time that the local
nonglobal minimizer of the Nesterov-Polyak subproblem could be founded by
solving a generalized eigenvalue problem.

The rest of this paper is organized as follows. In Section 2, we review some
properties for local minimizers of (TRS) and ($p$-RS). In Section 3 and 4, we
first show that monotonicity properties hold for first-order stationary points
of (TRS) and ($p$-RS), respectively. As the main result of this paper, we
prove that the local nonglobal minimums of (TRS) and ($p$-RS) are the second
smallest objective function value among all first-order stationary points,
respectively. Conclusion and discussions are made in Section 5.

##  2 Preliminaries

In this section, we review some characterizations for local minimizers of
(TRS) and ($p$-RS).

###  2.1 Characterizations for local minimizers of (TRS)

In this subsection, we briefly review some properties of local minimizers of
nonconvex (TRS):

$\displaystyle(\rm TRS)~{}~{}$ $\displaystyle\min$ $\displaystyle
f(x)=\frac{1}{2}x^{T}Qx+c^{T}x$ $\displaystyle{\rm s.t.}$ $\displaystyle
x^{T}x-1\leq 0,$

where $Q\in\mathbb{R}^{n\times n}$ is not positive semidefinite and
$c\in\mathbb{R}^{n}$. Up to an eigenvalue decomposition, we can always assume
that $Q$ is a diagonal matrix with $n$ diagonal elements being
$(0>)~{}\alpha_{1}\leq\alpha_{2}\leq\cdots\leq\alpha_{n}$.

$(x;\lambda)\in\mathbb{R}^{n}\times\mathbb{R}$ is a KKT point of (TRS) if and
only if it satisfies the following KKT system:

$\displaystyle Qx+c=-\lambda x,$ (2.1) $\displaystyle\lambda(x^{T}x-1)=0,$
(2.2) $\displaystyle x^{T}x-1\leq 0,~{}~{}\lambda\geq 0.$ (2.3)

In case of $\lambda>0$, (2.2)-(2.3) is equivalent to

$x^{T}x-1=0.$ (2.4)

Suppose $-\lambda\not\in\\{\alpha_{1},\alpha_{2},\cdots,\alpha_{n}\\}$. Then
by substituting the $x$-solution of (2.1) into (2.4), we obtain that $\lambda$
must be a zero point of the following so-called secular function:

$\varphi(\lambda)=\sum_{i=1}^{n}\frac{c_{i}^{2}}{(\alpha_{i}+\lambda)^{2}}-1.$
(2.5)

The first and second derivatives of $\varphi(\lambda)$ can be written as
follows:

$\displaystyle\varphi^{\prime}(\lambda)=-2\sum_{i=1}^{n}\frac{c_{i}^{2}}{(\alpha_{i}+\lambda)^{3}},$
$\displaystyle\varphi^{\prime\prime}(\lambda)=6\sum_{i=1}^{n}\frac{c_{i}^{2}}{(\alpha_{i}+\lambda)^{4}}.$

If $\varphi(\lambda)$ has zero points, then $c\neq 0$ and hence
$\varphi^{\prime\prime}(\lambda)>0$ on each nonempty interval
$(-\alpha_{i+1},-\alpha_{i})$ for $i\in\\{1,2,\cdots,n-1\\}$. It follows that,
on each interval, $\varphi(\lambda)$ is strongly convex and thus has at most
two zero points.

In the early 1980s, necessary and sufficient conditions for global minimizers
of (TRS) are established by Gay G81 , Sorensen S82 , Moré and Sorensen M83 .

######  Lemma 2.1 (G81 , Theorem 2.1)

$x^{*}$ is a global minimizer of (TRS) if and only if ${x^{*}}^{T}x^{*}=1$ and
there is a $\lambda^{*}\in\mathbb{R}$ such that (2.1) holds and

$\lambda^{*}\geq-\alpha_{1}(>0).$

In 1994, Martínez M94  identified two cases where the local nonglobal
minimizer of (TRS) does not exist.

######  Lemma 2.2 (M94 , Lemmas 3.2, 3.3)

Suppose either $\alpha_{1}=\alpha_{2}$ or $c_{1}=0$, there is no local
nonglobal minimizer of (TRS).

In the same paper, Martínez M94  proved that (TRS) has at most one local
nonglobal minimizer based on the following detailed characterization.

######  Lemma 2.3 (M94 , Theorem 3.1(i))

If $\underline{x}$ is a local nonglobal minimizer of (TRS), then there is a
nonnegative $\underline{\lambda}\in(-\alpha_{2},-\alpha_{1})$ such that (2.1)
holds, and $\underline{\lambda}$ is a zero point of $\varphi(\lambda)$ and

$\varphi^{\prime}(\underline{\lambda})\geq 0.$

Recently, Lemma 2.3 has been updated by Wang and Xia WX20 .

######  Lemma 2.4 (WX20 , Theorem 3.1)

$\underline{x}$ is a local nonglobal minimizer of (TRS), if and only if there
is a nonnegative $\underline{\lambda}\in(-\alpha_{2},-\alpha_{1})$ such that
(2.1) holds, and $\underline{\lambda}$ is a zero point of $\varphi(\lambda)$
and

$\varphi^{\prime}(\underline{\lambda})>0.$

In 1998, Lucidi et al. LPR98  proved that the strict complementarity condition
holds at the local nonglobal minimizer of (TRS).

######  Lemma 2.5 (LPR98 , Proposition 3.5)

If $\underline{x}$ is a local nonglobal minimizer of (TRS), the corresponding
Lagrangian multiplier is $\underline{\lambda}$, then (2.2) holds with

$\underline{x}^{T}\underline{x}-1=0,~{}~{}\underline{\lambda}>0.$

In 2017, Adachi et al. A17  observed the following result.

######  Lemma 2.6 (A17 , Lemma 3.1)

If $(x;\lambda)$ is a KKT point of (TRS) satisfying (2.1) and (2.4), then we
have

${\rm det}~{}M_{1}(\lambda)=0,$

where

$\displaystyle M_{1}(\lambda):=\left[{\begin{array}[]{*{20}{c}}-I&Q\\\
Q&-cc^{T}\end{array}}\right]-\lambda\left[{\begin{array}[]{*{20}{c}}0&-I\\\
-I&0\end{array}}\right].$ (2.10)

That is, $\lambda$ is a real generalized eigenvalue of $M_{1}(\lambda)$.

Let $(x;\lambda)$ satisfying (2.1). We can rewrite $M_{1}(\lambda)$ as

$M_{1}(\lambda)=\left[{\begin{array}[]{*{20}{c}}I&{}\hfil\\\ {}\hfil&Q+\lambda
I\end{array}}\right]\left[{\begin{array}[]{*{20}{c}}-I&I\\\
I&-xx^{T}\end{array}}\right]\left[{\begin{array}[]{*{20}{c}}I&{}\hfil\\\
{}\hfil&Q+\lambda I\end{array}}\right],$

Then,

$\displaystyle{\rm det}~{}M_{1}(\lambda)=(-1)^{n}({\rm det}(Q+\lambda
I))^{2}(1-x^{T}x).$ (2.11)

The following result follows from Lemma 2.6 and (2.11).

######  Corollary 2.1 (S17 , Corollary 3)

The set of real generalized eigenvalues of $M_{1}(\lambda)$ is nonempty.
Moreover, if ${\rm det}~{}M_{1}(\lambda)=0$, then either $-\lambda$ is an
eigenvalue of $Q$ or $\lambda$ is a zero point of $\varphi(\lambda)$.

###  2.2 Characterizations for local minimizers of ($p$-RS)

In this subsection, we review some characterizations for local minimizers of
($p$-RS):

$\displaystyle(p\rm-RS)~{}~{}$
$\displaystyle\mathop{\min}\limits_{x\in\mathbb{R}^{n}}$
$\displaystyle\left\\{g(x)=\frac{1}{2}x^{T}Qx+c^{T}x+\frac{\sigma}{p}\|x\|^{p}\right\\},$

where $\sigma>0$, $p>2$, $\|\cdot\|$ is the Euclidean norm, and $Q,c$ have the
same definitions as in (TRS). The global minimum of ($p$-RS) is attainable,
even if $Q$ is not positive semidefinite. Actually, $g(x)$ is coercive, i.e.,

$\lim_{\|x\|\rightarrow+\infty}g(x)=+\infty,$

since $p>2$ and then

$\lim_{\|x\|\rightarrow+\infty}\frac{\frac{1}{2}x^{T}Qx+c^{T}x}{\|x\|^{p}}=0.$

We call $x\in\mathbb{R}^{n}$ a critical point of ($p$-RS) if and only if it
satisfies the following first order necessary condition for local minimizers
of ($p$-RS):

$\displaystyle\nabla g(x)=(Q+\sigma\|x\|^{p-2}I)x+c=0,$ (2.12)

where $\nabla g$ denotes the gradient of $g(\cdot)$. Let

$\displaystyle t=\|x\|^{p-2},$ (2.13)

then (2.12) can be rewritten as

$\displaystyle(Q+\sigma tI)x+c=0.$ (2.14)

Suppose $-\sigma t\not\in\\{\alpha_{1},\alpha_{2},\cdots,\alpha_{n}\\}$. Then
by substituting the $x$-solution of (2.14) into (2.13), we obtain that $t$
must be a zero point of the following secular function:

$h(t)=\sum_{i=1}^{n}\frac{c_{i}^{2}}{(\sigma
t+\alpha_{i})^{2}}-t^{\frac{2}{p-2}}.$ (2.15)

In the following, we make an assumption that $x\neq 0$, then it holds that
$c\neq 0$ from (2.12). Under this assumption, $h(t)$ has the same zeros as

$p(t)={\rm log}\left(\sum_{i=1}^{n}\frac{c_{i}^{2}}{(\sigma
t+\alpha_{i})^{2}}\right)-\frac{2}{p-2}{\rm log}(t).$

Hsia et al. Hsia17  proved that $p^{\prime\prime}(t)>0$ on each nonempty
interval
$\left(-\frac{\alpha_{i+1}}{\sigma},-\frac{\alpha_{i}}{\sigma}\right)$ for
$i\in\\{1,2,\cdots,n-1\\}$.

The first derivatives of $h(t)$ and $p(t)$ can be written as follows:

$\displaystyle h^{\prime}(t)=-\sum_{i=1}^{n}\frac{2\sigma c_{i}^{2}}{(\sigma
t+\alpha_{i})^{3}}-\frac{2}{(p-2)t}\cdot t^{\frac{2}{p-2}},$ (2.16)
$\displaystyle p^{\prime}(t)=-\frac{\sum_{i=1}^{n}\frac{2\sigma
c_{i}^{2}}{(\sigma t+\alpha_{i})^{3}}}{\sum_{i=1}^{n}\frac{c_{i}^{2}}{(\sigma
t+\alpha_{i})^{2}}}-\frac{2}{(p-2)t}.$ (2.17)

If $\hat{t}$ is a zero point of $h(t)$, according to (2.15), (2.16) and
(2.17), we obtain

$\displaystyle
p^{\prime}(\hat{t})=h^{\prime}(\hat{t})\cdot\hat{t}^{{\frac{-2}{p-2}}}.$
(2.18)

The following detailed characterizations of local minimizers of ($p$-RS) are
all due to Hsia et al. Hsia17 .

######  Lemma 2.7 (Hsia17 , Theorem 2.2)

$x^{*}$ is a global minimizer of ($p$-RS) for $p>2$ if and only if it is a
critical point satisfying (2.12) and

$\sigma\|x^{*}\|^{p-2}+\alpha_{1}\geq 0.$

Moreover, the $\ell_{2}$ norms of all the global minimizers are equal.

######  Lemma 2.8 (Hsia17 , Lemma 3.1)

Suppose $\underline{x}$ is a local non-global minimizer of ($p$-RS) for $p>2$.
It holds that $\underline{x}_{1}\neq 0$, $\alpha_{1}<\alpha_{2}$, and

$\sigma\|\underline{x}\|^{p-2}+\alpha_{2}>0.$

######  Lemma 2.9 (Hsia17 , Theorem 3.2)

The point $\underline{x}$ is a local nonglobal minimizer of ($p$-RS) for $p>2$
if and only if $\underline{x}$ is a critical point satisfying (2.12) and
$\underline{t}$ (=$\|\underline{x}\|^{p-2}$) is a root of

$h(t)=0,~{}t\in\left({\rm
max}\left\\{-\frac{\alpha_{2}}{\sigma},0\right\\},-\frac{\alpha_{1}}{\sigma}\right)$

such that $h^{\prime}(\underline{t})>0$.

Motivated by the observation that the Lagrange multiplier of each KKT point of
(TRS) satisfying (2.1) and (2.4) is a generalized eigenvalue of
$M_{1}(\lambda)$ A17 , Lieder L20  pointed out the following result holds for
the Nesterov-Polyak subproblem (i.e., ($p$-RS) for $p=3$).

######  Lemma 2.10 (L20 , Lemma 3.1)

If $x$ is a critical point of the cubic regularization, then $t~{}(=\|x\|)$ is
a real generalized eigenvalue of

$\displaystyle M_{2}(t)=\left[{\begin{array}[]{*{20}{c}}0&0&0&c^{T}\\\
0&-\sigma I&0&Q\\\ 0&0&-\sigma&0\\\
c&Q&0&0\end{array}}\right]-t\left[{\begin{array}[]{*{20}{c}}0&0&-\sigma&0\\\
0&0&0&-\sigma I\\\ -\sigma&0&0&0\\\ 0&-\sigma I&0&0\end{array}}\right].$
(2.27)

That is, $t$ is a real root of ${\rm det}~{}M_{2}(t)=0$.

Let $(x;t)$ satisfying (2.14). We can reformulate $M_{2}(t)$ as

$M_{2}(t)=\left[{\begin{array}[]{*{20}{c}}I&{}\hfil&{}\hfil&{}\hfil\\\
{}\hfil&I&{}\hfil&{}\hfil\\\ {}\hfil&{}\hfil&1&{}\hfil\\\
{}\hfil&{}\hfil&{}\hfil&Q+\sigma
tI\end{array}}\right]\left[{\begin{array}[]{*{20}{c}}0&0&\sigma t&-x^{T}\\\
0&-\sigma I&0&I\\\ \sigma t&0&-\sigma&0\\\
-x&I&0&0\end{array}}\right]\left[{\begin{array}[]{*{20}{c}}I&{}\hfil&{}\hfil&{}\hfil\\\
{}\hfil&I&{}\hfil&{}\hfil\\\ {}\hfil&{}\hfil&1&{}\hfil\\\
{}\hfil&{}\hfil&{}\hfil&Q+\sigma tI\end{array}}\right],$

Then, by Schur complement,

$\displaystyle{\rm det}~{}M_{2}(t)=(-1)^{n+1}\sigma^{2}\cdot({\rm
det}(Q+\sigma tI))^{2}(t^{2}-x^{T}x).$ (2.28)

The following result follows from Lemma 2.10 and (2.28).

###### Corollary 2.2

The set of real generalized eigenvalues of $M_{2}(t)$ is nonempty. Moreover,
if ${\rm det}~{}M_{2}(t)=0$, then one of the following cases happens: (1)
$t=0$; (2) $-\sigma t$ is an eigenvalue of $Q$; (3) $t$ is a zero point of
$p(t)$.

##  3 On local nonglobal minimum of (TRS)

We first prove a monotonicity property for KKT points of (TRS).

###### Lemma 3.1

Let $(x_{1};\lambda_{1})$ and $(x_{2};\lambda_{2})$ be two KKT points for
(TRS). Then

$\lambda_{2}\geq(>)\lambda_{1}~{}\Longrightarrow~{}f(x_{2})\leq(<)f(x_{1}).$

###### Proof

According to (2.1) and the definitions of $(x_{1};\lambda_{1})$ and
$(x_{2};\lambda_{2})$, we have

$\displaystyle Qx_{1}+c=-\lambda_{1}x_{1},~{}Qx_{2}+c=-\lambda_{2}x_{2}.$
(3.1)

By multiplying $x_{1}^{T}$ and $x_{2}^{T}$ to both sides of the two equalities
in (3.1), respectively, we have

$\displaystyle
x_{1}^{T}Qx_{1}+x_{1}^{T}c=-\lambda_{1}x_{1}^{T}x_{1},~{}~{}x_{1}^{T}Qx_{2}+x_{1}^{T}c=-\lambda_{2}x_{1}^{T}x_{2},$
$\displaystyle
x_{2}^{T}Qx_{2}+x_{2}^{T}c=-\lambda_{2}x_{2}^{T}x_{2},~{}~{}x_{2}^{T}Qx_{1}+x_{2}^{T}c=-\lambda_{1}x_{2}^{T}x_{1}.$

Then, it holds that

$\displaystyle 2(f(x_{2})-f(x_{1}))$ $\displaystyle=$
$\displaystyle(x_{2}^{T}Qx_{2}+x_{2}^{T}c)-(x_{1}^{T}Qx_{1}+x_{1}^{T}c)+x_{2}^{T}c-x_{1}^{T}c$
(3.2) $\displaystyle=$
$\displaystyle(-\lambda_{2}x_{2}^{T}x_{2})-(-\lambda_{1}x_{1}^{T}x_{1})+(-\lambda_{1}x_{2}^{T}x_{1})-(-\lambda_{2}x_{1}^{T}x_{2})$
$\displaystyle=$
$\displaystyle(\lambda_{2}-\lambda_{1})x_{1}^{T}x_{2}+\lambda_{1}x_{1}^{T}x_{1}-\lambda_{2}x_{2}^{T}x_{2}$
$\displaystyle\leq$ $\displaystyle(\lambda_{2}-\lambda_{1})\cdot
1+\lambda_{1}x_{1}^{T}x_{1}-\lambda_{2}x_{2}^{T}x_{2}$ (3.3) $\displaystyle=$
$\displaystyle\lambda_{1}(x_{1}^{T}x_{1}-1)-\lambda_{2}(x_{2}^{T}x_{2}-1)=0,$
(3.4)

where (3.3) follows from Cauchy-Schwarz inequality and (2.3), more precisely,

$x_{1}^{T}x_{2}\leq\|x_{1}\|\cdot\|x_{2}\|\leq 1,$ (3.5)

the last equality in (3.4) holds due to the complementarity condition (2.2).

The equality $f(x_{2})=f(x_{1})$ holds true if and only if the inequality in
(3.3) holds as an equality, that is

$(\lambda_{2}-\lambda_{1})(x_{1}^{T}x_{2}-1)=0.$ (3.6)

In case of $\lambda_{2}>\lambda_{1}$, (3.6) is equivalent to

$x_{1}^{T}x_{2}=1.$ (3.7)

It implies from (3.5) and (3.7) that

$x_{1}=x_{2}\neq 0.$ (3.8)

Combining the assumption $\lambda_{2}>\lambda_{1}$ with (2.1) and (3.8) yields
a contradiction. Therefore, $f(x_{2})=f(x_{1})$ if and only if
$\lambda_{2}=\lambda_{1}$. The proof is complete.

###### Remark 3.1

Let $(x_{1};\lambda_{1})$ and $(x_{2};\lambda_{2})$ be two KKT points for
(TRS) satisfying $\|x_{1}\|=\|x_{2}\|$, since

$x_{1}^{T}x_{1}-x_{1}^{T}x_{2}=x_{2}^{T}x_{2}-x_{1}^{T}x_{2}=\frac{1}{2}\left\|x_{1}-x_{2}\right\|^{2},$

it follows from (3.2) that

$f(x_{2})-f(x_{1})=\frac{\lambda_{1}-\lambda_{2}}{4}\|x_{1}-x_{2}\|^{2},$

which is presented in G80  for the relaxed least square problem.

As a corollary of Lemma 3.1, the following result holds.

###### Corollary 3.1

If $(x^{*};\lambda^{*})$ is KKT point of (TRS), $x^{*}$ is a global minimizer,
then $\lambda^{*}$ is the largest Lagrangian multiplier among all KKT points.
Moreover, all the global minimizers of (TRS) share the same Lagrangian
multiplier $\lambda^{*}$.

Besides, according to Lemma 3.1, Lemma 2.1, Lemma 2.6 and Corollary 2.1, we
can recover the following result for the global minimizer of (TRS).

######  Corollary 3.2 (A17 , Theorem 3.2)

If $x^{*}$ is a global minimizer of (TRS), $\lambda^{*}$ is the corresponding
Lagrangian multiplier, then $\lambda^{*}$ is equal to the largest real
generalized eigenvalue of $M_{1}(\lambda)$.

We now present the main result of this section.

###### Theorem 3.1

Let $(x_{1};\lambda_{1})$ and $(x_{2};\lambda_{2})$ be any two different KKT
points for (TRS). If $x_{1}$ is a local nonglobal minimizer of (TRS) and
$x_{2}$ is not a global minimizer of (TRS), then $f(x_{1})<f(x_{2})$.

###### Proof

According to Lemma 3.1, it is sufficient to prove $\lambda_{2}<\lambda_{1}$.

If $x_{2}^{T}x_{2}<1$, then (2.2) implies that $\lambda_{2}=0$. According to
Lemma 2.5, it holds that $\lambda_{1}>0$. Then, we obtain

$\lambda_{2}=0<\lambda_{1}.$

If $x_{2}^{T}x_{2}=1$, then either
$-\lambda_{2}\in\\{\alpha_{1},\alpha_{2},\cdots,\alpha_{n}\\}$ or
$\lambda_{2}$ is a zero point of $\varphi(\lambda)$ (2.5). Since $x_{2}$ is
not a global minimizer of (TRS), it follows from Lemma 2.1 that

$\lambda_{2}<-\alpha_{1}.$

According to Lemma 2.3, $\lambda_{1}$ is the largest zero point of
$\varphi(\lambda)$ in $(-\alpha_{2},-\alpha_{1})$. Then we have

$\lambda_{2}<\lambda_{1}<-\alpha_{1}.$

The proof is complete.

According to Lemma 2.2, if there is a local nonglobal minimizer of (TRS), then
we have $\alpha_{1}\neq\alpha_{2}$ and $c_{1}\neq 0$. Then, ${\rm
det}~{}M_{1}(-\alpha_{1})\neq 0$. Otherwise, it is not difficult to verify
that $c_{1}=0$, which leads to a contradiction. Combining ${\rm
det}~{}M_{1}(-\alpha_{1})\neq 0$, the proof of Theorem 3.1, Lemma 2.1, Lemma
2.6 and Corollary 2.1 together, we can recover the following result.

######  Corollary 3.3 (S17 , Theorem 5)

If $\underline{x}$ is a local nonglobal minimizer of (TRS),
$\underline{\lambda}$ is the corresponding Lagrangian multiplier, then
$\underline{\lambda}$ is equal to the second largest real generalized
eigenvalue of $M_{1}(\lambda)$.

##  4 On local nonglobal minimum of ($p$-RS)

We first prove a monotonicity property for critical points of ($p$-RS).

###### Lemma 4.1

Let $x_{1}$ and $x_{2}$ be two critical points of ($p$-RS) for $p>2$. Then

$\|x_{2}\|\geq(>)\|x_{1}\|~{}\Longrightarrow~{}g(x_{2})\leq(<)g(x_{1}).$

###### Proof

According to (2.12) and the definitions of $x_{1}$ and $x_{2}$, we have

$\displaystyle(Q+\sigma\|x_{1}\|^{p-2}I)x_{1}+c=0,~{}~{}(Q+\sigma\|x_{2}\|^{p-2}I)x_{2}+c=0.$
(4.1)

By multiplying $x_{1}^{T}$ and $x_{2}^{T}$ to both sides of the two equalities
in (4.1), respectively, we have

$\displaystyle
x_{1}^{T}Qx_{1}+x_{1}^{T}c+\sigma\|x_{1}\|^{p}=0,~{}~{}x_{1}^{T}Qx_{2}+x_{1}^{T}c+\sigma\|x_{2}\|^{p-2}x_{1}^{T}x_{2}=0,$
$\displaystyle
x_{2}^{T}Qx_{2}+x_{2}^{T}c+\sigma\|x_{2}\|^{p}=0,~{}~{}x_{2}^{T}Qx_{1}+x_{2}^{T}c+\sigma\|x_{1}\|^{p-2}x_{1}^{T}x_{2}=0.$

Then, $g(x_{2})-g(x_{1})$ is equal to

$\displaystyle\frac{1}{2}((x_{2}^{T}Qx_{2}+x_{2}^{T}c)-(x_{1}^{T}Qx_{1}+x_{1}^{T}c))+\frac{1}{2}(x_{2}^{T}c-x_{1}^{T}c)+\frac{\sigma}{p}(\|x_{2}\|^{p}-\|x_{1}\|^{p})$
(4.2) $\displaystyle=$
$\displaystyle(\frac{\sigma}{p}-\frac{\sigma}{2})(\|x_{2}\|^{p}-\|x_{1}\|^{p})+\frac{\sigma}{2}(\|x_{2}\|^{p-2}-\|x_{1}\|^{p-2})x_{1}^{T}x_{2}.$

Then by (4.2), we have

$\displaystyle\|x_{2}\|=\|x_{1}\|~{}\Longrightarrow~{}g(x_{2})=g(x_{1}).$
(4.3)

Next, we show that

$\displaystyle\|x_{2}\|>\|x_{1}\|~{}\Longrightarrow~{}g(x_{2})<g(x_{1}).$
(4.4)

If $x_{1}=0$ and $x_{2}\neq 0$, it follows from (4.2) and $p>2$ that

$\displaystyle g(x_{2})-g(x_{1})=$
$\displaystyle(\frac{\sigma}{p}-\frac{\sigma}{2})\|x_{2}\|^{p}<0.$ (4.5)

If $x_{1}\neq 0$, let $\|x_{2}\|=\alpha\|x_{1}\|$, where $\alpha\geq 1$.
Substituting $\|x_{2}\|=\alpha\|x_{1}\|$ to (4.2) yields that

$\displaystyle g(x_{2})-g(x_{1})$ $\displaystyle=$
$\displaystyle(\frac{\sigma}{p}-\frac{\sigma}{2})(\alpha^{p}-1)\|x_{1}\|^{p}+\frac{\sigma}{2}(\alpha^{p-2}-1)\|x_{1}\|^{p-2}\cdot
x_{1}^{T}x_{2}$ (4.6) $\displaystyle\leq$
$\displaystyle\frac{\sigma}{2}\|x_{1}\|^{p}\cdot\left(\frac{2-p}{p}(\alpha^{p}-1)+\alpha^{p-1}-\alpha\right).$

where (4.6) follows from Cauchy-Schwarz inequality, i.e.,

$x_{1}^{T}x_{2}\leq\|x_{1}\|\cdot\|x_{2}\|=\alpha\|x_{1}\|^{2}.$

Define

$\phi(\alpha)=\frac{2-p}{p}(\alpha^{p}-1)+\alpha^{p-1}-\alpha.$

The first and second derivatives of $\phi(\alpha)$ can be written as follows:

$\displaystyle\phi^{\prime}(\alpha)$ $\displaystyle=$
$\displaystyle(p-1)\alpha^{p-2}-(p-2)\alpha^{p-1}-1,$
$\displaystyle\phi^{\prime\prime}(\alpha)$ $\displaystyle=$
$\displaystyle(p-1)(p-2)(1-\alpha)\alpha^{p-3}.$ (4.7)

If $\alpha=1$, we can verify that

$\displaystyle\phi(1)=\phi^{\prime}(1)=\phi^{\prime\prime}(1)=0.$ (4.8)

If $\alpha>1$, according to (4.7) and (4.8), we have

$\displaystyle\phi^{\prime\prime}(\alpha)<\phi^{\prime\prime}(1)=0~{}\Longrightarrow~{}\phi^{\prime}(\alpha)<\phi^{\prime}(1)=0\Longrightarrow~{}\phi(\alpha)<\phi(1)=0.$
(4.9)

Then (4.4) holds due to (4.5) and (4.9). Combining (4.3) with (4.4), the proof
is complete.

If $p=3$, according to Lemma 4.1, Lemma 2.7, Lemma 2.10 and Corollary 2.2, we
recover the following result.

######  Corollary 4.1 (L20 , Corollary 3.2)

If $x^{*}$ is a global minimizer of ($p$-RS) for $p=3$, then $t^{*}=\|x^{*}\|$
is equal to the largest real generalized eigenvalue of $M_{2}(t)$.

Now we present the main result of this section.

###### Theorem 4.1

Let $x_{1}$ and $x_{2}$ be any two different critical points of ($p$-RS) for
$p>2$. If $x_{1}$ is a local nonglobal minimizer of ($p$-RS) and $x_{2}$ is
not a global minimizer, then $g(x_{1})<g(x_{2})$.

###### Proof

Let $t_{i}=\|x_{i}\|^{p-2}$ for $i=1,2$. According to Lemma 4.1, it is
sufficient to prove $\|x_{2}\|<\|x_{1}\|$ or equivalently, $t_{2}<t_{1}$.

We first prove that $c\neq 0$. As $x_{1}$ is a local nonglobal minimizer of
($p$-RS), according to Lemma 2.7, we have

$\displaystyle\sigma\|x_{1}\|^{p-2}+\alpha_{1}<0.$ (4.10)

Then, it follows from (4.10), (2.12) and Lemma 2.8 that $c_{1}\neq 0$, and
hence $c\neq 0$.

Next, as $c\neq 0$, (2.12) implies that $x_{1},x_{2}\neq 0$. Then, $p(t)$ is
well defined. We can conclude that either $-\sigma
t_{2}\in\\{\alpha_{1},\alpha_{2},\cdots,\alpha_{n}\\}$ or $t_{2}$ is a zero
point of $p(t)$. Since $x_{2}$ is not a global minimizer of ($p$-RS), it
follows from Lemma 2.7 that

$t_{2}<-\frac{\alpha_{1}}{\sigma}.$

Since $p(t)$ is strongly convex for $t\in\left({\rm
max}\left\\{-\frac{\alpha_{2}}{\sigma},0\right\\},-\frac{\alpha_{1}}{\sigma}\right)$,
it has at most two zero points in this interval. According to Lemma 2.9 and
(2.18), we have

$p^{\prime}(t_{1})>0.$

Thus $t_{1}$ is the largest zero point of $p(t)$ in the interval $\left({\rm
max}\left\\{-\frac{\alpha_{2}}{\sigma},0\right\\},-\frac{\alpha_{1}}{\sigma}\right)$.
Then we have

$t_{2}<t_{1}<-\frac{\alpha_{1}}{\sigma}.$

The proof is complete.

In Section 3, we have shown that if there is a local nonglobal minimizer of
(TRS), then ${\rm det}~{}M_{1}(-\alpha_{1})\neq 0$. Similarly, we can show
that if there is a local nonglobal minimizer of ($p$-RS) for $p=3$, then ${\rm
det}~{}M_{2}\left(-\frac{\alpha_{1}}{\sigma}\right)\neq 0$. Combing this
result with the proof of Theorem 4.1, Lemmas 2.7, 2.10 and Corollary 2.2, we
give the following result.

###### Corollary 4.2

If $\underline{x}$ is the local nonglobal minimizer of the Nesterov-Polyak
subproblem (i.e., ($p$-RS) for $p=3$), then $\underline{t}=\|\underline{x}\|$
is equal to the second largest real generalized eigenvalue of $M_{2}(t)$.

##  5 Conclusion and discussions

This paper shows that there are monotonicity properties for first-order
stationary points of trust-region subproblem (TRS) and $p$-regularized
subproblem ($p$-RS). Based on these properties, we point out that all the
global minimizers of (TRS) share the same Lagrangian multiplier and recover
the result that finding a global minimizer of (TRS) (or ($p$-RS)) corresponds
to a generalized eigenvalue problem. As the main contribution, we prove that
the local nonglobal minimizer of (TRS) (or ($p$-RS)), if it exists, has the
second smallest objective function value among all KKT points (critical
points). For the Nesterov-Polyak subproblem (i.e., ($p$-RS) with $p=3$), we
show for the first time that the local nonglobal minimizer, if it exists,
could also be founded by solving a generalized eigenvalue problem.

Our main results for ($p$-RS) or (TRS) may fail to hold for more general
optimization problems. Consider the problem of minimizing a univariate sextic
function:

$s(x)=\frac{1}{6}x^{6}-\frac{21}{10}x^{5}+\frac{57}{8}x^{4}-\frac{1}{12}x^{3}-\frac{435}{32}x^{2}-\frac{297}{32}x.$

There are four critical points, denoted by $x_{i}$, $i=1,2,3,4$. As
demonstrated in Fig. 1, $x_{1}=-0.5$ is neither a local minimizer nor a local
maximizer, $x_{2}=1.5$ is a global minimizer, $x_{3}=4.5$ is a local
maximizer, and $x_{4}=5.5$ is a local nonglobal minimizer. It is observed from
Fig. 1 that, the unique local nonglobal minimum is not the second smallest
objective function value among all critical points.

Figure 1: Variation of $s(x)$ with four critical points satisfying
$s(x_{2})<s(x_{1})<s(x_{4})<s(x_{3})$.

Moreover, the above unconstrained minimization problem can be reformulated as
the following nonconvex quadratic optimization with a quintic constraint:

$\displaystyle\mathop{\min}\limits_{x,y\in\mathbb{R}}$ $\displaystyle xy$
$\displaystyle{\rm s.t.}$ $\displaystyle
y=\frac{1}{6}x^{5}-\frac{21}{10}x^{4}+\frac{57}{8}x^{3}-\frac{1}{12}x^{2}-\frac{435}{32}x-\frac{297}{32}.$

One can verify that the constrained problem has four KKT points, corresponding
the four critical points of $s(x)$.

It seems that there is still a little room for possible extension. We conclude
this paper with the following conjecture:

###### Conjecture 5.1

For ${\rm(CDT)}$ subproblem or nonconvex quadratic optimization with a quartic
constraint, the smallest local nonglobal minimum has the second smallest
objective function value among all KKT points.

## References

  * (1) Y. Yuan, Recent advances in trust region algorithms, Math. Program., 151(1) (2015), pp. 249–281.
  * (2) C. Buchheim, R. Hübner, A. Schöbel, Ellipsoid bounds for convex quadratic integer programming, SIAM J. Optim., 25(2) (2015), pp. 741–769.
  * (3) A. H. Phan, M. Yamagishi, D. Mandic, A. Cichocki, Quadratic programming over ellipsoids with applications to constrained linear regression and tensor decomposition, Neural Comput. Appl., 32(2020), pp. 7097–7120.
  * (4) Y. Xia, A survey of hidden convex optimization, J. Oper. Res. Soc. China., 8(1) (2020), pp. 1–28.
  * (5) D. M. Gay, Computing optimal locally constrained steps, SIAM J. Sci. Stat. Comput., 2(2) (1981), pp. 186–197.
  * (6) D. C. Sorensen, Newton’s method with a model trust region modification, SIAM J. Numer. Anal., 19(2) (1982), pp. 409–426.
  * (7) J. J. Moré, D. C. Sorensen, Computing a trust region step, SIAM J. Sci. Statist. Comput., 4(3) (1983), pp. 553–572.
  * (8) J. M. Martínez, Local minimizers of quadratic functions on Euclidean balls and spheres, SIAM J. Optim., 4(1) (1994), pp. 159–176.
  * (9) J. Wang, Y. Xia, Closing the gap between necessary and sufficient conditions for local nonglobal minimizer of trust region subproblem, SIAM J. Optim., 30(3) (2020), pp. 1980–1995.
  * (10) S. Lucidi, L. Palagi, M. Roma, On some properties of quadratic programs with a convex quadratic constraint, SIAM J. Optim., 8(1) (1998), pp. 105–122.
  * (11) W. Gander, Least squares with a quadratic constraint, Numer. Math., 36(3) (1980), pp. 291–307.
  * (12) M. R. Celis, J. E. Dennis, R. A. Tapia, A trust region strategy for nonlinear inequality constrained optimization, in Numerical Optimization, R.T. Boggs, R.H. Byrd, and R.B. Schnabel, eds. (1984), SIAM, Philadelphia, PA, (1985), pp. 71–82.
  * (13) X. Chen, Y. Yuan, On local solutions of the Celis-Dennis-Tapia subproblem, SIAM J. Optim., 10(2) (2000), pp. 359–383.
  * (14) N. I. M. Gould, D. P. Robinson, H. Sue Thorne, On solving trust-region and other regularised subproblems in optimization, Math. Program. Comput., 2(1) (2010), pp. 21–57.
  * (15) Y. Hsia, R. L. Sheu, Y. Yuan, Theory and application of $p$-regularized subproblems for $p>2$, Optim. Method Softw., 32(5) (2017), pp. 1059-1077.
  * (16) Y. Nesterov, B. T. Polyak, Cubic regularization of Newton method and its global performance, Math. Program., 108(1) (2006), pp. 177–205.
  * (17) M. Weiser, P. Deuflhard, B. Erdmann, Affine conjugate adaptive Newton methods for nonlinear elastomechanics, Optim. Methods Softw., 22(3) (2007), pp. 413–431.
  * (18) C. Cartis, N. I. M. Gould, Ph. L. Toint, Adaptive cubic regularisation methods for unconstrained optimization. Part I: motivation, convergence and numerical results, Math. Program., 127(2) (2011), pp. 245–295.
  * (19) C. Cartis, N. I. M. Gould, Ph. L. Toint, Adaptive cubic regularisation methods for unconstrained optimization. Part II: worst-case function-and derivative-evaluation complexity, Math. Program., 130(2) (2011), pp. 295–319.
  * (20) F. Lieder, Solving large-Scale cubic regularization by a generalized eigenvalue problem, SIAM J. Optim., 30(4) 2020, pp. 3345–3358.
  * (21) S. C. Fang, D. Gao, G. X. Lin, R. L. Sheu, W. Xing, , Double well potential function and its optimization in the n-dimensional real space: part I, J. Ind. Manag. Optim., 13(3) (2017), pp. 1291–1305.
  * (22) Y. Xia, R. L. Sheu, S. C. Fang, W. Xing, Double well potential function and its optimization in the n-dimensional real space: part II, J. Ind. Manag. Optim., 13(3) (2017), pp. 1307–1328.
  * (23) S. Adachi, S. Iwata, Y. Nakatsukasa, A. Takeda, Solving the trust-region subproblem by a generalized eigenvalue problem, SIAM J. Optim., 27(1) (2017), pp. 269–291.
  * (24) M. Salahi, A. Taati, H. Wolkowicz, Local nonglobal minima for solving large-scale extended trust-region subproblems, Comput. Optim. Appl., 66(2) (2017), pp. 223–244.


# Constraints on cosmological coupling from the accretion history of
supermassive black holes

Mark Lacy National Radio Astronomy Observatory, Charlottesville, VA, USA
Athena Engholm Department of Physics and Astronomy, University of Hawai‘i at
Mãnoa, 2505 Correa Rd., Honolulu, HI 96822, USA  Institute for Astronomy,
University of Hawai‘i, 2680 Woodlawn Dr., Honolulu, HI 96822, USA  Duncan
Farrah Department of Physics and Astronomy, University of Hawai‘i at Mãnoa,
2505 Correa Rd., Honolulu, HI 96822, USA  Institute for Astronomy, University
of Hawai‘i, 2680 Woodlawn Dr., Honolulu, HI 96822, USA  Kiana Ejercito
Department of Physics and Astronomy, University of Hawai‘i at Mãnoa, 2505
Correa Rd., Honolulu, HI 96822, USA  Institute for Astronomy, University of
Hawai‘i, 2680 Woodlawn Dr., Honolulu, HI 96822, USA

###### Abstract

Coupling of black hole mass to the cosmic expansion has been suggested as a
possible path to understanding the dark energy content of the Universe. We
test this hypothesis by comparing the supermassive black hole (SMBH) mass
density at $z=0$ to the total mass accreted in AGN since $z=6$, to constrain
how much of the SMBH mass density can arise from cosmologically-coupled
growth, as opposed to growth by accretion. Using an estimate of the local SMBH
mass density of $\approx 1.0\times 10^{6}\,$M${}_{\odot}\,$Mpc-1, a radiative
accretion efficiency, $\eta$: $0.05<\eta<0.3$, and the observed AGN luminosity
density at $z\approx 4$, we constrain the value of the coupling constant
between the scale size of the Universe and the black hole mass, $k$, to lie in
the range $0<k\stackrel{{\scriptstyle<}}{{{}_{\sim}}}2$, below the value of
$k=3$ needed for black holes to be the source term for dark energy. Initial
estimates of the gravitational wave background using pulsar timing arrays,
however, favor a higher SMBH mass density at $z=0$. We show that if we adopt
such a mass density at $z=0$ of $\approx 7.4\times
10^{6}\,$M${}_{\odot}\,$Mpc-1, this makes $k=3$ viable even for low radiative
efficiencies, and may exclude non-zero cosmological coupling. We conclude
that, although current estimates of the SMBH mass density based on the black
hole mass – bulge mass relation probably exclude $k=3$, the possibility
remains open that, if the GWB is due to SMBH mergers, $k>2$ is preferred.

Supermassive black holes (1663): Luminosity function (942): Active Galactic
Nuclei (16)

##  1 Introduction

###  1.1 The Soltan Argument

There exists an elegant argument for constraining the growth of supermassive
black hole (SMBH) mass in galaxies via accretion. Soltan (1982) showed that
the integral of the AGN luminosity density over cosmic history gives an
estimate of the total accreted mass onto SMBHs, which can be compared to the
local mass density of SMBHs to constrain how they grow. This relation, between
two independently derived quantities, can give excellent constraints on both
AGN physics and the nature of SMBHs themselves. For example, using then
available data, Soltan (1982) showed that the history of accreted mass is
close to the local SMBH mass density, if a radiative efficiency of
$\eta\approx 0.1$ is assumed. This range in $\eta$ is in good agreement with
other estimates. Theoretical arguments give a range of
$0.05\lesssim\eta\lesssim 0.43$, depending on the spin of the black hole (BH)
(Bardeen, 1970; Thorne, 1974). This range is itself in good agreement with
most other studies, which find values in the range $0.1\lesssim\eta\lesssim
0.3$ (e.g. Yu & Tremaine, 2002; Elvis et al., 2002; Hopkins et al., 2007;
Shankar et al., 2009; Lacy et al., 2015; Zhang & Lu, 2017; Shankar et al.,
2020; Shen et al., 2020; Farrah et al., 2022). This consistency argues for the
increase in SMBH mass density with cosmic time to be driven by accretion.

###  1.2 Cosmologically-coupled black holes

For several decades it has been recognized that BHs with non-singular
interiors (including ‘dark energy’-like interiors) may have exteriors
identical to BHs with interior singularities (e.g. Gliner, 1966; Dymnikova,
1992; Faraoni & Jacques, 2007; Beltracchi & Gondolo, 2019; Cadoni et al.,
2023a). Recently, however, it was realized that the masses of non-singular BHs
could be affected by the cosmic expansion (Croker & Weiner, 2019). Indeed, in
a recent paper, Cadoni et al. (2023b) show that such a coupling is necessary
in General Relativity (GR) if BHs do not contain singularities. Thus the
detection of a non-zero cosmological coupling would have important
implications for our understanding of BHs.

The mass coupling between BHs and the cosmic expansion can be parameterized
via:

$M_{\bullet}(a)=M_{\bullet}(a_{i})\left(\frac{a}{a_{i}}\right)^{k}$ (1)

(Croker et al., 2020, 2021) in which $M_{\bullet}(a)$ and $M_{\bullet}(a_{i})$
are the BH mass at some later and earlier times respectively, $a$ and $a_{i}$
are the scale factors at those times ($a=(1+z)^{-1}$ at redshift $z$), and $k$
is the cosmological coupling strength, with $-3\leq k\leq 3$ from a causality
constraint.

Croker & Weiner (2019) show that the value of $k$ is related to the equation
of state of the medium internal to the BH. Expressing the equation of state as
$p=w\rho$, where $p$ is the pressure and $\rho$ the density, Croker et al.
(2021) show that $k=-3w$. Thus, for vacuum energy where $p=-\rho$, $k=3$. In
this case, the mass density of BHs is conserved as the Universe expands: their
number density decreases with $a^{-3}$ but their mass increases as $a^{3}$ and
so BHs could account for the cosmological constant (Croker & Weiner, 2019;
Farrah et al., 2023a) (see also Prescod-Weinstein et al., 2009, for an
alternative way BHs could be the origin of dark energy). In contrast, GR
solutions for non-singular BHs prefer $k=1$ (Cadoni et al., 2023b). In this
case, the BH mass growth would not compensate for cosmological dilution of
their density, so they would not contribute to dark energy.

Observationally, the evidence is mixed. Farrah et al. (2023b, a) find, based
on an analysis of SMBH versus host galaxy masses out to $z\approx 2$, that a
value of $k\approx 3$ is suggested and $k=0$ can be ruled out at $\approx
3\sigma$. On the other hand, Lei et al. (2023) argue, based on a much smaller
sample of AGN at $z>4.5$ with much lower stellar masses, that $k=3$ can be
ruled out at $\approx 2\sigma$. For stellar mass black holes, Gao & Li (2023)
argue that cosmological coupling of black holes can explain the distribution
of compact object masses in low mass X-ray binaries, in particular the “mass
gap” between the highest mass neutron stars and lowest mass BHs. However,
Rodriguez (2023a) disfavor coupling based on observations of a BH binary in
the globular cluster NGC 3201 and Andrae & El-Badry (2023) argue that two BH
binaries found by Gaia are likely to have formed with BH masses that were too
small if $k=3$.111It is worth noting that these constraints from BH binaries
assume that the Tolman-Oppenheimer-Volkov limit for the minimum mass of a BH
of $\sim 2.2M_{\odot}$ holds. Also, Ghodla et al. (2023) argue that $k=3$
would result in too many BH mergers compared to observations from
gravitational wave detectors.

###  1.3 Scope of this paper

As the observational picture seems far from settled, we decided to investigate
an independent constraint on BH growth. For cosmologically-coupled mass growth
of BHs to be viable, the mass increase must fit within the constraint provided
by the Soltan argument - that is, the local SMBH mass density must accommodate
both the mass increase due to accretion, and the mass increase due to
cosmological coupling. In this paper, we examine whether cosmologically-
coupled BH growth can be accommodated within the constraints from the Soltan
argument, and, if so, what constraints this can set on the value of $k$. We
adopt the latest measures of both the AGN luminosity function (LF) over cosmic
history, and of the local SMBH mass density, and explore values for the
cosmological coupling strength in the range $0<k<3$. We discuss uncertainties
on our results arising from uncertainties in both $\eta$ and the Eddington
ratio $\lambda$. §2 presents our methods and §3 our sources of data. §4
presents our results and §5 discusses these results in the context of
uncertainties on observed parameters. §6 summarizes our conclusions. We assume
a cosmology with $H_{0}=70\,{\rm kms^{-1}Mpc}$, $\Omega_{\rm M}=0.3$ and
$\Omega_{\Lambda}=0.7$.

##  2 Methods

###  2.1 Generalization of the Soltan Argument for cosmologically-coupled
black holes

The SMBH mass and luminosity are connected for each AGN through

$L=\eta\dot{M}_{\rm acc}c^{2},$ (2)

where $L$ is the bolometric luminosity of the AGN, $\eta$ is the radiative
efficiency, and $\dot{M}_{\rm acc}$ is the mass accretion rate of the BH.
Mass/energy that is not radiated is accreted onto the BH, so the rate of mass
increase of the BH, $\dot{M}_{\bullet}=(1-\eta)L_{\rm{bol}}/(\eta c^{2})$.

In the absence of cosmological coupling (and ignoring the mass density of any
high redshift seed population, which we assume to be negligible throughout
this paper, the mass density in SMBHs today is equal to the integral of the
AGN luminosity density, $U(z)$ over Cosmic Time, modulated by $\eta$:

$\rho_{\bullet}(0)=\int_{0}^{\infty}\left(\frac{1-\eta}{\eta}\right)\frac{U(z)}{c^{2}}\frac{\mathrm{d}t}{\mathrm{d}z}\mathrm{d}z,$
(3)

where $U$:

$\displaystyle U(z)$ $\displaystyle=$
$\displaystyle\int_{\mathrm{L_{min}}}^{\mathrm{L_{max}}}L\frac{\partial\phi(L,z)}{\partial\mathrm{(log_{10}}L)}\,\mathrm{d(log_{10}}L)$
(4) $\displaystyle=$
$\displaystyle\frac{1}{\mathrm{ln}10}\int_{\mathrm{L_{min}}}^{\mathrm{L_{max}}}\frac{\partial\phi(L,z)}{\partial\mathrm{(log_{10}}L)}\,\mathrm{d}L,$
(5)

$\frac{\partial\phi(L,z)}{\partial\mathrm{(log_{10}}L)}$is the bolometric
luminosity function of AGN and

$\frac{\mathrm{d}t}{\mathrm{d}z}=\frac{-1}{(1+z)H(z)},$ (6)

where

$H(z)=H_{0}\sqrt{(1+z)^{3}\Omega_{\mathrm{m}}+\Omega_{\Lambda}}\ ,$ (7)

and $H_{0}$ is the Hubble Constant at the current epoch. If, in addition, BH
mass is increasing with the scale factor of the Universe as per Equation 1,
then the contribution to the BH mass density at redshift $z^{\prime}$,
$\delta\rho_{\bullet}(z^{\prime})$, to that at redshift $z$ is

$\delta\rho_{\bullet}(z)(1+z)^{k}=\delta\rho_{\bullet}(z^{\prime})(1+z^{\prime})^{k}$
(8)

and the version of Equation 3 with cosmologically-coupled BHs is thus:

$\rho_{\bullet}(z)=\frac{1}{(1+z)^{k}}\int_{\infty}^{z}\left(\frac{1-\eta}{\eta}\right)\frac{U(z^{\prime})(1+z^{\prime})^{k-1}}{H(z^{\prime})c^{2}}\,\mathrm{d}z^{\prime}$
(9)

###  2.2 Constraints from the AGN luminosity function at high redshifts

The high redshift AGN luminosity function also provides a constraint on SMBH
growth. If the growth of SMBHs by accretion at high redshifts is too slow,
there is insufficient mass density in SMBHs to power the observed AGN
luminosity density at early times
($z\stackrel{{\scriptstyle>}}{{{}_{\sim}}}2$), when a large fraction of the
SMBH population was actively accreting. We relate the luminosity density in
AGN to the mass density in SMBHs via a mean Eddington ratio, $\lambda=L/L_{\rm
Edd}$, where $L_{\rm Edd}$ is the Eddington luminosity. (We note that both
this analysis and the Soltan analysis above could be made more constraining by
considering the observed distribution of SMBH masses and AGN luminosities and
using a continuity-equation approach (e.g. Raimundo et al., 2012; Tucci &
Volonteri, 2017). However, given the uncertainties in our understanding of the
possible mass/luminosity dependencies of variables such as $\lambda$ and
$\eta$, and uncertainties in the exact forms of the mass and luminosity
functions and SMBH merger histories, we retain an integral approach for this
study.)

##  3 Data

Two measurements are needed to use the formalism in §2. First is a measurement
of the local SMBH mass density. Second is a measure of the accretion history
of SMBHs since the first quasars, via a bolometric AGN luminosity function. We
describe both these sources of data in this section.

###  3.1 The local mass density of SMBHs

Graham & Driver (2007, hereafter GD07) compared several estimates of the local
SMBH mass density available in the literature at the time(see Figure 1, left-
hand panel). These were derived from the SMBH mass – bulge velocity dispersion
($M_{\bullet}-\sigma$) and the SMBH mass – bulge mass ($M_{\bullet}-M_{\rm
bulge}$) relations, and adjusted to a common $H_{0}=70{\rm kms^{-1}Mpc^{-1}}$.
GD07 determined a local SMBH mass density in the range 4.4-5.9 $\times
10^{5}M_{\odot}$Mpc-3. Subsequently, Kormendy & Ho (2013) presented a
recalibration of the $M_{\bullet}-\sigma$ and $M_{\bullet}-M_{\rm bulge}$
relations by excluding pseudobulges and galaxy-galaxy mergers from the fitting
procedure. This resulted in an increase of the $M_{\bullet}/M_{\rm bulge}$
ratio from $\approx 0.2$% to $\approx 0.5$%. Thus, the SMBH mass density
calculated from a stellar mass function and the Kormendy & Ho (2013)
$M_{\bullet}-M_{\rm bulge}$ relation results in a local SMBH mass density
substantially higher than the Graham & Driver (2007) estimate (we note also
that the work of McConnell & Ma (2013) also suggests a relatively high
normalization, $\approx 0.3$%, however, to be conservative, we adopt the
higher Kormendy & Ho (2013) normalization). Furthermore, the difference in the
shapes of the stellar mass functions for disks and bulges means that the
assumption of a single bulge fraction independent of stellar mass is not
strictly correct.

We therefore construct an updated range for the local SMBH mass density as
follows. We use the study of Thanjavur et al. (2016), who split the stellar
mass function into spheroid and disk components. Using their spheroid mass
function and the $M_{\bullet}-M_{\rm bulge}$ relation of Kormendy & Ho (2013)
results in an estimate of $\rho_{\bullet}(0)=1.0\times 10^{6}M_{\odot}$Mpc-3,
higher than that of Graham & Driver (2007), but lower than that obtained by
assuming a fixed bulge fraction and using an all-galaxy stellar mass function.

The uncertainty in this estimate is substantial. It may be too low – there is
evidence that a significant fraction of SMBH mass may be missed by using
spheroid mass of velocity dispersion as a proxy. For example, Voggel et al.
(2019) argue for the existence of a population of ‘stripped’ galaxy nuclei in
the haloes of giant ellipticals. These objects have modest stellar masses, but
can harbor SMBHs of up to $\sim 10^{7}$M⊙. Voggel et al. (2019) further argue
that this population may increase the total SMBH mass density by up to $\sim
30\%$. At high velocity dispersions (or equivalently, high stellar masses) the
$M_{\bullet}-\sigma$ relation may also be biased low, as there are several
examples of ‘overmassive’ SMBHs for their hosts (e.g. Dullo et al., 2021;
Nightingale et al., 2023). Another plausible source of error is ejected SMBHs
(Postman et al., 2012; Chu et al., 2023). Conversely, it is possible our
estimate is too high. Selection biases (Schulze & Wisotzki, 2011; Shankar et
al., 2016) could lead to an overestimate of the SMBH mass density by a factor
$\approx 3$, though the existence of this bias has not been confirmed
observationally. Also, the spheroid mass function still includes pseudobulges,
which, as Kormendy & Ho (2013) emphasize, typically fall below the
$M_{\bullet}-M_{\rm bulge}$ relation (Hopkins et al., 2008, though
pseudobulges are typically less massive than classical bulges, making up
$\stackrel{{\scriptstyle<}}{{{}_{\sim}}}10$% of the stellar mass density in
spheroids). It is beyond the scope of this paper to synthesize these disparate
results into a firm range for $\rho_{\bullet}$. However, a maximum value of
$\approx 1.6\times 10^{6}$M⊙Mpc-1 does not appear to be ruled out (assuming
$\approx 30$% of the mass density is stripped and $\approx 30$% ejected), and
a minimum value of $3\times 10^{5}$M⊙Mpc-1 also seems plausible, based on the
most extreme estimates of bias. We therefore adopt an estimate and range of
${\rm log_{10}}(\rho_{\bullet}(0)/$M⊙Mpc${}^{-1})=6.0^{+0.2}_{-0.5}$ based on
the $M_{\bullet}-M_{\rm bulge}$ relation.

Finally, we note that the recent detection of a gravitational wave background
(GWB) by pulsar timing arrays, if ascribed purely to the mergers of SMBHs,
results in a much higher SMBH mass density estimate than those discussed
above. Using the results of Agazie et al. (2023) we find that the peaks of
their posterior distributions of the most relevant parameters from the fits
based on their Phenom+Astro priors (the galaxy stellar mass function
characteristic density and cutoff mass, the $M_{\bullet}/M_{\rm bulge}$ ratio,
and the logarithmic width of the $M_{\bullet}/M_{\rm bulge}$ distribution)
lead to an estimate of $7.4\times 10^{6}M_{\odot}$Mpc-3. Using the posterior
distribution of $\rho_{\bullet}$ constructed from their Monte-Carlo
simulations results in a 95% confidence range of $2-1400\times
10^{6}M_{\odot}$Mpc-3 and an average of $\approx 1\times 10^{8}M_{\odot}$Mpc-3
(we note also that Casey-Clyde et al., 2022, found a similar result using the
earlier 12.5 yr NanoGrav data release). Given the large uncertainties in the
SMBH merger estimates it is probably too early to be concerned with the
tension between the GWB versus the $M_{\bullet}-M_{\rm bulge}$ relation (see
Afzal et al., 2023, for a discussion of the strength of the tension between
the SMBH binary model and the measured GWB in the context of new physics),
nevertheless we consider the implications of a GWB-based estimate of
$\rho_{\bullet}(0)=7.4\times 10^{6}M_{\odot}$Mpc-3 in our analysis as well.

###  3.2 The bolometric AGN luminosity function

To compute the comoving AGN luminosity density, we need an estimate of the
bolometric AGN luminosity function as a function of redshift. Integrating
under this luminosity function as a function of redshift then gives the
comoving AGN luminosity density. To perform this calculation, we adopt the
Shen et al. (2020) (their “global fit B”) luminosity function (Equation 10),
which agrees well with, but is somewhat more refined than other estimates
(e.g. Hopkins et al., 2007; Lacy et al., 2015; Runburg et al., 2022):

$\frac{\partial\phi(L,z)}{\partial\mathrm{(log_{10}}L)}=\frac{\phi(z)^{*}}{((L/L(z)^{*})^{\gamma_{1}(z)}+(L/L(z)^{*})^{\gamma_{2}(z)})},$
(10)

where the redshift dependencies of $\phi^{*}$, $L^{*}$, $\gamma_{1}$ and
$\gamma_{2}$ are detailed in Shen et al. (2020).

Figure 1:  Middle and right panels: constraints on the value of the
cosmological coupling strength $k$ (Equation 1) based on comparing the growth
of SMBH density via accretion to the local density of SMBHs. Middle: assuming
$\eta=0.1$ (which agrees well with the GD07 estimate of $\rho_{\bullet}(0)$ if
$k=0$) and right: $\eta$=0.3, the maximum theoretical value for a maximally-
spinning BH (Thorne, 1974). In both, the red triangle and error bar represents
our estimate of the local SMBH mass density from the $M_{\bullet}-M_{\rm
bulge}$ relation and the blue square and error bar (which extends well above
the plotted range) the preliminary estimate from the GWB, as discussed in
§3.2. Also in both, the black lines represent the minimum mass density in
SMBHs needed to account for the observed luminosity function for different
mean values of $\lambda$ and the cyan-solid, magenta-dotted and red-dashed
lines the growth in BH mass density for different values of $k$ assuming the
AGN luminosity function of Shen et al. (2020). The shaded region around the
$k=0$ (cyan-solid) line corresponds to the uncertainty in Shen et al. (2020)
(uncertainty ranges for other values of $k$ are similar). Far left a
compilation of estimates of the local mass density in black holes from the
literature: S04: Shankar et al. (2004), M04: Marconi et al. (2004), H07:
Hopkins et al. (2007), G07: Graham et al. (2007), GD07: Graham & Driver
(2007), YL08: Yu & Lu (2008), S09: Shankar et al. (2009), with the cyan shaded
area indicating the range of uncertainty at $z=0$ for the Shen et al. (2020)
luminosity function, assuming $\eta=0.1$.  Figure 2: The allowed values of the
cosmological coupling strength $k$ (Equation 1) and the Salpeter time
(Equation 11) $t_{S}$ for quasars assuming accretion-dominated growth of the
supermassive black hole density. Left: assuming the GD07 estimate of
$\rho_{\bullet}(0)$, middle assuming the best estimate from this paper and
right assuming the value from the peaks of the GWB posterior fits. The hatched
blue regions represent the parameter space where the Salpeter time is too long
to grow enough black holes at $z=4$ to explain the observed luminosity
function and the red regions are those where the Eddington ratio needs to
exceed unity to do this.The yellow hatched regions are excluded based on the
plausible range of $0.05<\eta<0.3$.

##  4 Results

Figure 1 (middle and right panels) plots the SMBH mass densities from §3.1
with the result of integrating Equation 9 from $z=6$ to $z=0$ for fixed $\eta$
values of $0.1$ and $0.3$ (Thorne, 1974). We show curves for $k=0$ (no
cosmological coupling), $k=1$ (e.g. the exact solutions of Faraoni & Jacques
2007) and $k=3$ (Farrah et al., 2023a). This indicates that, if we assume the
GD07 value for the local SMBH mass density and $\langle\eta\rangle\simeq 0.3$,
then an upper limit can be set of $k\stackrel{{\scriptstyle<}}{{{}_{\sim}}}2$.
Alternatively, using the GWB value for the local SMBH mass density, then a
non-zero value of $k$ is required for any reasonable value of $\eta$.

To generalize this argument, we note that, for a given value of $k$, $\eta$
can uniquely determined by matching the observed local black hole mass
density, so if we assume a value of the Eddington ratio, $\lambda$ between the
minimum value needed to explain the AGN luminosity function at high redshift
(assuming all SMBHs are actively-accreting) and $\lambda=1$ then we can define
an allowed region in AGN accretion parameter space as a function of $k$. In
Figure 2, we use the Salpeter time:

$t_{S}=4\times 10^{8}\frac{\eta}{\lambda(1-\eta)}~{}{\rm yr}$ (11)

as a proxy for accretion rate as a function of $\eta$, and plot $k$ vs $t_{S}$
to define an allowed region of parameter space in which quasars can (on
average) lie (assuming no evolution in the mean $\eta$ or $\lambda$ and that
the increase in the SMBH mass density is dominated by accretion). Here, a
similar result is seen. With the GD07 local density we constrain $k$ to $<2$,
with our derived density range then $0<k<3$ can (just) be accommodated (though
only with very high values of $\eta\approx 0.55$ and $\lambda\approx 1$), but
with the GWB local density range then an approximate lower limit can be set of
$k>2$, assuming standard accretion disk dominated growth.

##  5 Discussion

Conventional assumptions about SMBH accretion, combined with estimates of the
local SMBH mass density based on the $M_{\bullet}-M_{\rm bulge}$ relation,
imply an allowed range of $0<k\stackrel{{\scriptstyle<}}{{{}_{\sim}}}2$
(Figure 2). This range is in slight tension ($\sim 90\%$ confidence interval)
with the results in Farrah et al. (2023a), though it is consistent with some
of the recent constraints discussed in §1.2 (e.g. Lei et al., 2023; Rodriguez,
2023b; Andrae & El-Badry, 2023; Ghodla et al., 2023).

Our constraints on $k$ from the Soltan argument are, however, dependent on
four variables: the local SMBH mass density, the integrated luminosity density
from the AGN luminosity function and the mean values of $\lambda$ and $\eta$,
which regulate the total accretion derived from integrating the AGN luminosity
function. We first briefly review the uncertainties on the luminosity density,
$\eta$ and $\lambda$ (the uncertainties on the local SMBH density are
discussed in §3.1), before discussing the effect of these uncertainties on our
limit on $k$.

###  5.1 Uncertainties in the luminosity density

The bolometric luminosity function of AGN is still uncertain, but most of the
uncertainty lies at the faint end at high redshifts, whereas the dominant
contribution to the luminosity density arises around the break in the
luminosity function, where AGN surveys are relatively complete.

As briefly discussed in §3.2, we adopt the luminosity function of Shen et al.
(2020), however, if, for example, we use the mid-infrared luminosity function
of Glikman et al. (2018) and convert to the bolometric luminosity function
using the prescription of Lacy et al. (2015) we recover a black hole mass
density 60% higher than that from the Shen et al. (2020) luminosity function.
(Using this value would be even more constraining on the value of $k$.) Based
on this small difference between two independently-derived luminosity density
estimates we believe that the uncertainties from the AGN luminosity function,
although not negligible, are smaller than those from the local black hole mass
density.

###  5.2 Uncertainties in accretion efficiency

The typical accretion efficiency of an SMBH is poorly constrained (e.g.
Raimundo et al., 2012). The canonical value of $\eta\approx 0.1$ (e.g. Soltan,
1982) is not measured directly, but obtained by matching early estimates of
accretion luminosity to local SMBH mass density. The theoretical limit is
thought to be $\eta\approx 0.3$ (Thorne, 1974) Observational estimates of
$\eta$ vary from $\stackrel{{\scriptstyle<}}{{{}_{\sim}}}0.1$ to potentially
$>0.5$ (e.g. Trakhtenbrot, 2014; Jana et al., 2021; Farrah et al., 2022; Lai
et al., 2023). Several factors are thought to affect $\eta$; for example,
magnetized discs may have $\eta>0.3$ (Avara et al., 2016; Kinch et al., 2021),
but there is no consensus on how these factors interrelate. For Figure 2 we
have assumed the conventional range of $0.05<\eta<0.3$, but a larger range is
certainly not excluded. We note though that high values of $\eta$ reduce the
available time for black hole growth by accretion in the early Universe, which
is already problematic (e.g. Bañados et al., 2018).

###  5.3 Uncertainties in Eddington ratio

Directly observed Eddington ratios in quasars are consistent with $\lambda\sim
0.05$, with a broad tail towards higher Eddington ratios (e.g. Farrah et al.,
2023b). This ratio may however be biased low. Most SMBH mass is likely
accreted in obscured phases (Martínez-Sansigre et al., 2005) during which
Eddington ratios can be higher (Teng et al., 2014; Farrah et al., 2022).
Furthermore, quasar surveys may find AGN whose SMBH masses and Eddington
ratios are biased low and high, respectively (e.g. Shen et al., 2008).

##  6 Conclusions

By analysing the growth of the mass density in SMBHs via accretion we show
that high values of cosmological coupling ($k\approx 3$) needed by Farrah et
al. (2023a) to account for dark energy via BHs is disfavored if conventional
assumptions regarding the local SMBH mass density and accretion mechanisms are
made. However, we cannot rule out scenarios in which $k=3$ is possible, and
indeed, recent estimates of the GWB from pulsar timing arrays favor a much
higher mass density of SMBHs than estimates based on the $M_{\bullet}-M_{\rm
bulge}$ relation and would comfortably allow $k=3$ (Figure 2, right panel).
Better constraints on the local mass function of SMBHs are thus essential both
to our understanding of the GWB and for more reliably constraining any
cosmological coupling of black holes.

We are grateful to L.Z. Kelley and the NANOGrav team for supplying the data
from their simulations of supermassive black hole merger statistics. The
National Radio Astronomy Observatory is a facility of the National Science
Foundation operated under cooperative agreement by Associated Universities,
Inc.

## References

  * Afzal et al. (2023) Afzal, A., Agazie, G., Anumarlapudi, A., et al. 2023, ApJ, 951, L11, doi: 10.3847/2041-8213/acdc91
  * Agazie et al. (2023) Agazie, G., Anumarlapudi, A., Archibald, A. M., et al. 2023, arXiv e-prints, arXiv:2306.16220, doi: 10.48550/arXiv.2306.16220
  * Andrae & El-Badry (2023) Andrae, R., & El-Badry, K. 2023, A&A, 673, L10, doi: 10.1051/0004-6361/202346350
  * Avara et al. (2016) Avara, M. J., McKinney, J. C., & Reynolds, C. S. 2016, MNRAS, 462, 636, doi: 10.1093/mnras/stw1643
  * Bañados et al. (2018) Bañados, E., Venemans, B. P., Mazzucchelli, C., et al. 2018, Nature, 553, 473, doi: 10.1038/nature25180
  * Bardeen (1970) Bardeen, J. M. 1970, Nature, 226, 64, doi: 10.1038/226064a0
  * Beltracchi & Gondolo (2019) Beltracchi, P., & Gondolo, P. 2019, Phys. Rev. D, 99, 044037, doi: 10.1103/PhysRevD.99.044037
  * Cadoni et al. (2023a) Cadoni, M., Murgia, R., Pitzalis, M., & Sanna, A. P. 2023a, arXiv e-prints, arXiv:2309.16444, doi: 10.48550/arXiv.2309.16444
  * Cadoni et al. (2023b) Cadoni, M., Sanna, A. P., Pitzalis, M., et al. 2023b, J. Cosmology Astropart. Phys, 2023, 007, doi: 10.1088/1475-7516/2023/11/007
  * Casey-Clyde et al. (2022) Casey-Clyde, J. A., Mingarelli, C. M. F., Greene, J. E., et al. 2022, ApJ, 924, 93, doi: 10.3847/1538-4357/ac32de
  * Chu et al. (2023) Chu, A., Boldrini, P., & Silk, J. 2023, MNRAS, 522, 948, doi: 10.1093/mnras/stad1033
  * Croker et al. (2020) Croker, K. S., Nishimura, K. A., & Farrah, D. 2020, ApJ, 889, 115, doi: 10.3847/1538-4357/ab5aff
  * Croker & Weiner (2019) Croker, K. S., & Weiner, J. L. 2019, ApJ, 882, 19, doi: 10.3847/1538-4357/ab32da
  * Croker et al. (2021) Croker, K. S., Zevin, M., Farrah, D., Nishimura, K. A., & Tarlé, G. 2021, ApJ, 921, L22, doi: 10.3847/2041-8213/ac2fad
  * Dullo et al. (2021) Dullo, B. T., Gil de Paz, A., & Knapen, J. H. 2021, ApJ, 908, 134, doi: 10.3847/1538-4357/abceae
  * Dymnikova (1992) Dymnikova, I. 1992, General Relativity and Gravitation, 24, 235, doi: 10.1007/BF00760226
  * Elvis et al. (2002) Elvis, M., Risaliti, G., & Zamorani, G. 2002, ApJ, 565, L75, doi: 10.1086/339197
  * Faraoni & Jacques (2007) Faraoni, V., & Jacques, A. 2007, Phys. Rev. D, 76, 063510, doi: 10.1103/PhysRevD.76.063510
  * Farrah et al. (2022) Farrah, D., Efstathiou, A., Afonso, J., et al. 2022, MNRAS, 513, 4770, doi: 10.1093/mnras/stac980
  * Farrah et al. (2023a) Farrah, D., Croker, K. S., Zevin, M., et al. 2023a, ApJ, 944, L31, doi: 10.3847/2041-8213/acb704
  * Farrah et al. (2023b) Farrah, D., Petty, S., Croker, K. S., et al. 2023b, ApJ, 943, 133, doi: 10.3847/1538-4357/acac2e
  * Gao & Li (2023) Gao, S.-J., & Li, X.-D. 2023, ApJ, 956, 128, doi: 10.3847/1538-4357/ace890
  * Ghodla et al. (2023) Ghodla, S., Easther, R., Briel, M. M., & Eldridge, J. J. 2023, The Open Journal of Astrophysics, 6, 25, doi: 10.21105/astro.2306.08199
  * Glikman et al. (2018) Glikman, E., Lacy, M., LaMassa, S., et al. 2018, ApJ, 861, 37, doi: 10.3847/1538-4357/aac5d8
  * Gliner (1966) Gliner, E. B. 1966, Soviet Journal of Experimental and Theoretical Physics, 22, 378
  * Graham & Driver (2007) Graham, A. W., & Driver, S. P. 2007, MNRAS, 380, L15, doi: 10.1111/j.1745-3933.2007.00340.x
  * Graham et al. (2007) Graham, A. W., Driver, S. P., Allen, P. D., & Liske, J. 2007, MNRAS, 378, 198, doi: 10.1111/j.1365-2966.2007.11770.x
  * Hopkins et al. (2008) Hopkins, P. F., Hernquist, L., Cox, T. J., & Kereš, D. 2008, ApJS, 175, 356, doi: 10.1086/524362
  * Hopkins et al. (2007) Hopkins, P. F., Richards, G. T., & Hernquist, L. 2007, ApJ, 654, 731, doi: 10.1086/509629
  * Jana et al. (2021) Jana, A., Naik, S., Chatterjee, D., & Jaisawal, G. K. 2021, MNRAS, 507, 4779, doi: 10.1093/mnras/stab2448
  * Kinch et al. (2021) Kinch, B. E., Schnittman, J. D., Noble, S. C., Kallman, T. R., & Krolik, J. H. 2021, The Astrophysical Journal, 922, 270, doi: 10.3847/1538-4357/ac2b9a
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511, doi: 10.1146/annurev-astro-082708-101811
  * Lacy et al. (2015) Lacy, M., Ridgway, S. E., Sajina, A., et al. 2015, ApJ, 802, 102, doi: 10.1088/0004-637X/802/2/102
  * Lai et al. (2023) Lai, S., Wolf, C., Onken, C. A., & Bian, F. 2023, MNRAS, 521, 3682, doi: 10.1093/mnras/stad651
  * Lei et al. (2023) Lei, L., Zu, L., Yuan, G.-W., et al. 2023, arXiv e-prints, arXiv:2305.03408, doi: 10.48550/arXiv.2305.03408
  * Marconi et al. (2004) Marconi, A., Risaliti, G., Gilli, R., et al. 2004, MNRAS, 351, 169, doi: 10.1111/j.1365-2966.2004.07765.x
  * Martínez-Sansigre et al. (2005) Martínez-Sansigre, A., Rawlings, S., Lacy, M., et al. 2005, Nature, 436, 666, doi: 10.1038/nature03829
  * McConnell & Ma (2013) McConnell, N. J., & Ma, C.-P. 2013, ApJ, 764, 184, doi: 10.1088/0004-637X/764/2/184
  * Nightingale et al. (2023) Nightingale, J. W., Smith, R. J., He, Q., et al. 2023, MNRAS, 521, 3298, doi: 10.1093/mnras/stad587
  * Postman et al. (2012) Postman, M., Lauer, T. R., Donahue, M., et al. 2012, ApJ, 756, 159, doi: 10.1088/0004-637X/756/2/159
  * Prescod-Weinstein et al. (2009) Prescod-Weinstein, C., Afshordi, N., & Balogh, M. L. 2009, Phys. Rev. D, 80, 043513, doi: 10.1103/PhysRevD.80.043513
  * Raimundo et al. (2012) Raimundo, S. I., Fabian, A. C., Vasudevan, R. V., Gandhi, P., & Wu, J. 2012, MNRAS, 419, 2529, doi: 10.1111/j.1365-2966.2011.19904.x
  * Rodriguez (2023a) Rodriguez, C. L. 2023a, ApJ, 947, L12, doi: 10.3847/2041-8213/acc9b6
  * Rodriguez (2023b) —. 2023b, ApJ, 947, L12, doi: 10.3847/2041-8213/acc9b6
  * Runburg et al. (2022) Runburg, J., Farrah, D., Sajina, A., et al. 2022, ApJ, 924, 133, doi: 10.3847/1538-4357/ac37b8
  * Schulze & Wisotzki (2011) Schulze, A., & Wisotzki, L. 2011, A&A, 535, A87, doi: 10.1051/0004-6361/201117564
  * Shankar et al. (2004) Shankar, F., Salucci, P., Granato, G. L., De Zotti, G., & Danese, L. 2004, MNRAS, 354, 1020, doi: 10.1111/j.1365-2966.2004.08261.x
  * Shankar et al. (2009) Shankar, F., Weinberg, D. H., & Miralda-Escudé, J. 2009, ApJ, 690, 20, doi: 10.1088/0004-637X/690/1/20
  * Shankar et al. (2016) Shankar, F., Bernardi, M., Sheth, R. K., et al. 2016, MNRAS, 460, 3119, doi: 10.1093/mnras/stw678
  * Shankar et al. (2020) Shankar, F., Weinberg, D. H., Marsden, C., et al. 2020, MNRAS, 493, 1500, doi: 10.1093/mnras/stz3522
  * Shen et al. (2020) Shen, X., Hopkins, P. F., Faucher-Giguère, C.-A., et al. 2020, MNRAS, 495, 3252, doi: 10.1093/mnras/staa1381
  * Shen et al. (2008) Shen, Y., Greene, J. E., Strauss, M. A., Richards, G. T., & Schneider, D. P. 2008, ApJ, 680, 169, doi: 10.1086/587475
  * Soltan (1982) Soltan, A. 1982, MNRAS, 200, 115, doi: 10.1093/mnras/200.1.115
  * Teng et al. (2014) Teng, S. H., Brandt, W. N., Harrison, F. A., et al. 2014, ApJ, 785, 19, doi: 10.1088/0004-637X/785/1/19
  * Thanjavur et al. (2016) Thanjavur, K., Simard, L., Bluck, A. F. L., & Mendel, T. 2016, MNRAS, 459, 44, doi: 10.1093/mnras/stw495
  * Thorne (1974) Thorne, K. S. 1974, ApJ, 191, 507, doi: 10.1086/152991
  * Trakhtenbrot (2014) Trakhtenbrot, B. 2014, ApJ, 789, L9, doi: 10.1088/2041-8205/789/1/L9
  * Tucci & Volonteri (2017) Tucci, M., & Volonteri, M. 2017, A&A, 600, A64, doi: 10.1051/0004-6361/201628419
  * Voggel et al. (2019) Voggel, K. T., Seth, A. C., Baumgardt, H., et al. 2019, ApJ, 871, 159, doi: 10.3847/1538-4357/aaf735
  * Yu & Lu (2008) Yu, Q., & Lu, Y. 2008, ApJ, 689, 732, doi: 10.1086/592770
  * Yu & Tremaine (2002) Yu, Q., & Tremaine, S. 2002, MNRAS, 335, 965, doi: 10.1046/j.1365-8711.2002.05532.x
  * Zhang & Lu (2017) Zhang, X., & Lu, Y. 2017, Science China Physics, Mechanics, and Astronomy, 60, 109511, doi: 10.1007/s11433-017-9062-1