# Yang-Mills as a Liouville Theory

Stephan Stieberger1, Tomasz R. Taylor2,3, Bin Zhu4
1Max–Planck–Institut für Physik, Werner–Heisenberg–Institut, 80805 München,
Germany
2Department of Physics, Northeastern University, Boston, MA 02115, USA
3Faculty of Physics, University of Warsaw, ul. Pasteura 5, 02-093 Warsaw,
Poland
4School of Mathematics and Maxwell Institute for Mathematical Sciences,
University of Edinburgh, EH9 3FD, UK
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We propose a description of the gluon scattering amplitudes as the inverse
Mellin transforms of the conformal correlators of light operators in two-
dimensional Liouville theory tensored with WZW-like chiral currents on the
celestial sphere. The dimensions of operators are Mellin dual to gluon light
cone energies while their positions are determined by the gluon momentum
directions. Tree-level approximation in Yang-Mills theory corresponds to the
semiclassical limit of Liouville theory. By comparing subleading corrections,
we find $b^{2}=(8\pi^{2})^{-1}\beta_{0}\,g^{2}(M)$, where $b$ is the Liouville
coupling constant, $g(M)$ is the Yang Mills coupling at the renormalization
scale $M$ and $\beta_{0}$ is the one-loop coefficient of the Yang-Mills beta
function.

##

In a recent paper Stieberger:2022zyk , we established an intriguing connection
between the tree-level gluon scattering amplitudes and the correlators of two-
dimensional Liouville theory on the celestial sphere. The gluon amplitudes
were evaluated in the presence of a dilaton source and transformed into
“celestial” amplitudes Pasterski:2017kqt ; Pasterski:2017ylz  by taking Mellin
transforms with respect to the light cone energies of scattered gluons. The
dimensions of Liouville operators were Mellin duals of such energies. Their
positions were determined by the celestial map between the directions of
light-like momenta and points on two-dimensional celestial sphere.111See
reviews of celestial holography in Refs.Strominger:2017zoo ; Raclariu:2021zjz
; Pasterski:2021rjz ; Pasterski:2021raf . Most of the recent work has focused
on extracting CFT data of the putative celestial CFT from scattering
amplitudes in four dimensions, e.g., celestial OPEs Fan:2019emx ; Pate:2019lpp
, infinite-dimensional algebras Fotopoulos:2019vac ; Guevara:2021abz ;
Strominger:2021lvk , differential equations Banerjee:2020zlg ;
Banerjee:2020vnt ; Hu:2021lrx , and connections to twistor theory
Adamo:2019ipt ; Adamo:2021lrv ; Adamo:2021zpw ; Costello:2022wso ;
Bittleston:2023bzp . The celestial amplitudes matched the Liouville
correlators evaluated in the limit of small Liouville coupling, $b\to 0$,
which corresponds to the infinite central charge limit. This construction has
been recently generalized in Ref.Taylor:2023bzj  to celestial amplitudes in
$\mathcal{N}=1$ supersymmetric Yang-Mills theory coupled to dilatons.

In the present work we proceed in the opposite direction. We start from the
operators associated with gluons, constructed as the products of holomorphic
Wess-Zumino-Witten (WZW) curents times the so-called light Liouville
operators. The current part carries the information about gluon gauge charges
and spins. The Liouville part determines their dimensions. The three-point
correlation functions of such operators factorize into a relatively simple,
exactly known WZW correlators times the three-point correlators of light
Liouville operators. The latter ones are known exactly from DOZZ formula
Dorn:1994xn ; Zamolodchikov:1995aa  and can be expressed in terms of
Zamolodchikovs’ $\Upsilon$ function Zamolodchikov:1995aa . We perform inverse
Mellin transformations on the two-dimensional correlators. By using the
celestial map, we construct the corresponding gluon scattering amplitudes. We
can recover the gluon amplitudes, at the tree level and beyond, without the
dilaton background, by taking the limit of inverse Mellin transforms in which
the dilatons decouple. This procedure can be performed exactly at the leading
order in the Liouville coupling ($b\to 0$), corresponding to the tree level
approximation in Yang-Mills theory. We also go beyond the leading order and
identifiy some corrections pointing towards a direct relation between the
Liouville and Yang-Mills couplings.

The Lagrangian density of two-dimensional Liouville theory is given by

${\cal L}=\frac{1}{\pi}{\partial\phi\over\partial
z}{\partial\phi\over\partial\bar{z}}+\mu e^{2b\phi}\ ,$ (1)

where $z$ and $\bar{z}$ are the complexified (Euclidean) spacetime
coordinates, $b$ is the dimensionless Liouville coupling constant and $\mu$ is
the “cosmological constant” scale parameter. The theory has a “background
charge at infinity,”

$Q=b+\frac{1}{b},$ (2)

which is related to the central charge by

$c=1+6Q^{2}.$ (3)

The “light” primary field operators have the form:

$V_{\sigma}(z,\bar{z})=e^{2\sigma b\phi(z,\bar{z})},$ (4)

with the exponents parametrized by $b$-independent parameters $\sigma$. Their
conformal dimensions are given by

$d(\sigma)=2\sigma+2b^{2}\sigma(1-\sigma)\ .$ (5)

We introduce spin and gauge charges into the two-dimensional system by
including a WZW-like holomorphic sector. The WZW current $J^{a}(z)$, with $a$
labeling the adjoint representation of the Lie group, has chiral weights
$(h,\bar{h})=(1,0)$. We also include another operator in the adjoint
representation, $\widehat{J}^{a}(z)$, with $(h,\bar{h})=(-1,0)$. The only
property of this chiral system222For a more detailed discussions of this
chiral system, see Ref.Costello:2022wso . relevant to our discussion is the
form of the three-point correlator

$\big{\langle}\widehat{J}^{a_{1}}(z_{1})\widehat{J}^{a_{2}}(z_{2})J^{a_{3}}(z_{3})\big{\rangle}=f^{a_{1}a_{2}a_{3}}\frac{z_{12}^{3}}{z_{23}z_{31}},$
(6)

where $z_{ij}=z_{i}-z_{j}$ and $f^{a_{1}a_{2}a_{3}}$ are the structure
constants.

We construct the operators associated with the positive helicity gluons in the
following way:

$O^{+a}_{\Delta}(z,\bar{z})=F_{+}(\Delta,\mu,b)\,J^{a}(z)e^{2\sigma(\Delta-1)b\phi(z,\bar{z})}\,,$
(7)

where $F_{+}(\Delta,\mu,b)$ is a normalization factor and $2\sigma(\Delta-1)$
ensures dimension $\Delta-1$ of the Liouville operator. At the leading order
${\cal O}(b^{0})$, $2\sigma(\Delta-1)=\Delta-1$. Similarly, for the negative
helicity gluon,

$O^{-a}_{\Delta}(z,\bar{z})=F_{-}(\Delta,\mu,b)\,J^{a}(z)e^{2\sigma(\Delta+1)b\phi(z,\bar{z})}\,,$
(8)

Note that the normalization factors $F_{\pm}(\Delta,\mu,b)$ depend on the
dimensions $\Delta$, therefore they contribute to inverse Mellin transforms in
a nontrivial way.

We are interested in the “MHV” correlator

$\displaystyle\hskip-38.74988pt\big{\langle}O^{-a_{1}}_{\Delta_{1}}(z_{1},\bar{z}_{1})O^{-a_{2}}_{\Delta_{2}}(z_{2},\bar{z}_{2})O^{+a_{3}}_{\Delta_{3}}(z_{3},\bar{z}_{3})\big{\rangle}=f^{a_{1}a_{2}a_{3}}\displaystyle\frac{z_{12}^{3}}{z_{23}z_{31}}F_{1-}F_{2-}F_{3+}\\!\times~{}~{}~{}~{}~{}~{}~{}~{}$
$\displaystyle\times(z_{12}\bar{z}_{12})^{\frac{\Delta_{3}-\Delta_{1}-\Delta_{2}-3}{2}}(z_{23}\bar{z}_{23})^{\frac{\Delta_{1}-\Delta_{2}-\Delta_{3}+1}{2}}(z_{13}\bar{z}_{13})^{\frac{\Delta_{2}-\Delta_{1}-\Delta_{3}+1}{2}}\times
C(\alpha_{1},\alpha_{2},\alpha_{3})\ ,$

where the three-point Liouville coefficient is given by the famous DOZZ
formula Dorn:1994xn ; Zamolodchikov:1995aa :

$\displaystyle C(\alpha_{1},\alpha_{2},\alpha_{3})=$
$\displaystyle\Big{[}\pi\mu\gamma(b^{2})b^{2-2b^{2}}\Big{]}^{(Q-\sum\alpha_{i})/b}\times$
(10)
$\displaystyle\frac{\Upsilon_{0}\Upsilon(2\alpha_{1})\Upsilon(2\alpha_{2})\Upsilon(2\alpha_{3})}{\Upsilon(\alpha_{1}+\alpha_{2}+\alpha_{3}-Q)\Upsilon(\alpha_{1}+\alpha_{2}-\alpha_{3})\Upsilon(\alpha_{2}+\alpha_{3}-\alpha_{1})\Upsilon(\alpha_{3}+\alpha_{1}-\alpha_{2})}\,,$

in our case specified to the case of light operators with
$\alpha_{i}=\sigma_{i}b$. Here, $\Upsilon$ is the function defined in
Zamolodchikovs’ Ref.Zamolodchikov:1995aa .

The semiclassical ($b\to 0$) limit of the three-point correlator of light
Liouville fields has been studied before by Harlow, Maltz and Witten
Harlow:2011ny . We use the following formulas from Ref.Harlow:2011ny :

$\displaystyle\Upsilon(x-1/b)$
$\displaystyle=\gamma(x/b-1/b^{2})^{-1}b^{1+2/b^{2}-2x/b}\,\Upsilon(x)\,,$
(11) $\displaystyle\Upsilon_{0}$ $\displaystyle=\frac{\cal
C}{b^{1/2}}\exp\left(-\frac{1}{4b^{2}}\log b+\dots\right)\,,$ (12)
$\displaystyle\Upsilon_{b}(\sigma b)$ $\displaystyle=\frac{{\cal
C}b^{1/2-\sigma}}{\Gamma(\sigma)}\exp\left(-\frac{1}{4b^{2}}\log
b+\dots\right)\,,$ (13)

where $\cal C$ is a constant and $\gamma(x)=\Gamma(x)/\Gamma(1-x)$. In this
way, we find

$\displaystyle C(\sigma_{1}b,\sigma_{2}b,\sigma_{3}b)=$
$\displaystyle\frac{\pi\tilde{\mu}\gamma(1/b^{2})\,\gamma(\sum\sigma_{i}-1-1/b^{2})\,[\pi\mu\gamma(b^{2})b^{-2b^{2}}]^{1-\sum\sigma_{i}}}{b^{5}}$
(14)
$\displaystyle\times\frac{\Gamma(\sigma_{1}+\sigma_{2}+\sigma_{3}-1)\Gamma(\sigma_{1}+\sigma_{2}-\sigma_{3})\Gamma(\sigma_{2}+\sigma_{3}-\sigma_{1})\Gamma(\sigma_{3}+\sigma_{1}-\sigma_{2})}{\Gamma(2\sigma_{1})\Gamma(2\sigma_{2})\Gamma(2\sigma_{3})}\,,$

where the “dual” cosmological constant $\tilde{\mu}$ is related to $\mu$ as
follows

$\displaystyle\pi\tilde{\mu}\gamma(1/b^{2})=(\pi\mu\gamma(b^{2}))^{1/b^{2}}\,.$
(15)

Our goal is to apply the celestial map to the inverse Mellin transform,

$\displaystyle{\cal A}_{3G}(\omega_{i},z_{i},\bar{z}_{i})=\displaystyle
M^{\Delta_{1}+\Delta_{2}+\Delta_{3}-3}\left(\frac{1}{2\pi i}\right)^{3}$
$\displaystyle\displaystyle\\!\\!\int_{c-i\infty}^{c+i\infty}d\Delta_{1}d\Delta_{2}d\Delta_{3}\,\omega_{1}^{-\Delta_{1}}\,\omega_{2}^{-\Delta_{2}}\,\omega_{3}^{-\Delta_{3}}~{}~{}~{}~{}~{}$
$\displaystyle\;\;\times\,\big{\langle}O^{-a_{1}}_{\Delta_{1}}(z_{1},\bar{z}_{1})O^{-a_{2}}_{\Delta_{2}}(z_{2},\bar{z}_{2})O^{+a_{3}}_{\Delta_{3}}(z_{3},\bar{z}_{3})\big{\rangle}\
,$

where the integrations are performed on the complex plane along the lines of
real constant $c>0$; at the end, we will take the limit of $c\to 0^{+}$. Note
that connecting two to four dimensions necessitates introducing a
“renormalization” scale $M$ in order to ensure the correct mass dimension $-3$
of the three-gluon amplitude. As mentioned before, the integrands depend on
the normalization constants $F_{\pm}$. We will see below that the following
choice leads to the desired result in the semiclassical limit:

$\displaystyle F_{+}(\Delta,\mu,b)=$
$\displaystyle[\pi\mu\gamma(b^{2})b^{-2b^{2}}]^{\sigma(\Delta-1)}\;\Gamma[2\sigma(\Delta-1)]\
,$ (17) $\displaystyle F_{-}(\Delta,\mu,b)=$
$\displaystyle[\pi\mu\gamma(b^{2})b^{-2b^{2}}]^{\sigma(\Delta+1)-1/2}\;\Gamma[2\sigma(\Delta+1)]\
.$ (18)

Then as $b\to 0$, when $2\sigma=\Delta-1$ for positive helicity gluon and
$2\sigma=\Delta+1$ for negative helicity gluon, the leading term becomes

${\cal
A}_{3G}^{(0)}(\omega_{i},z_{i},\bar{z}_{i})=\frac{\pi\,\tilde{\mu}}{b\,M^{2}}f^{a_{1}a_{2}a_{3}}\displaystyle\frac{z_{12}^{3}}{z_{23}z_{31}}I^{(0)}(\omega_{1},\omega_{2},\omega_{3})\
,$ (19)

where

$\displaystyle I^{(0)}(\omega_{1},\omega_{2},\omega_{3})=$
$\displaystyle\displaystyle\\!\\!\\!\hskip-38.74988pt~{}~{}~{}~{}~{}\left(\frac{1}{2\pi
i}\right)^{3}\,\int_{c-i\infty}^{c+i\infty}d\Delta_{1}d\Delta_{2}d\Delta_{3}\,M^{\Delta_{1}+\Delta_{2}+\Delta_{3}-1}\omega_{1}^{-\Delta_{1}}\,\omega_{2}^{-\Delta_{2}}\,\omega_{3}^{-\Delta_{3}}~{}~{}~{}~{}~{}~{}$
$\displaystyle\times\,\Gamma\left(\frac{\Delta_{1}+\Delta_{2}+\Delta_{3}-1}{2}\right)\Gamma\left(\frac{\Delta_{1}+\Delta_{3}-\Delta_{2}-1}{2}\right)\Gamma\left(\frac{\Delta_{2}+\Delta_{3}-\Delta_{1}-1}{2}\right)\Gamma\left(\frac{\Delta_{1}+\Delta_{2}-\Delta_{3}+3}{2}\right)$
$\displaystyle\times\,(z_{12}\bar{z}_{12})^{\frac{\Delta_{3}-\Delta_{1}-\Delta_{2}-3}{2}}(z_{23}\bar{z}_{23})^{\frac{\Delta_{1}-\Delta_{2}-\Delta_{3}+1}{2}}(z_{13}\bar{z}_{13})^{\frac{\Delta_{2}-\Delta_{1}-\Delta_{3}+1}{2}}$

It is convenient to use the integral representation

$\Gamma(z)=\int_{0}^{+\infty}dt\,e^{-t}\,t^{z-1}\,.$ (21)

to rewrite the inverse Mellin transform as

$\displaystyle I^{(0)}(\omega_{1},\omega_{2},\omega_{3})=$
$\displaystyle\frac{1}{M}\,\left(\frac{1}{2\pi
i}\right)^{3}\,\int_{c-i\infty}^{c+i\infty}d\Delta_{1}d\Delta_{2}d\Delta_{3}\,\int_{0}^{+\infty}dt_{0}dt_{1}dt_{2}dt_{3}\,e^{\Delta_{1}x_{1}}e^{\Delta_{2}x_{2}}e^{\Delta_{3}x_{3}}$
(22)
$\displaystyle\times\,e^{-t_{0}-t_{1}-t_{2}-t_{3}}\,\frac{t_{0}^{-\frac{1}{2}}\,t_{1}^{-\frac{1}{2}}\,t_{2}^{-\frac{1}{2}}\,t_{3}^{\frac{3}{2}}}{t_{0}\,t_{1}\,t_{2}\,t_{3}}(z_{12}\bar{z}_{12})^{-\frac{3}{2}}(z_{23}\bar{z}_{23})^{\frac{1}{2}}(z_{13}\bar{z}_{13})^{\frac{1}{2}}$

where

$\displaystyle x_{1}$
$\displaystyle=\frac{1}{2}\ln\left(\frac{M^{2}t_{0}\,t_{1}\,t_{3}\,z_{23}\bar{z}_{23}}{\omega_{1}^{2}\,t_{2}\,z_{12}\bar{z}_{12}\,z_{13}\bar{z}_{13}}\right)$
(23) $\displaystyle x_{2}$
$\displaystyle=\frac{1}{2}\ln\left(\frac{M^{2}t_{0}\,t_{2}\,t_{3}\,z_{13}\bar{z}_{13}}{\omega_{2}^{2}\,t_{1}\,z_{12}\bar{z}_{12}\,z_{23}\bar{z}_{23}}\right)$
(24) $\displaystyle x_{3}$
$\displaystyle=\frac{1}{2}\ln\left(\frac{M^{2}t_{0}\,t_{1}\,t_{2}\,z_{12}\bar{z}_{12}}{\omega_{3}^{2}\,t_{3}\,z_{23}\bar{z}_{23}\,z_{13}\bar{z}_{13}}\right)\,.$
(25)

In terms of these variables,

$\displaystyle
t_{1}=\frac{\omega_{1}\,\omega_{3}\,e^{x_{1}+x_{3}}\,|z_{13}|^{2}}{M^{2}t_{0}}\
,~{}t_{2}=\frac{\omega_{2}\,\omega_{3}\,e^{x_{2}+x_{3}}\,|z_{23}|^{2}}{M^{2}t_{0}}\
,~{}t_{3}=\frac{\omega_{1}\,\omega_{2}\,e^{x_{1}+x_{2}}\,|z_{12}|^{2}}{M^{2}t_{0}}\
.$ (26)

After changing the integration variables from $t_{1},t_{2},t_{3}$ to
$x_{1},x_{2},x_{3}$ and performing inverse Mellin transforms, we obtain

$I^{(0)}(\omega_{1},\omega_{2},\omega_{3})=\frac{2\,\omega_{1}\omega_{2}}{\omega_{3}M^{2}}\int_{0}^{+\infty}\,dt_{0}\,e^{-t_{0}-\frac{Q^{2}}{M^{2}t_{0}}}\,t_{0}^{-2}\
.$ (27)

where

$Q^{2}=\omega_{1}\omega_{2}|z_{12}|^{2}+\omega_{1}\omega_{3}|z_{13}|^{2}+\omega_{2}\omega_{3}|z_{23}|^{2}\,.$
(28)

According to the celestial map,

$\sqrt{\omega_{i}\omega_{j}}z_{ij}=\langle ij\rangle\
,\qquad\omega_{i}\omega_{j}|z_{ij}|^{2}=2p_{i}\cdot p_{j}\ ,$ (29)

therefore

$Q^{2}=(p_{1}+p_{2}+p_{3})^{2}\ ,$ (30)

and $Q=p_{1}+p_{2}+p_{3}$ can be identified as the total momentum of the gluon
system. After inserting the result (27) into Eq.(19) and using

$\int_{0}^{+\infty}\,dt_{0}\,e^{-t_{0}-\frac{Q^{2}}{M^{2}t_{0}}}\,t_{0}^{-2}=2\sqrt{\frac{M^{2}}{Q^{2}}}\,K_{1}\left(2\sqrt{\frac{Q^{2}}{M^{2}}}\right)\,,$
(31)

where $K_{1}$ is a modified Bessel function, we obtain

${\cal
A}_{3G}^{(0)}(\omega_{i},z_{i},\bar{z}_{i})=\frac{4\pi\,\tilde{\mu}}{b\,M^{4}}f^{a_{1}a_{2}a_{3}}\displaystyle\frac{\langle
12\rangle^{3}}{\langle 23\rangle\langle
31\rangle}\sqrt{\frac{M^{2}}{Q^{2}}}\,K_{1}\left(2\sqrt{\frac{Q^{2}}{M^{2}}}\right).$
(32)

Note that Bessel integrals (31) had already appeared in AdS amplitudes
Penedones:2010ue . Here they appear in the inverse Mellin transform of the
WZW-Liouville correlator (32), which at this point seems to be different from
the three-gluon amplitude of Ref.Stieberger:2022zyk  evaluated in Minkowski
space. In the latter case, the amplitude was evaluated in the presence of a
dilaton background, which was taken into account by one insertion of the
dilaton source. It contained the pole $(Q^{2})^{-1}$ originating from the
massless dilaton propagator connecting the source to the gluon system. The
single source approximations, however, can be justified only in the limit of
small $Q^{2}$. In this limit, the Bessel function can be expanded as

$\displaystyle
2\sqrt{\frac{M^{2}}{Q^{2}}}\,K_{1}\left(2\sqrt{\frac{Q^{2}}{M^{2}}}\right)=\frac{M^{2}}{Q^{2}}+\dots\,,$
(33)

therefore

${\cal
A}_{3G}^{(0)}(\omega_{i},z_{i},\bar{z}_{i})=\frac{2\pi\,\tilde{\mu}}{bM^{2}Q^{2}}f^{a_{1}a_{2}a_{3}}\displaystyle\frac{\langle
12\rangle^{3}}{\langle 23\rangle\langle 31\rangle}+\dots$ (34)

We want to match this correlator with the tree-level amplitude

${\cal
A}_{3G}^{(0^{\prime})}(\omega_{i},z_{i},\bar{z}_{i})=\frac{g}{\Lambda\Lambda^{\prime}}f^{a_{1}a_{2}a_{3}}{1\over
Q^{2}}\displaystyle\frac{\langle 12\rangle^{3}}{\langle 23\rangle\langle
31\rangle}+\dots\ ,$ (35)

where $g$ is the Yang-Mills coupling constant, $\Lambda^{-1}$ is the canonical
coupling of the dilaton to the gauge field strength and $\Lambda^{\prime}$
determines the strength of the point-like dilaton source, ${\cal
J}(x)=\delta^{(4)}(x)/\Lambda^{\prime}$. The semiclassical limit of the
Liouville correlator is equal to the tree-level amplitude provided that the
Yang-Mills and dilaton parameters are related to the Liouville parameters and
the renormalization scale in the following way:

$\frac{gM^{2}}{\Lambda\Lambda^{\prime}}=\frac{2\pi\,\tilde{\mu}}{b}\ .$ (36)

The relation between Liouville correlators and Yang-Mills amplitudes can be
extended beyond the semiclassical limit. The limit of $Q^{2}\to 0$ singles out
gluon amplitudes with one insertion of the dilaton source. These amplitudes
contain the dilaton propagator and the coupling of the off-shell dilaton to
the gluon system. It is well known, however, that the dilaton decouples in the
zero-momentum limit Ademollo:1975pf ; Shapiro:1975cz ; DiVecchia:2015jaq .
Namely, the Feynman matrix element with one zero momentum dilaton is given by
the Feynman matrix element evaluated in the absence of dilatons – in our case
in pure Yang-Mills theory. This observation leads to
Proposition:

$\displaystyle{\cal
M}_{3G}(\omega_{i},z_{i},\bar{z}_{i})=\displaystyle\lim_{Q\to 0}$
$\displaystyle\displaystyle\displaystyle\frac{Q^{2}}{(2\pi
i)^{3}}\displaystyle\\!\\!\int_{c-i\infty}^{c+i\infty}d\Delta_{1}d\Delta_{2}d\Delta_{3}\,M^{\Delta_{1}+\Delta_{2}+\Delta_{3}-1}\omega_{1}^{-\Delta_{1}}\,\omega_{2}^{-\Delta_{2}}\,\omega_{3}^{-\Delta_{3}}~{}~{}~{}~{}~{}$
(37)
$\displaystyle\;\;\times\,\big{\langle}O^{-a_{1}}_{\Delta_{1}}(z_{1},\bar{z}_{1})O^{-a_{2}}_{\Delta_{2}}(z_{2},\bar{z}_{2})O^{+a_{3}}_{\Delta_{3}}(z_{3},\bar{z}_{3})\big{\rangle}\
,$

where ${\cal M}_{3G}$ is the exact three-gluon MHV Feynman matrix element (of
mass dimension 1) in Yang-Mills theory. The equation should be supplemented
with a prescription how to replace two-dimensional Liouville parameters on the
r.h.s. by four-dimensional Yang-Mills parameters on the l.h.s. All what we can
extract at the leading perturbative order is written in Eq.(36). We need an
exact and more direct relation, however, between Liouville and Yang-Mills
couplings. It can be extracted by going beyond the leading order on the Yang-
Mills and Liouville sides of Eq.(37).

In Yang-Mills theory, next-to-leading corrections originate from one-loop
diagrams and are of order ${\cal O}(g^{2})$ as compared to the tree level. In
Liouville theory, they are of order ${\cal O}(b^{2})$ and originate from
various sources. First of all, the $\Upsilon$ function has been expanded in
Ref.Harlow:2011ny  to the order ${\cal O}(b^{2}\ln b^{2})$ only, and more work
is needed to reach higher precision. Furthermore, there is a similar
uncertainity in the normalization factors $F_{\pm}$. In addition, DOZZ formula
is written in terms of the exponents $\sigma_{i}$ while the inverse Mellin
transforms involve integrations over the dimensions $\Delta_{i}$. Eq.(5)
implies that at the subleading order

$\displaystyle\sigma_{1}$
$\displaystyle=\frac{\Delta_{1}+1}{2}+\frac{b^{2}}{4}(\Delta_{1}+1)(\Delta_{1}-1)$
$\displaystyle\sigma_{2}$
$\displaystyle=\frac{\Delta_{2}+1}{2}+\frac{b^{2}}{4}(\Delta_{2}+1)(\Delta_{2}-1)$
$\displaystyle\sigma_{3}$
$\displaystyle=\frac{\Delta_{3}-1}{2}+\frac{b^{2}}{4}(\Delta_{3}-1)(\Delta_{3}-3)$
(38)

We leave full analysis of subleading Liouville corrections to future work,
nevertheless already at this point, we can get a preliminary insight by
discussing some consequences of Eq.(38).

After repeating the steps leading to Eq.(22), but now with the exponents
related to dimensions by Eq.(38), we obtain

$\displaystyle I^{(1)}(\omega_{1},\omega_{2},\omega_{3})=$
$\displaystyle\frac{1}{M}\,\left(\frac{1}{2\pi
i}\right)^{3}\,\int_{c-i\infty}^{c+i\infty}d\Delta_{1}d\Delta_{2}d\Delta_{3}\,\int_{0}^{+\infty}dt_{0}dt_{1}dt_{2}dt_{3}\,e^{\Delta_{1}x_{1}}e^{\Delta_{2}x_{2}}e^{\Delta_{3}x_{3}}$
$\displaystyle\times\,e^{-t_{0}-t_{1}-t_{2}-t_{3}}\,\frac{t_{0}^{-\frac{1}{2}}\,t_{1}^{-\frac{1}{2}}\,t_{2}^{-\frac{1}{2}}\,t_{3}^{\frac{3}{2}}}{t_{0}\,t_{1}\,t_{2}\,t_{3}}(z_{12}\bar{z}_{12})^{-\frac{3}{2}}(z_{23}\bar{z}_{23})^{\frac{1}{2}}(z_{13}\bar{z}_{13})^{\frac{1}{2}}$
(39)
$\displaystyle\hskip-25.83325pt\times\left(\frac{t_{0}\,t_{1}\,t_{3}}{t_{2}}\right)^{\frac{b^{2}}{4}(\Delta_{1}+1)(\Delta_{1}-1)}\left(\frac{t_{0}\,t_{2}\,t_{3}}{t_{1}}\right)^{\frac{b^{2}}{4}(\Delta_{2}+1)(\Delta_{2}-1)}\left(\frac{t_{0}\,t_{1}\,t_{2}}{t_{3}}\right)^{\frac{b^{2}}{4}(\Delta_{3}-1)(\Delta_{3}-3)}$

The difference between the present case and Eq.(22) is that the integrals over
dimensions $\Delta_{i}$ become Gaussian instead of delta functions. After
performing these integrals and changing the variables from $t_{1},t_{2},t_{3}$
to $x_{1},x_{2},x_{3}$, we obtain

$\displaystyle
I^{(1)}(\omega_{1},\omega_{2},\omega_{3})=\frac{2\,\omega_{1}\omega_{2}}{M^{2}\,\omega_{3}}\;e^{\big{[}-\frac{b^{2}}{4}\ln(\frac{2\,p_{1}\cdot
p_{2}}{M^{2}})-\frac{b^{2}}{4}\ln(\frac{2\,p_{2}\cdot
p_{3}}{M^{2}})-\frac{b^{2}}{4}\ln(\frac{2\,p_{1}\cdot p_{3}}{M^{2}})\big{]}}$
$\displaystyle~{}~{}~{}~{}\times\int_{-\infty}^{+\infty}dx_{1}dx_{2}dx_{3}\,\prod_{i=1}^{3}\frac{1}{\sqrt{\epsilon_{i}(x_{i})}}e^{\frac{-x_{i}^{2}}{\epsilon_{i}(x_{i})}+x_{i}(1-\frac{b^{2}}{2})}$
(40)
$\displaystyle~{}~{}~{}~{}~{}\times\int_{0}^{+\infty}dt_{0}\,t_{0}^{-2}\,e^{-t_{0}-\frac{e^{x_{1}+x_{3}}\,\omega_{1}\,\omega_{3}\,|z_{13}|^{2}}{M^{2}\,t_{0}}-\frac{e^{x_{2}+x_{3}}\,\omega_{2}\,\omega_{3}\,|z_{23}|^{2}}{M^{2}\,t_{0}}-\frac{e^{x_{1}+x_{2}}\,\omega_{1}\,\omega_{2}\,|z_{12}|^{2}}{M^{2}t_{0}}}\
,$ (41)

where

$\displaystyle\epsilon_{1}(x_{1})=$ $\displaystyle\pi
b^{2}\Big{[}\ln\Big{(}\frac{2\,p_{1}\cdot p_{3}\,p_{1}\cdot
p_{2}}{M^{2}\,p_{2}\cdot p_{3}}\Big{)}+2x_{1}\Big{]},$
$\displaystyle\epsilon_{2}(x_{2})=$ $\displaystyle\pi
b^{2}\Big{[}\ln\Big{(}\frac{2\,p_{2}\cdot p_{3}\,p_{1}\cdot
p_{2}}{M^{2}\,p_{1}\cdot p_{3}}\Big{)}+2x_{2}\Big{]},$ (42)
$\displaystyle\epsilon_{3}(x_{3})=$ $\displaystyle\pi
b^{2}\Big{[}\ln\Big{(}\frac{2\,p_{1}\cdot p_{3}\,p_{2}\cdot
p_{3}}{M^{2}\,p_{1}\cdot p_{2}}\Big{)}+2x_{3}\Big{]}.$

Since $\epsilon_{i}(x_{i})\sim b^{2}$, we can use the expansion

$\frac{1}{\sqrt{4\pi\epsilon}}e^{\frac{-x^{2}}{4\epsilon}}=e^{\epsilon\partial_{x}^{2}}\,\delta(x)\,.$
(43)

which yields the delta functions fixing $x_{i}=0$ at the leading order. After
expanding the remaining factors, we obtain

$\displaystyle
I^{(1)}(\omega_{1},\omega_{2},\omega_{3})=I^{(0)}(\omega_{1},\omega_{2},\omega_{3})$
$\displaystyle\displaystyle\hskip-91.04872pt\times\Big{[}1-\frac{b^{2}}{4}\ln\Big{(}\frac{2\,p_{1}\cdot
p_{2}}{M^{2}}\Big{)}-\frac{b^{2}}{4}\ln\Big{(}\frac{2\,p_{2}\cdot
p_{3}}{M^{2}}\Big{)}-\frac{b^{2}}{4}\ln\Big{(}\frac{2\,p_{1}\cdot
p_{3}}{M^{2}}\Big{)}\Big{]}+\dots\ .$

The presence of logarithmic corrections in Liouville theory indicates that the
arbitrary mass scale $M$, introduced as a parameter linking Liouville and
Yang-Mills theories, plays the role of renormalization scale in four
dimensions. Assuming that this is indeed the case, we can extract a more
precise relation between Liouville and Yang-Mills couplings by comparing Eq.(
‣ Yang-Mills as a Liouville Theory) with the one-loop correction to the
scattering amplitude of one dilaton with three gluons.

The one-loop corrections to the dilaton-gluon amplitudes have been studied
before in Ref.Inami:1982xt ; Badger:2007si ; Pereira:2022mne . We are
interested in the ultraviolet divergent part only, which after renormalization
leads to the logarithmic running of the gauge coupling $g(Q^{2})$ and of the
dilaton coupling $1/\Lambda$. For three gluons Inami:1982xt ; Badger:2007si ;
Pereira:2022mne :

$\frac{g}{\Lambda}(Q^{2})=\frac{g}{\Lambda}(M^{2})\left[1-\frac{3g^{2}(M^{2})}{2(4\pi)^{2}}\beta_{0}\,\ln\Big{(}\frac{Q^{2}}{M^{2}}\Big{)}\right]\,.$
(45)

where $\beta_{0}=11c_{A}/3$ ($c_{A}$ is the Casimir operator in the adjoint
representation of the gauge group) is the one-loop coefficient of the Yang-
Mills beta function. By comparing the renormalization scale dependence of Eqs.
( ‣ Yang-Mills as a Liouville Theory) and (45), we find

$b^{2}=\frac{\beta_{0}\,g^{2}(M)}{8\pi^{2}}\ .$ (46)

This relation should be taken with a grain of salt though, because it is based
on a partial analysis only of the subleading Liouville corrections.

We admit that the Proposition (37), together with the relation (46) contain
very strong statements. Does it make sense talking about $\it exact$ gluon
scattering amplitudes at all? Evidently, Yang-Mills theory confines gluons and
has a mass gap. Nevertheless, gluon-like states (jets) are physically
observable and according to our proposal, they are described by light
Liouville operators. Massive glueballs are probably described by some other
type of operators and their amplitudes have more string-like character.

## Acknowledgements

TRT is supported by the National Science Foundation under Grants Number
PHY-1913328 and PHY-2209903, by the NAWA Grant “Celestial Holography of
Fundamental Interactions” and by the Simons Foundation Grant MP-
SCMPS-00001550-05. Any opinions, findings, and conclusions or recommendations
expressed in this material are those of the authors and do not necessarily
reflect the views of the National Science Foundation. BZ is supported by the
Royal Society.

## References

  * (1) S. Stieberger, T. R. Taylor and B. Zhu, “Celestial Liouville theory for Yang-Mills amplitudes,” Phys. Lett. B 836, 137588 (2023) doi:10.1016/j.physletb.2022.137588 [arXiv:2209.02724 [hep-th]].
  * (2) S. Pasterski and S. H. Shao, “Conformal basis for flat space amplitudes,” Phys. Rev. D 96, no.6, 065022 (2017) doi:10.1103/PhysRevD.96.065022 [arXiv:1705.01027 [hep-th]].
  * (3) S. Pasterski, S. H. Shao and A. Strominger, “Gluon Amplitudes as 2d Conformal Correlators,” Phys. Rev. D 96, no.8, 085006 (2017) doi:10.1103/PhysRevD.96.085006 [arXiv:1706.03917 [hep-th]].
  * (4) A. Strominger, “Lectures on the Infrared Structure of Gravity and Gauge Theory,” [arXiv:1703.05448 [hep-th]].
  * (5) A. M. Raclariu, “Lectures on Celestial Holography,” [arXiv:2107.02075 [hep-th]].
  * (6) S. Pasterski, “Lectures on celestial amplitudes,” Eur. Phys. J. C 81, no.12, 1062 (2021) doi:10.1140/epjc/s10052-021-09846-7 [arXiv:2108.04801 [hep-th]].
  * (7) S. Pasterski, M. Pate and A. M. Raclariu, “Celestial Holography,” [arXiv:2111.11392 [hep-th]].
  * (8) W. Fan, A. Fotopoulos and T. R. Taylor, “Soft Limits of Yang-Mills Amplitudes and Conformal Correlators,” JHEP 05, 121 (2019) doi:10.1007/JHEP05(2019)121 [arXiv:1903.01676 [hep-th]].
  * (9) M. Pate, A. M. Raclariu, A. Strominger and E. Y. Yuan, “Celestial operator products of gluons and gravitons,” Rev. Math. Phys. 33, no.09, 2140003 (2021) doi:10.1142/S0129055X21400031 [arXiv:1910.07424 [hep-th]].
  * (10) A. Fotopoulos, S. Stieberger, T. R. Taylor and B. Zhu, “Extended BMS Algebra of Celestial CFT,” JHEP 03, 130 (2020) doi:10.1007/JHEP03(2020)130 [arXiv:1912.10973 [hep-th]].
  * (11) A. Guevara, E. Himwich, M. Pate and A. Strominger, “Holographic symmetry algebras for gauge theory and gravity,” JHEP 11, 152 (2021) doi:10.1007/JHEP11(2021)152 [arXiv:2103.03961 [hep-th]].
  * (12) A. Strominger, “w(1+infinity) and the Celestial Sphere,” [arXiv:2105.14346 [hep-th]].
  * (13) S. Banerjee, S. Ghosh and P. Paul, “MHV graviton scattering amplitudes and current algebra on the celestial sphere,” JHEP 02, 176 (2021) doi:10.1007/JHEP02(2021)176 [arXiv:2008.04330 [hep-th]].
  * (14) S. Banerjee and S. Ghosh, “MHV gluon scattering amplitudes from celestial current algebras,” JHEP 10, 111 (2021) doi:10.1007/JHEP10(2021)111 [arXiv:2011.00017 [hep-th]].
  * (15) Y. Hu, L. Ren, A. Y. Srikant and A. Volovich, “Celestial dual superconformal symmetry, MHV amplitudes and differential equations,” JHEP 12, 171 (2021) doi:10.1007/JHEP12(2021)171 [arXiv:2106.16111 [hep-th]].
  * (16) T. Adamo, L. Mason and A. Sharma, “Celestial amplitudes and conformal soft theorems,” Class. Quant. Grav. 36, no.20, 205018 (2019) doi:10.1088/1361-6382/ab42ce [arXiv:1905.09224 [hep-th]].
  * (17) T. Adamo, L. Mason and A. Sharma, “Celestial $w_{1+\infty}$ Symmetries from Twistor Space,” SIGMA 18, 016 (2022) doi:10.3842/SIGMA.2022.016 [arXiv:2110.06066 [hep-th]].
  * (18) T. Adamo, W. Bu, E. Casali and A. Sharma, “Celestial operator products from the worldsheet,” JHEP 06, 052 (2022) doi:10.1007/JHEP06(2022)052 [arXiv:2111.02279 [hep-th]].
  * (19) K. Costello and N. M. Paquette, “Celestial holography meets twisted holography: 4d amplitudes from chiral correlators,” JHEP 10, 193 (2022) doi:10.1007/JHEP10(2022)193 [arXiv:2201.02595 [hep-th]].
  * (20) R. Bittleston, S. Heuveline and D. Skinner, “The Celestial Chiral Algebra of Self-Dual Gravity on Eguchi-Hanson Space,” [arXiv:2305.09451 [hep-th]].
  * (21) T. R. Taylor and B. Zhu, “Celestial Supersymmetry,” JHEP 06, 210 (2023) doi:10.1007/JHEP06(2023)210 [arXiv:2302.12830 [hep-th]].
  * (22) H. Dorn and H. J. Otto, “Two and three point functions in Liouville theory,” Nucl. Phys. B 429, 375-388 (1994) doi:10.1016/0550-3213(94)00352-1 [arXiv:hep-th/9403141 [hep-th]].
  * (23) A. B. Zamolodchikov and A. B. Zamolodchikov, “Structure constants and conformal bootstrap in Liouville field theory,” Nucl. Phys. B 477, 577-605 (1996) doi:10.1016/0550-3213(96)00351-3 [arXiv:hep-th/9506136 [hep-th]].
  * (24) D. Harlow, J. Maltz and E. Witten, “Analytic Continuation of Liouville Theory,” JHEP 12, 071 (2011) doi:10.1007/JHEP12(2011)071 [arXiv:1108.4417 [hep-th]].
  * (25) J. Penedones, “Writing CFT correlation functions as AdS scattering amplitudes,” JHEP 03, 025 (2011) doi:10.1007/JHEP03(2011)025 [arXiv:1011.1485 [hep-th]].
  * (26) M. Ademollo, A. D’Adda, R. D’Auria, F. Gliozzi, E. Napolitano, S. Sciuto and P. Di Vecchia, “Soft Dilations and Scale Renormalization in Dual Theories,” Nucl. Phys. B 94, 221-259 (1975) doi:10.1016/0550-3213(75)90491-5
  * (27) J. A. Shapiro, “On the Renormalization of Dual Models,” Phys. Rev. D 11, 2937 (1975) doi:10.1103/PhysRevD.11.2937
  * (28) P. Di Vecchia, R. Marotta, M. Mojaza and J. Nohle, “New soft theorems for the gravity dilaton and the Nambu-Goldstone dilaton at subsubleading order,” Phys. Rev. D 93, no.8, 085015 (2016) doi:10.1103/PhysRevD.93.085015 [arXiv:1512.03316 [hep-th]].
  * (29) T. Inami, T. Kubota and Y. Okada, “Effective Gauge Theory and the Effect of Heavy Quarks in Higgs Boson Decays,” Z. Phys. C 18, 69-80 (1983) doi:10.1007/BF01571710
  * (30) S. D. Badger, E. W. N. Glover and K. Risager, “One-loop phi-MHV amplitudes using the unitarity bootstrap,” JHEP 07, 066 (2007) doi:10.1088/1126-6708/2007/07/066 [arXiv:0704.3914 [hep-ph]].
  * (31) A. Pereira, A. Cherchiglia, M. Sampaio and B. Hiller, “Higgs boson decay into gluons in a 4D regularization: IR cancellation without evanescent fields to NLO,” Eur. Phys. J. C 83, no.1, 73 (2023) doi:10.1140/epjc/s10052-023-11173-y [arXiv:2207.05187 [hep-ph]].


# Weak Cosmic Censorship Conjecture of Hairy Black Holes in Einstein Gravity

Aofei Sang<EMAIL_ADDRESS>College of Education for the Future,
Beijing Normal University, Zhuhai 519087, China  Department of Physics,
Beijing Normal University, Beijing 100875, China  Jie Jiang  Corresponding
author<EMAIL_ADDRESS>College of Education for the Future, Beijing
Normal University, Zhuhai 519087, China  Department of Physics, Beijing Normal
University, Beijing 100875, China

###### Abstract

Apart from the Kerr-Newman black holes, the hairy black holes in general
relativity have been widely investigated in the gravity and cosmology. In this
paper, we extend the Sorce-Wald method to prove the weak cosmic censorship
conjecture (WCCC) of the static and spherically symmetric hairy black holes in
the Einstein gravity without using the explicit expressions of the metric and
lagrangian for matter fields. We examine the WCCC in the collision process
without requiring the spherically symmetry of the perturbation matter fields.
After assuming the stability condition of spacetime and applying the Gaussian
null coordinates into the variational identities, we derive the first two
order perturbation inequalities which reflects the null energy condition of
the matter fields. As a result, we find that nearly extremal static hairy
black hole cannot be destroyed in the above perturbation process under the
second-order approximation. Our result implies that the valid of the WCCC is
universal for the hairy black hole in the Einstein gravity as long as the
matter fields satisfy the null energy condition under the perturbation level.

##  I Introduction

The weak cosmic censorship conjecture (WCCC) Penrose:1969pc ; Wald:1997wa  is
one of the unproved conjectures in classical general relativity. It is
proposed by Penrose Penrose:1969pc  and states that there is no naked
singularity in our universe. This conjecture ensures that the spacetime
singularity caused by the gravitational collapsing body is always hidden
inside the event horizon of the black hole. It also means that any physical
process cannot destroy the black hole and cause the naked singularity. To
examine the WCCC, Wald proposed a gedanken experiment in 1974 Wald:1974wl , in
which they considered a test particle throwing into an extremal Kerr-Newman
(KN) black hole and showed that any particles which can enter the horizon
cannot destroy the black hole under first-order approximation. After that,
Hubeny Hubeny:1998ga  started with a nearly extremal black hole and found that
there exists some situation under the second-order approximation where the
black hole can be destroyed if the test particle satisfies some certain
conditions. Since then, Hubeny-type violations were found in many black holes
and gravitational theories Jacobson:2010iu ; Chirco:2010rq ; Saa:2011wq ;
Gao:2012ca . However, this method viewed the black hole as a background and
does not consider the second-order effects, such as the self-force effect,
radiative and finite-size effects Jacobson:2010iu . Later, numerical
simulation considering the self-force effect on particles indicates the
Hubeny-type violation will not happen Zimmerman:2012zu ; Colleoni:2015ena ;
Barausse:2011vx ; Barausse:2010ka .

In 2017, Sorce and Wald Sorce:2017dst  put forward a new version of gedanken
experiment. By assuming the stability condition and the null energy condition
of the matter field, they derived the first-order and second-order
perturbation inequalities using the variational identities. Based on these
inequalities, they proved that a nearly extremal KN black hole cannot be over-
charged or over-spun under the second-order approximation of perturbation. In
this method, they consider the full dynamical process of both the spacetime
and the matter fields, and therefore self-force effect, finite-size effect,
and other second-order corrections are automatically taken into account.
Recently, Sorce-Wald gedanken experiments are widely investigated in various
black hole models and gravitational theories.

The no-hair theorem states that a stationary black hole in Einstein-Maxwell
theory is only determined by the mass, electric charge, and angular momentum
Israel:1967wq ; Carter:1971zc ; Ruffini:1971bza , i.e., it is described by the
KN black hole solutions. However, the hairy black holes can be found when the
gravity couples to some extra field degree of freedom, such as the Yang-Mills
fields MSDV ; Bizon:1990sr ; Greene:1992fw , dilaton filed Kanti:1995vq , and
Skyrme hairs Luckock:1986tr ; Droz:1991cx . Most recently, by considering the
complex scalar field without the static symmetry, Herdeiro and Radu Herdeiro1
; Herdeiro2 ; Herdeiro:2015waa  presented a family of spinning scalar-hairy
black hole solutions in Einstein gravity minimally coupled to a massive
complex scalar field. Moreover, Hong, Suzuki, and Yamada Hong:2020miv  also
showed that there exists static and spherical charged scalar hair after the
scalar mass term is considered. Additionally, the Astro observation
Planck:2015fie  implies that there exists dark matter Bertone:2004pz  and dark
energy Peebles:2002gy  in real universe. It is natural to ask whether the WCCC
is valid in these hairy black holes. In this paper, we would like to extend
the Sorce-Wald method to prove the WCCC of static hairy black holes for the
Einstein gravity coupled to some additional matter fields without using the
explicit expressions of the metric and Lagrangian for matter fields.

The outline of this paper is as follows. In Sec. II, we consider the collision
process in the static and spherically symmetric hairy black in Einstein
gravity couples to some additional matter fields. Then, we assume that the
spacetime satisfies the stability condition and introduce two coordinate
systems to describe this assumption. In Sec. III, we review the Noether charge
method of Iyer and Wald and derive the first two order off-shell variational
identities in Einstein gravity. In Secs. IV and V, after applying the Gaussian
null coordinate of the hypersurface into the variational identities, we
derived the first- and second-order perturbation inequalities from the null
energy condition. Then, we prove that the nearly extremal static hairy black
hole cannot be destroyed under the above collision process. Finally, the
conclusion and discussion are presented in Sec. VI.

##  II Perturbed static black hole geometries

In this paper, we consider a static hairy black hole in the Einstein gravity
coupled to some matter fields, such as the electromagnetic field, scalar
field, dilaton field, dark matter, and dark energy. The Lagrangian $n$-form of
this theory is given by

$\displaystyle\begin{aligned}
\bm{L}=\frac{\bm{\epsilon}}{16\pi}R+\bm{L}_{\text{mt}}\,,\end{aligned}$ (1)

in which $R$ is the Ricci scalar, $\bm{\epsilon}$ is the volume element of the
metric $g_{ab}$, and $\bm{L}_{\text{mt}}$ is the Lagrangian of the matter
fields. The equation of motion of this theory is given by

$\displaystyle G_{ab}=8\pi T_{ab}\,,$ (2)

where $G_{ab}=R_{ab}-1/2Rg_{ab}$ is the Einstein tensor, and $T_{ab}$ is the
total energy stress tensor of the matter fields.

Then, we take a general ansatz for the static and spherically symmetric
solution with the line element

$\displaystyle
ds^{2}=-f(r)dv^{2}+2\chi(r)dvdr+r^{2}h_{ij}d\theta^{i}d\theta^{j}\,,$ (3)

with $i,j=1,2\dots,n-2$ in the ingoing Eddington’s coordinate
$O:\\{v,r,\theta\\}$, in which $h_{ij}d\theta^{i}d\theta^{j}$ is the line
element of the unit sphere, and $\chi(r)$ is a positive function of $r$. If
there exists a positive $r_{h}$ such that $f(r_{h})=0$, the above solution
describes a black hole. In this case, the black hole possess an event horizon
at $r=r_{h}$. If there is no positive root of $f(r)=0$, the geometry describes
a naked singularity.

Next, to test whether the black hole can be destroyed, we consider a process
where some matter falls into the black hole from a finite portion of the
horizon and finally, the spacetime settles down to a static state, i.e. we
assume that the black hole satisfies the stability condition Hollands:2012sf ;
Sorce:2017dst . This means, at a very early time, the spacetime geometry is
the same as the background spacetime, i.e. there exists a coordinate system
$O_{1}:\\{v_{1},r_{1},\theta_{1}\\}$ such that the metric has the same form
with Eq. (3), i.e.,

$\displaystyle
ds^{2}=-f(r_{1})dv^{2}_{1}+2\chi(r_{1})dv_{1}dr_{1}+r^{2}_{1}h_{ij}d\theta^{i}_{1}d\theta^{j}_{1}\,.$
(4)

The stability condition implies that there exists another coordinate system
$O_{2}:\\{v_{2},r_{2},\theta_{2}\\}$ such that the late-time geometry can be
described by the line element

$\displaystyle
ds^{2}=-f_{2}(r_{2})dv^{2}_{2}+2\chi_{2}(r_{2})dv_{2}dr_{2}+r^{2}_{2}h_{ij}d\theta^{i}_{2}d\theta^{j}_{2}\,.$
(5)

The form of $f_{2}(r)$ and $\chi_{2}(r)$ are determined by the collision
process. Note that $O_{1}$ and $O_{2}$ are two separated coordinate systems
and both of them cannot cover the whole spacetime. With these setups,
examining the WCCC converts to checking whether $f_{2}(r)$ possesses a
positive root.

Now, we consider a one-parameter family $\phi(\lambda)$ of the field
configurations Sorce:2017dst , in which every element $\phi(\lambda)$
represents a process described above and $\phi(0)$ is the background static
geometry with the line element (3). When $\lambda$ is a small parameter, the
configuration $\phi(\lambda)$ describes a perturbation process on $\phi(0)$.
In the background geometry, both the coordinate systems $O_{1}$ and $O_{2}$
become $O$. The early-time geometry is given by Eq. (4) in the coordinate
$O_{1}$ and the corresponding metric component is independent on $\lambda$.
Noting that The dependence on $\lambda$ of the components show the different
influences of the collision process, the late-time geometric function
$f_{2}(r)$ and $\chi_{2}(r)$ are determined by $\lambda$. Therefore, the line
element of the late-time geometry in $O_{2}$ can be shown by

$\displaystyle
ds^{2}(\lambda)=-f(r_{2},\lambda)dv^{2}_{2}+2\chi(r_{2},\lambda)dv_{2}dr_{2}+r^{2}_{2}h_{ij}d\theta^{i}_{2}d\theta^{j}_{2}\,.\quad$
(6)

Under the background spacetime with $\lambda=0$, we have $f(r,0)=f(r)$ and
$\chi(r,0)=\chi(r)$. The above setups indicate examining the WCCC is
equivalent to finding the root of $f(r,\lambda)$.

##  III Off-shell variational identity in Einstein gravity

In this section, we would like to review the Noether current method proposed
by Iyer and Wald Iyer:1994ys  and derive the off-shell variational identities.
In Einstein gravity, the gravitational part of the Lagrangian $n$-form is

$\displaystyle\begin{aligned}
\bm{L}_{\text{grav}}=\frac{\bm{\epsilon}}{16\pi}R\,.\end{aligned}$ (7)

Variation of the Lagrangian gives

$\delta\bm{L}_{\text{grav}}=\bm{E}_{g}^{ab}\delta g_{ab}+d\bm{\Theta}(g,\delta
g)\,,$ (8)

where

$\displaystyle\begin{aligned}
\bm{E}_{\text{grav}}^{ab}=-\frac{\bm{\epsilon}}{16\pi}G^{ab}=-\frac{\bm{\epsilon}}{2}T^{ab}\,.\end{aligned}$
(9)

describes the gravitational part of the equation of motion, and

$\displaystyle\begin{aligned} \bm{\Theta}_{a_{2}\dots a_{n}}(g,\delta
g)=\frac{1}{16\pi}\bm{\epsilon}_{da_{2}\dots
a_{n}}g^{de}g^{fg}\left(\nabla_{g}\delta g_{ef}-\nabla_{e}\delta
g_{fg}\right)\,.\quad\end{aligned}$ (10)

is the symplectic potential $(n-1)$-form. Here we define

$\displaystyle\begin{aligned}
\delta^{k}\eta(x)=\left.\frac{\partial^{k}\eta(x,\lambda)}{\partial\lambda^{k}}\right|_{\lambda=0}\,.\end{aligned}$
(11)

to denote the $k$th-order variation of the quantity $\eta(x,\lambda)$. Note
that the symbol $\delta$ always denotes the partial derivative respect to
$\lambda$ when we fix the variable in the other slot of $\eta$. We can define
the symplectic current $(n-1)$-form as

$\bm{\omega}(g,\delta_{1}g,\delta_{2}g)=\delta_{1}\bm{\Theta}(g,\delta_{2}g)-\delta_{2}\bm{\Theta}(g,\delta_{1}g)\,,$
(12)

which can be explicitly written as

$\omega_{a_{2}\dots a_{n}}=\frac{1}{16\pi}\bm{\epsilon}_{da_{2}\dots
a_{n}}w^{d},$ (13)

where

$w^{a}=P^{abcdef}\left(\delta_{2}g_{bc}\nabla_{d}\delta_{1}g_{ef}-\delta_{1}g_{bc}\nabla_{d}\delta_{2}g_{ef}\right),$
(14)

with

$\displaystyle\begin{aligned}
P^{abcdef}&=g^{ae}g^{fb}g^{cd}-\frac{1}{2}g^{ad}g^{be}g^{fc}-\frac{1}{2}g^{ab}g^{cd}g^{ef}\\\
&-\frac{1}{2}g^{bc}g^{ae}g^{fd}+\frac{1}{2}g^{bc}g^{ad}g^{ef}.\end{aligned}$
(15)

The Noether current $(n-1)$-form $\bm{J}_{\zeta}$ associated with the vector
field $\zeta^{a}$ is given by

$\bm{J}_{\zeta}=\bm{\Theta}(g,\mathcal{L}_{\zeta}g)-\zeta\cdot\bm{L}.$ (16)

On the other hand, the Noether current can also be written as Iyer:1994ys

$\bm{J}_{\zeta}=\bm{C}_{\zeta}+d\bm{Q}_{\zeta}$ (17)

with the Noether charge

$\displaystyle\begin{aligned} \left(\bm{Q}_{\zeta}\right)_{a_{3}\dots
a_{n}}=-\frac{1}{16\pi}\bm{\epsilon}_{a_{3}\dots
a_{n}cd}\nabla^{c}\zeta^{d}\end{aligned}$ (18)

and the constraints

$\displaystyle\begin{aligned} (\bm{C}_{\zeta})_{a_{2}\dots
a_{n}}=2\bm{\epsilon}_{ea_{2}\dots a_{n}}\zeta^{a}G_{a}^{e}\,.\end{aligned}$
(19)

Then, combining the above expressions and replacing $\zeta^{a}$ by the Killing
vector filed $\xi^{a}$ in the background geometry $\phi(0)$, when $\xi^{a}$ is
fixed under the variation, we can further obtain the first-order variational
identitySorce:2017dst

$\displaystyle\begin{aligned} &\quad
d\left[\delta\bm{Q}_{\xi}-\zeta\cdot\bm{\Theta}\left(g,\delta
g\right)\right]+\xi\cdot\bm{E}_{g}^{ab}\delta
g_{ab}+\delta\bm{C}_{\xi}=0\end{aligned}$ (20)

and second-order variational identity

$\displaystyle\begin{aligned} &\quad
d\left[\delta^{2}\bm{Q}_{\xi}-\zeta\cdot\delta\bm{\Theta}\left(g,\delta
g\right)\right]\\\ &=\bm{\omega}\left(g,\delta g,\mathcal{L}_{\xi}\delta
g\right)-\xi\cdot\bm{E}_{g}^{ab}\delta
g_{ab}-\delta\bm{C}_{\xi}\,.\end{aligned}$ (21)

It is worth noting that the validity of these variational identities is
independent of the gravitational theory imposed in spacetime.

##  IV First-order Sorce-Wald Gedanken experiment

In this section, we consider the first-order Sorce-Wald gedanken experiment in
the above perturbation process. To test the WCCC, we need to know whether
there exists an event horizon after perturbation, which is equivalent to
checking whether the root of $f(r,\lambda)=0$ exists. Then, it is natural to
convert this question to check whether the minimum value of $f(r,\lambda)$ is
negative. Thus, we define a function

$\displaystyle h(\lambda)=f(r_{m}(\lambda),\lambda),$ (22)

where $r_{m}(\lambda)$ denotes the radius where $f(r_{m}(\lambda),\lambda)$ is
the minimum value of $f(r,\lambda)$ and $r_{m}(\lambda)$ can be determined by
$\partial_{r}f(r_{m}(\lambda),\lambda)=0$. For a nearly extremal black hole,
we can define a small parameter $\epsilon$ such that $r_{m}=(1-\epsilon)r_{h}$
with $r_{h}$ satisfied $f(r_{h})=f(r_{h},0)=0$. Then $h(\lambda)$ can be
expanded as

$\displaystyle h(\lambda)$ $\displaystyle=f(r_{m})+\lambda\delta f(r_{m})$
(23) $\displaystyle+\frac{\lambda^{2}}{2}\left[\delta^{2}f(r_{m})+2\delta
r_{m}\delta f^{\prime}(r_{m})+\delta
r_{m}^{2}f^{\prime\prime}(r_{m})\right]+O(\lambda^{3})\,,$

Considering $f^{\prime}(r_{m})=0$, we have

$\displaystyle f^{\prime}(r_{h})=\epsilon
r_{h}f^{\prime\prime}(r_{h})+\mathcal{O}(\epsilon^{2})$ (24)

For the first term of Eq. (23), replacing $r_{m}$ by $(1-\epsilon)r_{h}$,
expanding to the second order of $\epsilon$ and , we can find

$\displaystyle
f(r_{m})=-\frac{1}{2}\epsilon^{2}r_{h}^{2}f^{\prime\prime}(r_{h})+O(\epsilon^{3})$
(25)

Considering these conditions, we can find that the final expression of
$h(\lambda)$ should be

$\displaystyle\begin{aligned} h(\lambda)=&\lambda\delta
f(r_{h})-\frac{1}{2}\epsilon^{2}r_{h}^{2}f^{\prime\prime}(r_{h})-\lambda\epsilon
r_{h}\delta f^{\prime}(r_{h})\\\
&+\frac{\lambda^{2}}{2}\left[\delta^{2}f(r_{h})+2\delta r_{m}\delta
f^{\prime}(r_{h})+\delta r_{m}^{2}f^{\prime\prime}(r_{h})\right]\\\
&+\mathcal{O}(\lambda^{3},\epsilon^{3},\dots)\end{aligned}$ (26)

If the order of $\lambda$ is higher than $\epsilon^{2}$, the leading term of
$h(\lambda)$ is $-1/2\epsilon^{2}r_{h}^{2}f^{\prime\prime}(r_{h})$, which is
negative. Thus the black hole cannot be destroyed in this case. If $\lambda$
has the same order with $\epsilon^{2}$, we can find that the black hole is
much easier to be destroyed when
$-1/2\epsilon^{2}r_{h}^{2}f^{\prime\prime}(r_{h})$ is small. Thus, we consider
the case where $\lambda$ and $\epsilon$ have the same order.

Note that the coordinate $O_{1}$ and $O_{2}$ are two separated coordinate
systems and therefore they cannot be used to describe the collision process.
To obtain a coordinate system to describe the process, we first construct a
Gaussian null coordinate near a hypersurface $\Sigma$ which is the horizon
$r=r_{h}$ in the background geometry. Moreover, it is given by $r_{1}=r_{h}$
in $O_{1}$ at the early time as well as $r_{2}=r_{h}$ in $O_{2}$ at the late
time. For now, the choice of its middle part outside $O_{1}$ and $O_{2}$ is
arbitrary. Here $r_{h}$ is only the horizon radius of the background geometry
and it is independent of the variational parameter $\lambda$, i.e., it is not
the root of the blackening factor $f(r,\lambda)$.

Next, we would like to introduce the Gaussian null (GN) coordinate
$O_{\text{GN}}$, $\\{u,z,\theta^{i}\\}$, near the hypersurface $\Sigma$. In
Gaussian null coordinate system,

$\displaystyle k^{a}=\left(\frac{\partial}{\partial u}\right)^{a}$ (27)

is the tangent vector of $\Sigma$. In the background, the hypersurface is the
horizon, which is a null hypersurface and thus $k^{a}k_{a}=0$ on the
background. In this coordinate system, coordinates $\theta^{i}$ on $\Sigma$
are carried out by $k^{a}$ from one cross-section to another on $\Sigma$. And
the coordinates outside the surface $\Sigma$ is carried out by
$l^{a}=(\partial/\partial z)^{a}$, which is a null vector field satisfying
$l^{a}k_{a}=1$ on $\Sigma$ and $z$ is an affine parameter. i.e., we have

$\displaystyle\begin{aligned} l^{a}\nabla_{a}l^{b}=0\,.\end{aligned}$ (28)

Moreover, we also assume that $l^{a}$ is normal to the cross-section of
$\Sigma$, i.e., we have

$\displaystyle l^{a}\left(\frac{\partial}{\partial\theta^{i}}\right)_{a}=0\,.$
(29)

With the above preparations, in the GN coordinate, the line element in the
neighborhood of $\Sigma$ can always be written as Hollands:2012sf

$\displaystyle\begin{aligned} ds^{2}=2(dz-\alpha
du-\beta_{i}d\theta^{i})du+\gamma_{ij}d\theta^{i}d\theta^{j},\end{aligned}$
(30)

Next, we would like to find the line element in the Gaussian null coordinate
system at the late time. We assume there is a point $p$ on $\Sigma$ with
coordinates $(v_{0},r_{h},\theta^{i}_{0})$ in coordinate system $O_{2}$, in
which $x^{\mu}_{2}$ denotes the coordinates of $O_{2}$. Denote
$x^{\mu}_{2}=\tilde{x}^{\mu}(z)$ is the integral curve of $l^{a}$. This
integral curve connects $p$ with another point $p^{\prime}$ with coordinates
$(v_{0},z,\theta^{i}_{0})$ in $O_{\text{GN}}$ near $\Sigma$. Meanwhile, for a
point near $\Sigma$, the coordinates of $p^{\prime}$ in $O_{2}$ can be
expressed as

$\displaystyle\begin{aligned}
x{}^{\mu}_{2}=x^{\mu}_{0}+z\left.\frac{d\tilde{x}^{\mu}}{dz}\right|_{p}+\frac{z^{2}}{2}\left.\frac{d^{2}\tilde{x}^{\mu}}{dz^{2}}\right|_{p}+\mathcal{O}(z^{3})\,.\end{aligned}$

Since $l^{a}$ is a null vector and satisfies $l^{a}k_{a}|_{\Sigma}=0$, we can
find

$\displaystyle\begin{aligned}
l^{a}|_{\Sigma}=\frac{1}{\chi(r_{h},\lambda)}\left(\frac{\partial}{\partial
r}\right)^{a}\end{aligned}$ (32)

Using $l^{a}\nabla_{a}l^{b}=0$, we can further obtain

$\displaystyle\left.\frac{d^{2}\tilde{x}^{\mu}}{dz^{2}}\right|_{p}=-\left.\Gamma^{\mu}{}_{\nu\sigma}l^{\nu}l^{\sigma}\right|_{p}\,.$
(33)

Then, the coordinate transformation between $O_{2}$ and $O_{\text{GN}}$ can be
written as

$\displaystyle\begin{aligned}
r_{2}&=r_{h}+\frac{z}{\chi(r_{h},\lambda)}-\frac{z^{2}\chi^{\prime}(r_{h},\lambda)}{2\chi(r_{h},\lambda)^{3}}\\\
v_{2}&=v_{0}\,,\quad\quad\theta_{2}^{i}=\theta_{0}^{i}\,,\end{aligned}$ (34)

where we denote

$\displaystyle\begin{aligned}
\chi^{\prime}(r_{h},\lambda)=\left.\frac{d\chi(r,\lambda)}{dr}\right|_{r_{h}}\end{aligned}$
(35)

We can further obtain the line element in $O_{\text{GN}}$

$\displaystyle\begin{aligned} ds^{2}=2(dz-\alpha
du)du+\gamma_{ij}d\theta^{i}d\theta^{j}\,,\end{aligned}$ (36)

where

$\displaystyle\begin{aligned}
\alpha&=\frac{1}{2}f(r_{h},\lambda)+\frac{zf^{\prime}(r_{h},\lambda)}{2\chi(r_{h},\lambda)}\\\
\gamma_{ij}&=r_{h}\left(r_{h}+\frac{2z}{\chi(r_{h},\lambda)}\right)h_{ij}\,,\end{aligned}$
(37)

Next, we choose a gauge such that the GN coordinate system is fixed under the
variation. Taking $\zeta^{a}=\left({\partial}/{\partial u}\right)^{a}$ in the
first-order variational identity (LABEL:variation1) and integrating it over
$\Sigma$, we have

$\displaystyle\begin{aligned}
\int_{B_{1}}\left[\delta\bm{Q}_{\zeta}-\zeta\cdot\bm{\Theta}(g,\delta
g)\right]=-\int_{\Sigma}\delta\bm{C}_{\zeta},\end{aligned}$ (38)

After considering $f(r_{h})=0$ on background and using the line element (36)
of late-time geometry in GN coordinate, we can find the explicit expression of
the left side of Eq. (38),

$\displaystyle\begin{aligned}
\int_{B_{1}}\left[\delta\bm{Q}_{\zeta}-\zeta\cdot\bm{\Theta}(g,\delta
g)\right]=\frac{-\delta f(r_{h})}{8\pi r_{h}\chi(r_{h})}A_{H}\,,\end{aligned}$
(39)

where $A_{H}$ is the area of the horizon in the background geometry. The right
side of Eq. (38) is

$\displaystyle\begin{aligned}
\int_{\Sigma}-\delta\bm{C}_{\zeta}=2\int_{\Sigma}\hat{\bm{\epsilon}}\left(dz\right)_{e}\left(\frac{\partial}{\partial
u}\right)^{a}\delta T_{a}^{e}\,,\end{aligned}$ (40)

where $\hat{\bm{\epsilon}}$ is the induced volume element on $\Sigma$. We
assume the perturbation matter field satisfies the null energy condition ,i.e.
$T_{ab}(\lambda)k^{a}(\lambda)k^{b}(\lambda)\geq 0$ for any null vector
$k^{a}$. We choose the null vector

$\displaystyle\begin{aligned} k^{a}(\lambda)=\left(\frac{\partial}{\partial
u}\right)^{a}+\alpha(\lambda)\left(\frac{\partial}{\partial
z}\right)^{a}\end{aligned}$ (41)

Using the line element (30), we also have the matrix of the inverse metric is

$\displaystyle\begin{aligned} &g^{zu}=1\,,\qquad
g^{zi}=-\gamma^{ij}\beta_{j}\,,\\\
&g^{zz}=2\alpha+\gamma^{ij}\beta_{i}\beta_{j}\,,\qquad
g^{ij}=\gamma^{ij}\,,\end{aligned}$ (42)

in which $\gamma^{ij}$ is the inverse of $\gamma_{ij}$. Then, it is not hard
to verify

$\displaystyle\begin{aligned}
T_{ab}k^{a}k^{b}=T_{u}^{z}-\gamma^{ij}\beta_{i}\beta_{j}T_{uz}+\alpha^{2}T_{zz}+\gamma^{ij}\beta_{j}T_{ui}\,.\end{aligned}$
(43)

For simplification, we have neglected $\lambda$ in the above expression. From
the above relation, we can get

$\displaystyle\delta(T_{ab}k^{a}k^{b})=\delta T_{u}^{z}\,,$ (44)

Here we have used the fact that $\alpha=\beta_{i}=T_{ui}=0$ under the
background geometry. And therefore the null energy condition under the first-
order approximation implies

$\displaystyle\delta T_{u}^{z}\geq 0\,.$ (45)

Then, the variational identity reduces to

$\displaystyle\delta f(r_{h})\leq 0\,.$ (46)

Then, Eq. (26) reduces to

$\displaystyle\begin{aligned}
h(\lambda)\leq\mathcal{O}(\lambda^{2},\varepsilon^{2},\lambda\varepsilon)\,,\end{aligned}$
(47)

which means $f(r,\lambda)$ exists a positive root under the first-order
approximation and therefore the black hole cannot be destroyed under the
first-order approximation. However, there exists an optimal condition $\delta
f(r_{h})=0$ such that the leading term of $h(\lambda)$ is second order, which
means that we cannot determine the sign of $h(\lambda)$ under the first-order
approximation in this case and we need to consider the second-order
approximation. From the above derivation, it is not hard to see that the
optimal condition $\delta f(r_{h})=0$ also implies

$\displaystyle\begin{aligned} \delta T_{kk}=\delta
T_{u}^{z}=0\,.\end{aligned}$ (48)

##  V Second-order Sorce-Wald gedanken experiment

In this section, we consider the second-order Sorce-Wald gedanken experiment
when the optimal condition is satisfied. Note that the choice of the middle of
the hypersurface $\Sigma$ is arbitrary in the last section. In the appendixes,
we have shown that we can choose the hypersurface such that it is a null
hypersurface when $\varepsilon=0$ under the first-order approximation of
$\lambda$, i.e., we have

$\displaystyle\alpha=z\tilde{\alpha}+\mathcal{O}(\varepsilon\lambda,\lambda^{2},\varepsilon^{2})\,,\quad\beta_{i}=z\tilde{\beta}_{i}+\mathcal{O}(\varepsilon\lambda,\lambda^{2},\varepsilon^{2})\,,$
(49)

and the optimal condition can also give

$\displaystyle\delta T_{u}^{z}=\delta T_{uu}=0$ (50)

on this new hypersurface, in which the variation is evaluated in its GN
coordinate system.

Next, we consider the second-order variation identity

$\displaystyle
d\delta\left[\delta\bm{Q}_{\zeta}-\zeta\cdot\bm{\Theta}(g,\delta g)\right]$
(51) $\displaystyle=\bm{\omega}(g;\delta g,\mathcal{L}_{\zeta}\delta
g)-\zeta\cdot\delta(\bm{E}_{g}\delta g)-\delta^{2}\bm{C}_{\zeta}\,,$

and integrate it over the null hypersurface described above. Noting
$\zeta^{a}=\left({\partial}/{\partial u}\right)^{a}$ is tangent to $\Sigma$,
we can get

$\displaystyle\int_{B_{1}}\delta\left[\delta\bm{Q}_{\zeta}-\zeta\cdot\bm{\Theta}(g,\delta
g)\right]$ (52) $\displaystyle=\int_{\Sigma}\bm{\omega}(g;\delta
g,\mathcal{L}_{\zeta}\delta g)-\int_{\Sigma}\delta^{2}\bm{C}_{\zeta}\,,$

We can use the metric (30) to calculate the left side directly. After
considering the optimal condition $\delta f(r_{h})=0$, we can find

$\displaystyle\int_{B_{1}}\delta\left[\delta\bm{Q}_{\zeta}-\zeta\cdot\bm{\Theta}(g,\delta
g)\right]=-\frac{\delta^{2}f(r_{h})}{8\pi r_{h}h(r_{h})}A_{H}\,.$ (53)

Expanding (43) to the second order and using the first order optimal condition
$\delta(T_{ab}k^{a}k^{b})=\delta T_{u}^{z}=0$, (43) gives

$\displaystyle\delta^{2}(T_{ab}k^{a}k^{b})=\delta^{2}T_{u}^{z}+\mathcal{O}(\varepsilon)\geq
0$ (54)

Thus,

$\displaystyle-\int_{\Sigma}\delta^{2}\bm{C}_{\zeta}=\int_{\Sigma}\hat{\epsilon}\delta^{2}T_{u}^{z}\geq\mathcal{O}(\lambda)$
(55)

Next, we define the canonical energy of the perturbation $\delta g$ on
$\Sigma$ by

$\displaystyle\mathcal{E}_{\Sigma}(g,\delta
g)\equiv\int_{\Sigma}\bm{\omega}(g;\delta g,\mathcal{L}_{\xi}\delta g)\,.$
(56)

The canonical energy can be obtained as

$\displaystyle\begin{aligned} &\mathcal{E}_{\Sigma}(g,\delta g)\\\
&=\int_{\Sigma}\bm{\epsilon}_{aa_{1}\dots
a_{n-1}}P^{abcdef}\left(\mathcal{L}_{\xi}\delta g_{bc}\nabla_{d}\delta
g_{ef}-\delta g_{bc}\nabla_{d}\mathcal{L}_{\xi}\delta g_{ef}\right)\\\
&=-2\int_{\Sigma}\tilde{\bm{\epsilon}}P^{zbcdef}\mathcal{L}_{\xi}\delta
g_{bc}\nabla_{d}\delta
g_{ef}+\int_{\Sigma}\mathcal{L}_{\xi}(\tilde{\bm{\epsilon}}P^{zbcdef}\delta
g_{bc}\nabla_{d}\delta g_{ef})\end{aligned}$

Using the Stokes’ theorem, considering the assumption that the perturbation
vanishes on $B_{0}$ and $\delta
g_{ab}=\delta\gamma_{ab}+\mathcal{O}(\epsilon)$ on $\Sigma$ according to (49),
we can get

$\displaystyle\begin{aligned} \mathcal{E}_{\Sigma}&(g,\delta g)\\\
=&-2\int_{\Sigma}\tilde{\bm{\epsilon}}P^{zbcdef}\mathcal{L}_{\xi}\delta\gamma_{bc}\nabla_{d}\delta
g_{ef}+\int_{B_{1}}\hat{\bm{\epsilon}}P^{zbcdef}\delta\gamma_{bc}\nabla_{d}\delta
g_{ef}\\\ &+\mathcal{O}(\varepsilon)\end{aligned}$

Using the late-time metric (6), it is not hard to verify the boundary term
vanishes. Thus, we have

$\displaystyle\begin{aligned} \mathcal{E}_{\Sigma}(g,\delta
g)=-2\int_{\Sigma}\tilde{\bm{\epsilon}}P^{zbcdef}\mathcal{L}_{\xi}\delta\gamma_{bc}\nabla_{d}\delta
g_{ef}+\mathcal{O}(\varepsilon)\end{aligned}$ (59)

Using the optimal condition $\delta\theta=0$, we have
$2\delta\sigma_{ab}=\mathcal{L}_{\xi}\gamma_{ab}$. Then, we can get

$\displaystyle\begin{aligned} \mathcal{E}_{\Sigma}(g,\delta
g)=-4\int_{\Sigma}\tilde{\bm{\epsilon}}\delta\sigma_{bc}A^{bc}+\mathcal{O}(\varepsilon)\,.\end{aligned}$
(60)

with

$\displaystyle\begin{aligned} &A_{ij}=P^{z}{}_{ij}{}^{def}\nabla_{d}\delta
g_{ef}\\\ &=\delta g_{ui;j}-\frac{1}{2}\delta g_{ij;u}=\delta
g_{ui;j}-\delta\sigma_{ij}\end{aligned}$ (61)

Using the Gaussian null coordinate, it is not hard to verify $\delta
g_{ui;j}=0$. Therefore, we have $A_{ij}=-\delta\sigma_{ij}$, i.e.

$\displaystyle\begin{aligned} \mathcal{E}_{\Sigma}(g,\delta
g)=4\int_{\Sigma}\tilde{\bm{\epsilon}}\delta\sigma_{bc}\delta\sigma^{bc}+\mathcal{O}(\varepsilon)\geq\mathcal{O}(\varepsilon)\,.\end{aligned}$
(62)

Therefore, the second-order perturbation inequality is

$\displaystyle\begin{aligned}
\delta^{2}f(r_{h})\leq\mathcal{O}(\epsilon,\lambda)\,.\end{aligned}$ (63)

After considering the first-order optimal condition, (23) becomes

$\displaystyle\begin{aligned} f(\lambda)=&-\epsilon\lambda r_{h}\delta
f^{\prime}(r_{h})-\frac{1}{2}\epsilon^{2}r_{h}^{2}f^{\prime\prime}(r_{h})\\\
&+\frac{1}{2}\lambda^{2}\left[\delta^{2}f(r_{h})+2\delta r_{m}\delta
f^{\prime}(r_{h})+\delta
r_{m}^{2}f^{\prime\prime}(r_{h})\right]\,.\end{aligned}$ (64)

Considering $r_{m}(\lambda)$ is the minimum value of
$f(r_{m}(\lambda),\lambda)$, we have

$\displaystyle\begin{aligned} \delta r_{m}=\frac{\delta
f^{\prime}(r_{m})}{f^{\prime\prime}(r_{m})}=\frac{\delta
f^{\prime}(r_{h})}{f^{\prime\prime}(r_{h})}+\mathcal{O}(\varepsilon)\,,\end{aligned}$
(65)

and therefore, with the identity
$\delta^{2}f(r_{h})\leq\mathcal{O}(\varepsilon)$, we can find (64) can be
further expressed as

$\displaystyle\begin{aligned} f(\lambda)=-\frac{\left[\lambda\delta
f^{\prime}(r_{h})+\epsilon
r_{h}f^{\prime\prime}(r_{h})\right]^{2}}{2f^{\prime\prime}(r_{h})}\,,\end{aligned}$
(66)

under the second-order approximation of $\epsilon$ and $\lambda$. Therefore,
we have $f(\lambda)\leq 0$ under the second-order approximation, which means
the black hole cannot be destroyed.

##  VI Conclusion and Discussion

Apart from the KN black holes, the hairy black holes in general relativity
have been widely investigated in gravity and cosmology. In the current paper,
we extended the Sorce-Wald gedanken experiment to study the WCCC of the static
and spherically symmetric hairy black holes in Einstein gravity coupled to
some matter fields. We considered a collision process without requiring the
spherical symmetry of the perturbation matter fields. After assuming the
stability condition of the spacetime and based on the Noether charge method by
Iyer and Wald Iyer:1994ys , we derived the first- and second-order
perturbation inequalities which reflect the null energy condition of the
matter fields. Using these inequalities, we found that the nearly extremal
hairy black holes cannot be destroyed under the second-order approximation of
perturbation and therefore the WCCC is valid in these hairy black holes. Our
result is universal and independent of the explicit expressions of spacetime
metric and the Lagrangian for matter fields. Therefore, it implies that the
WCCC is generally valid for any black holes in Einstein gravity as long as the
matter fields satisfy the null energy condition.

In this paper, we only considered the cases where the background spacetime is
static and spherically symmetric. To give general proof of the WCCC in general
relativity under the perturbation level, it is necessary to consider the
perturbed spinning black holes in which the perturbation matter fields carry
the total angular momentum. Moreover, it is also interesting to extend the
discussion into the modified gravitational theories, such as the Einstein-
scalar-Gauss-Bonnet gravity and Lovelock gravities.

## Acknowledgement

This article is supported by the GuangDong Basic and Applied Basic Research
Foundation with Grant no. 217200003 and the Talents Introduction Foundation of
Beijing Normal University with Grant no. 310432102.

##  Appendix A Proof 1

In this appendix, we would like to proof that the first-order optimal
condition at an arbitrary hypersurface $H_{1}$, which is the horizon on the
background, is equivalent to $\delta T_{kk}=0$ at a null hypersurface $H_{2}$.

We assume the first order optimal condition is satisfied,i.e.

$\displaystyle\begin{aligned} \delta T_{u}^{z}=0\end{aligned}$ (67)

where $T_{u}^{z}=T_{u}^{z}(\tilde{x},\lambda)$ and we use
$O_{1}:\\{\tilde{u},\tilde{z},\tilde{\theta}\\}$ to represent the coordinate
system that we used to calculate the first-order inequality and $H_{1}$ is the
hypersurface given by $\tilde{z}=0$. We will use “$\sim$” to label the
quantities valued in $O_{1}$ later. Then, we assume $O_{2}:\\{u,z,\theta\\}$
is the Gaussian null coordinate around the null hypersurface $H_{2}$ with
$k^{a}=\left({\partial}/{\partial u}\right)^{a}$ the tangent null vector of
$H_{2}$. Taking a point $p_{1}$ on $H_{1}$ with coordinate
$\\{\tilde{x}_{p1}\\}$ and a point $p_{2}$ on $H_{2}$ with coordinate
$\\{\tilde{x}_{p2}\\}$ in $O_{1}$, here we use “$\sim$” to represent the
coordinate value in $O_{1}$. We also assume that the coordinate of $p_{2}$ in
$O_{2}$ is $\\{\hat{x}_{p2}\\}$. We associate $p_{1}$ and $p_{2}$ by demanding
$\hat{x}_{p2}=\tilde{x}_{p1}$. In $O_{1}$, the stress energy tensor can be
expressed as

$T_{ab}(\lambda)=\tilde{T}_{\mu\nu}(\tilde{x},\lambda)(d\tilde{x}^{\mu})_{a}(d\tilde{x}^{\nu})_{b}\,.$
(68)

The optimal condition valued at $p_{1}$ gives

$\displaystyle\delta\tilde{T}_{0}^{1}(\tilde{x}_{p1})=\delta\tilde{T}_{0}^{1}(\hat{x}_{p2})=0\,.$
(69)

Then, in the null hypersurface,

$\displaystyle T_{ab}k^{a}k^{b}|_{p_{2}}$
$\displaystyle=\left.T_{ab}\left(\frac{\partial}{\partial
u}\right)^{a}\left(\frac{\partial}{\partial u}\right)^{b}\right|_{p2}$ (70)
$\displaystyle=\tilde{T}_{\mu\nu}(\tilde{x}_{p2},\lambda)\left.\frac{\partial\tilde{x}^{\mu}}{\partial
u}\right|_{p2}\left.\frac{\partial\tilde{x}^{\nu}}{\partial u}\right|_{p2}\,.$

Then, we would like to calculate the variation of Eq. (70) in $O_{2}$. We
assume the coordinates transformation between $O_{1}$ and $O_{2}$ is
$\tilde{x}_{q}=\tilde{x}(\hat{x}_{q},\lambda)$ for any point $q$ and the
coordinates transformation satisfies $\tilde{x}(\hat{x},0)=\hat{x}=x$. Since
the equation (70) holds for any $p_{2}$ at $H_{2}$, therefore we will get rid
of $p_{2}$ for simplicity in the following calculation. Therefore (70) can be
written as

$\displaystyle\begin{aligned}
T_{kk}(\hat{x},\lambda)=\tilde{T}_{\mu\nu}(\tilde{x}(\hat{x},\lambda),\lambda)\frac{\partial\tilde{x}^{\mu}(\hat{x},\lambda)}{\partial
u}\frac{\partial\tilde{x}^{\nu}(\hat{x},\lambda)}{\partial u}\,.\end{aligned}$
(71)

When we calculate the variation, we fix $\hat{x}$ both side. Then, the
variation of $T_{kk}(\hat{x},\lambda)$ will be

$\displaystyle\begin{aligned} \delta
T_{kk}=&[\delta\tilde{T}_{\mu\nu}(x)+\partial_{\rho}\tilde{T}_{\mu\nu}(x)\delta\tilde{x}^{\rho}(x)]\frac{\partial\tilde{x}^{\mu}(\hat{x},0)}{\partial
u}\frac{\partial\tilde{x}^{\nu}(\hat{x},0)}{\partial u}\\\
&+2\tilde{T}_{\mu\nu}(x)\frac{\partial\tilde{x}^{\mu}(\hat{x},0)}{\partial
u}\frac{\partial\delta\tilde{x}^{\nu}(x)}{\partial u}\,.\end{aligned}$ (72)

Considering the coordinates of a point is the same in $O_{1}$ and $O_{2}$ when
$\lambda=0$, we can find

$\displaystyle\frac{\partial\tilde{x}^{\mu}(\hat{x},\lambda)}{\partial
u}|_{\lambda=0}=\delta_{0}^{\mu}\,.$ (73)

Then, considering $T_{0i}$ components vanish on background and $T_{00}=0$ on
the horizon $z=0$, we can obtain

$\displaystyle\begin{aligned} \delta
T_{kk}=\delta\tilde{T}_{00}(x)+\partial_{1}\tilde{T}_{00}(x)\delta\tilde{x}^{1}(x)+2\tilde{T}_{01}({x})\frac{\partial\delta\tilde{x}^{1}({x})}{\partial
u}\,.\quad\end{aligned}$ (74)

Using the fact that $g_{ab}(\frac{\partial}{\partial
u})^{a}(\frac{\partial}{\partial u})^{b}|_{H2}=0$, with a similar calculation,
we can find

$\displaystyle\begin{aligned} &\delta
g_{uu}=\delta\tilde{g}_{00}(x)+\partial_{1}\tilde{g}_{00}(x)\delta\tilde{x}^{1}(x)+2\tilde{g}_{01}({x})\frac{\partial\delta\tilde{x}^{1}({x})}{\partial
u}=0\,,\\\ &2\frac{\partial\delta\tilde{x}^{1}}{\partial
u}=-\delta\tilde{g}_{00}(x)-\partial_{1}\tilde{g}_{00}\delta\tilde{x}^{1}\,.\end{aligned}$
(75)

The zero component of the conservation law $\nabla_{a}T^{ab}=0$ gives

$\displaystyle\begin{aligned}
g^{ac}\partial_{c}T_{a0}-\Gamma^{da}{}_{a}T_{d0}-\Gamma^{da}{}_{0}T_{ad}=0\,.\end{aligned}$
(76)

On the background, the line element can be expressed as

$\displaystyle\begin{aligned} ds^{2}=2(dz-z\alpha
du-z\beta_{i}d\theta^{i})du+\gamma_{ij}d\theta^{i}d\theta^{j}\,.\end{aligned}$
(77)

It is early to find

$\displaystyle\begin{aligned}
g^{ac}\Gamma^{1}{}_{ac}=-2\alpha\,,\quad\Gamma^{\mu\nu}{}_{0}=-\Gamma^{\nu\mu}{}_{0}\,.\end{aligned}$
(78)

Therefore, we have

$\displaystyle\begin{aligned} \partial_{1}T_{00}+2\alpha
T_{01}=0\,,\end{aligned}$ (79)

i.e.

$\displaystyle\begin{aligned}
\partial_{1}{T}_{00}-{T}_{01}\partial_{1}{g}_{00}=0\,.\end{aligned}$ (80)

Then, we have

$\displaystyle\begin{aligned} \delta
T_{kk}&=\delta\tilde{T}_{00}-\tilde{T}_{01}\partial_{1}\tilde{g}_{00}\delta\tilde{x}^{1}+\tilde{T}_{01}(-\delta\tilde{g}_{00}-\partial_{1}\tilde{g}_{00}\delta\tilde{x}^{1})\\\
&=\delta\tilde{T}_{00}-\tilde{T}_{01}\delta\tilde{g}_{00}\\\
&=\delta\tilde{T}^{1}_{0}=0,\end{aligned}$ (81)

which is the optimal condition on the null hypersurface $H_{2}$.

##  Appendix B Proof 2

In this appendix, we would like to proof that there exists a null hypersurface
linking $\tilde{r}=r_{h}$ in
$O_{1}:\\{\tilde{v},\tilde{r},\tilde{\theta}^{i}\\}$ and $r=r_{h}$ in
$O_{2}:\\{v,r,\theta^{i}\\}$ for an extremal black hole when the optimal
condition is satisfied.

We assume a null hypersurface $\Sigma$ expressed as $\tilde{r}=r_{h}$ at early
time and is the horizon on background. Assuming there is a tangent vector
$\left(\frac{\partial}{\partial u}\right)^{a}$ of $\Sigma$ and the Gaussian
null coordinate of $\Sigma$ is $O_{\text{GN}}:\\{u,z,\tilde{\theta}^{i}\\}$,
then, the expression of $\Sigma$ in $O_{2}$ is

$\displaystyle\begin{aligned} v=u\,,\quad
r=r(u,\tilde{\theta}^{i},\lambda)\,,\quad\theta^{i}=\theta^{i}(u,\tilde{\theta}^{i},\lambda)\,,\end{aligned}$
(82)

where
$r(u,\tilde{\theta}^{i},0)=r_{h}\,,\theta^{i}=\theta^{i}(u,\tilde{\theta}^{i},0)=\tilde{\theta}^{i}$.
We can get the coordinate basis in Gaussian null coordinate system

$\displaystyle\begin{aligned} &k^{a}=\left(\frac{\partial}{\partial
u}\right)^{a}=\left(\frac{\partial}{\partial v}\right)^{a}+\frac{\partial
r}{\partial u}\left(\frac{\partial}{\partial
r}\right)^{a}+\frac{\partial\theta^{i}}{\partial
u}\left(\frac{\partial}{\partial\theta^{i}}\right)^{a}\,,\\\
&(e_{i})^{a}=\left(\frac{\partial}{\partial\tilde{\theta}^{i}}\right)^{a}=\frac{\partial
r}{\partial\tilde{\theta}^{i}}\left(\frac{\partial}{\partial
r}\right)^{a}+\frac{\partial\theta^{j}}{\partial\tilde{\theta}^{i}}\left(\frac{\partial}{\partial\theta^{j}}\right)^{a}\,.\end{aligned}$
(83)

Since the line element in $O_{2}$ is

$\displaystyle\begin{aligned}
ds^{2}=-f(r,\lambda)dv^{2}+2\chi(r,\lambda)dvdr+r^{2}h_{ij}d\theta^{i}d\theta^{j}\,,\end{aligned}$
(84)

and considering the Gaussian null coordinate of a null hypersurface satisfies

$\displaystyle\begin{aligned} k^{a}k_{a}=0\,,\quad
k^{a}(e_{i})_{a}=0\,,\end{aligned}$ (85)

we can find

$\displaystyle\begin{aligned}
&k^{a}k_{a}=-f(r,\lambda)+2\chi(r,\lambda)k^{r}+r^{2}h_{ij}k^{i}k^{j}=0\,,\\\
&k^{a}(e_{i})_{a}=2\chi(r,\lambda)(e_{i})^{r}+r^{2}h_{jl}k^{l}(e_{i})^{j}=0\,.\end{aligned}$
(86)

Considering $k^{r}=k^{i}=(e_{i})^{r}=0$ and $(e_{i})^{j}=\delta_{i}^{j}$ on
the background geometry $\lambda=0$, through the variation of
(LABEL:kappendixb), we can obtain

$\displaystyle\begin{aligned} \delta
k^{r}&=\frac{1}{2\chi(r_{h})}\left.\frac{d}{d\lambda}f\left(r(u,\tilde{\theta}^{i},\lambda\right),\lambda)\right|_{\lambda=0}\\\
&=\frac{1}{2\chi(r_{h})}[\partial_{r}f(r_{h})\delta r+\delta f]\,,\\\ \delta
k^{j}&=-\frac{2}{r_{h}}h^{ij}\chi(r_{h})\delta(e_{i})^{r}\,.\end{aligned}$
(87)

Noting that the variation is done fixing $\\{u,z,\tilde{\theta}^{i}\\}$. For
an extremal black hole, we have $\partial_{r}f(r_{h})=0$. When the first-order
optimal condition is satisfied,i.e. $\delta f=0$, we can find $\delta
k^{r}=0$, which gives

$\displaystyle\frac{\partial\delta r}{\partial u}=0\,.$ (88)

This means that $\delta r$ is a constant along $k^{a}$.

Then, according to the discussion in appendix B, it is obvious that the
optimal condition demands that $\delta T_{ab}k^{a}k^{b}=0$. From the
Raychaudhuri equation on $\Sigma$, we have

$\displaystyle\begin{aligned}
\frac{d\vartheta(\lambda)}{du}=&-\frac{1}{n-2}\vartheta(\lambda)^{2}\\\
&-\sigma_{ab}(\lambda)\sigma^{ab}(\lambda)-T_{ab}(\lambda)k^{a}k^{b}\,,\end{aligned}$
(90)

where

$\displaystyle\begin{aligned}
\vartheta&=\frac{1}{2}\gamma^{ab}\mathcal{L}_{k}\gamma_{ab}\\\
\sigma_{ab}&=\frac{1}{2}\mathcal{L}_{k}\gamma_{ab}-\frac{1}{n-2}\gamma_{ab}\vartheta\,,\end{aligned}$
(91)

Then, considering the fact that $\vartheta(0)=0\,,\sigma_{ab}(0)=0$ on
background, the optimal condition $\delta T_{kk}=0$ implies

$\displaystyle\begin{aligned} \frac{d\delta\vartheta}{du}=0\,,\end{aligned}$
(92)

which implies $\delta\vartheta=0$ on the entire horizon since
$\delta\vartheta=0$ at early-time cross-section $B_{0}$.

Since $\gamma_{ab}$ is independent with $u$ on the background, the variation
of $\vartheta$ can be written as

$\displaystyle\begin{aligned}
\delta\vartheta=\frac{1}{2}\gamma^{ij}\mathcal{L}_{k}\delta\gamma_{ij}\,,\end{aligned}$
(93)

where $\gamma_{ij}$ is the metric component in $O_{\text{GN}}$. Next, we would
like to calculate $\delta\vartheta$. Using the coordinate transformation
relation, we can express $\gamma_{ij}$ using the line element in $O_{2}$

$\displaystyle\begin{aligned} \gamma_{ij}(\tilde{x},\lambda)=\frac{\partial
x^{\mu}(\tilde{x},\lambda)}{\partial\tilde{\theta}^{i}}\frac{\partial
x^{\nu}(\tilde{x},\lambda)}{\partial\tilde{\theta}^{j}}g_{\mu\nu}(x(\tilde{x},\lambda),\lambda)\,,\end{aligned}$
(94)

where we used $\tilde{x}$ to denote the coordinates in $O_{\text{GN}}$. When
$\lambda=0$, there will be

$\displaystyle\begin{aligned}
\gamma_{ij}(\tilde{x})=r^{2}h_{ij}(x(\tilde{x},0),0)=r^{2}h_{ij}(\tilde{x})=\tilde{h}_{ij}(\tilde{x})\end{aligned}$
(95)

and the variation of $\gamma_{ij}$ will be

$\displaystyle\begin{aligned}
\delta\gamma_{ij}(\tilde{x})&=\left.\frac{d}{d\lambda}\left[\frac{\partial
x^{\mu}(\tilde{x},\lambda)}{\partial\tilde{\theta}^{i}}\frac{\partial
x^{\nu}(\tilde{x},\lambda)}{\partial\tilde{\theta}^{j}}g_{\mu\nu}(x(\tilde{x},\lambda),\lambda)\right]\right|_{\lambda=0}\\\
&=2\frac{\partial\delta
x^{\mu}(\tilde{x})}{\partial\tilde{\theta}^{i}}\frac{\partial
x^{\nu}(\tilde{x})}{\partial\tilde{\theta}^{j}}g_{\mu\nu}(\tilde{x})\\\
&\quad+\frac{\partial
x^{\mu}(\tilde{x})}{\partial\tilde{\theta}^{i}}\frac{\partial
x^{\nu}(\tilde{x})}{\partial\tilde{\theta}^{j}}(\delta
g_{\mu\nu}(\tilde{x})+\partial_{\sigma}g_{\mu\nu}(\tilde{x})\delta
x^{\sigma})\\\ &=2\frac{\partial\delta
x^{l}(\tilde{x})}{\partial\tilde{\theta}^{i}}\tilde{h}_{lj}(\tilde{x})+\delta
x^{\sigma}\partial_{\sigma}\tilde{h}_{ij}(\tilde{x})\,,\end{aligned}$ (96)

where we used $\tilde{h}_{ij}$ is independent with $\lambda$. And we can
further obtain

$\displaystyle\begin{aligned}
\mathcal{L}_{k}\delta\gamma_{ij}&=\frac{\partial}{\partial
u}\delta\gamma_{ij}=2\frac{\partial\delta\partial_{0}x^{l}}{\partial\tilde{\theta}^{i}}\tilde{h}_{lj}+\delta(\partial_{0}x^{\sigma})\partial_{\sigma}\tilde{h}_{ij}\\\
&=2\frac{\partial\delta
k^{l}}{\partial\tilde{\theta}^{i}}\tilde{h}_{lj}+\delta
k^{r}\partial_{r}\tilde{h}_{ij}+\delta k^{l}\partial_{l}\tilde{h}_{ij}\\\
&=2\frac{\partial\delta
k^{l}}{\partial\tilde{\theta}^{i}}\tilde{h}_{lj}+\delta
k^{l}\partial_{l}\tilde{h}_{ij}\end{aligned}$ (97)

Then, (93) becomes

$\displaystyle\begin{aligned}
\delta\vartheta&=\frac{1}{2}\gamma^{ij}\mathcal{L}_{k}\delta\gamma_{ij}\\\
&=\frac{\partial\delta
k^{l}}{\partial\tilde{\theta}^{i}}\tilde{h}_{lj}\tilde{h}^{ij}+\frac{1}{2}\delta
k^{l}\tilde{h}^{ij}\partial_{l}\tilde{h}_{ij}\\\
&=2\frac{1}{\sqrt{\gamma}}\partial_{i}(\sqrt{\gamma}\delta
k^{i})\end{aligned}$ (98)

Thus, the optimal condition gives

$\displaystyle\begin{aligned} \delta\vartheta&=D_{i}\delta
k^{i}=D_{i}[-\frac{2}{r_{h}^{2}}h^{ij}\chi(r_{h})\delta(e_{j})^{r}]\\\
&=-\frac{2}{r_{h}^{2}}\chi D_{a}D^{a}\delta r=0\,,\end{aligned}$ (99)

where we use $D_{i}$ to denote the derivative operator adopted to
$\tilde{h}_{ij}$. Therefore, we find that the optimal condition gives
$D_{a}D^{a}\delta r=0$. This is condition also gives

$\displaystyle\begin{aligned} 0&=\int_{B}d^{(n-2)}\theta\sqrt{\gamma}\delta
rD_{a}D^{a}\delta r\\\ &=\int_{B}d^{(n-2)}\theta\sqrt{\gamma}D_{a}(\delta
rD^{a}\delta r)-\int_{B}d^{(n-2)}\theta\sqrt{\gamma}D_{a}\delta rD^{a}\delta
r\,.\end{aligned}$ (100)

Since the sphere with no boundary, the boundary term vanishes. Considering
$D^{a}\delta r$ is spacelike, we can find

$\displaystyle\begin{aligned} D^{a}\delta r=0\,,\end{aligned}$ (101)

which means $\delta r$ is a constant. Thus, the hypersurface at late time
under the first-order approximation is given by $r=constant$. Considering
$\delta\vartheta=0$, the area of the cross-section is a constant under the
variation. Therefore, we find the null hypersurface which is $r=r_{h}$ at late
time.

## References

  * (1) R. Penrose, Gravitational collapse: The role of general relativity, Riv. Nuovo Cim. 1 , 252-276 (1969).
  * (2) R. M. Wald, doi:10.1007/978-94-017-0934-7_5 [arXiv:gr-qc/9710068 [gr-qc]].
  * (3) R. Wald, Gedanken experiments to destroy a black hole, Ann. Phys. (N.Y.) 82, 548 (1974).
  * (4) V. E. Hubeny, Phys. Rev. D 59, 064013 (1999) doi:10.1103/PhysRevD.59.064013 [arXiv:gr-qc/9808043 [gr-qc]].
  * (5) T. Jacobson and T. P. Sotiriou, J. Phys. Conf. Ser. 222, 012041 (2010) doi:10.1088/1742-6596/222/1/012041 [arXiv:1006.1764 [gr-qc]].
  * (6) G. Chirco, S. Liberati and T. P. Sotiriou, Phys. Rev. D 82, 104015 (2010) doi:10.1103/PhysRevD.82.104015 [arXiv:1006.3655 [gr-qc]].
  * (7) A. Saa and R. Santarelli, Phys. Rev. D 84, 027501 (2011) doi:10.1103/PhysRevD.84.027501 [arXiv:1105.3950 [gr-qc]].
  * (8) S. Gao and Y. Zhang, “Destroying extremal Kerr-Newman black holes with test particles,” Phys. Rev. D 87, no.4, 044028 (2013) doi:10.1103/PhysRevD.87.044028 [arXiv:1211.2631 [gr-qc]].
  * (9) P. Zimmerman, I. Vega, E. Poisson and R. Haas, Phys. Rev. D 87, no.4, 041501 (2013) doi:10.1103/PhysRevD.87.041501 [arXiv:1211.3889 [gr-qc]].
  * (10) E. Barausse, V. Cardoso and G. Khanna, Phys. Rev. Lett. 105, 261102 (2010) doi:10.1103/PhysRevLett.105.261102 [arXiv:1008.5159 [gr-qc]]. Barausse:2010ka
  * (11) E. Barausse, V. Cardoso and G. Khanna, Phys. Rev. D 84, 104006 (2011) doi:10.1103/PhysRevD.84.104006 [arXiv:1106.1692 [gr-qc]].
  * (12) M. Colleoni, L. Barack, A. G. Shah and M. van de Meent, Phys. Rev. D 92, no.8, 084044 (2015) doi:10.1103/PhysRevD.92.084044 [arXiv:1508.04031 [gr-qc]].
  * (13) J. Sorce and R. M. Wald, Phys. Rev. D 96, no.10, 104014 (2017) doi:10.1103/PhysRevD.96.104014 [arXiv:1707.05862 [gr-qc]].
  * (14) W. Israel, “Event horizons in static vacuum space-times,” Phys. Rev. 164, 1776-1779 (1967).
  * (15) B. Carter, “Axisymmetric Black Hole Has Only Two Degrees of Freedom,” Phys. Rev. Lett. 26, 331-333 (1971).
  * (16) R. Ruffini and J. A. Wheeler, “Introducing the black hole,” Phys. Today 24, no.1, 30 (1971)
  * (17) M. S. Volkov and D. V. Galtsov, “Non-Abelian Einstein Yang-Mills black holes”, JETP Lett. 50, 346 (1989).
  * (18) P. Bizon, “Colored black holes,” Phys. Rev. Lett. 64, 2844-2847 (1990).
  * (19) B. R. Greene, S. D. Mathur and C. M. O’Neill, “Eluding the no hair conjecture: Black holes in spontaneously broken gauge theories,” Phys. Rev. D 47, 2242-2259 (1993).
  * (20) P. Kanti, N. E. Mavromatos, J. Rizos, K. Tamvakis and E. Winstanley, “Dilatonic black holes in higher curvature string gravity,” Phys. Rev. D 54, 5049-5058 (1996).
  * (21) H. Luckock and I. Moss, “BLACK HOLES HAVE SKYRMION HAIR,” Phys. Lett. B 176, 341-345 (1986).
  * (22) S. Droz, M. Heusler and N. Straumann, “New black hole solutions with hair,” Phys. Lett. B 268, 371-376 (1991).
  * (23) C. Herdeiro and E. Radu, Phys. Rev. Lett. 112, 221101 (2014).
  * (24) C. Herdeiro and E. Radu, Classical Quantum Gravity 32, 144001 (2015).
  * (25) C. A. R. Herdeiro and E. Radu, “Asymptotically flat black holes with scalar hair: a review,” Int. J. Mod. Phys. D 24, no.09, 1542014 (2015).
  * (26) J. P. Hong, M. Suzuki and M. Yamada, “Spherically Symmetric Scalar Hair for Charged Black Holes,” Phys. Rev. Lett. 125, no.11, 111104 (2020).
  * (27) P. A. R. Ade et al. [Planck], Astron. Astrophys. 594, A13 (2016) doi:10.1051/0004-6361/201525830 [arXiv:1502.01589 [astro-ph.CO]].
  * (28) G. Bertone, D. Hooper and J. Silk, Phys. Rept. 405, 279-390 (2005) doi:10.1016/j.physrep.2004.08.031 [arXiv:hep-ph/0404175 [hep-ph]].
  * (29) P. J. E. Peebles and B. Ratra, Rev. Mod. Phys. 75, 559-606 (2003) doi:10.1103/RevModPhys.75.559 [arXiv:astro-ph/0207347 [astro-ph]].
  * (30) S. Hollands and R. M. Wald, Commun. Math. Phys. 321, 629-680 (2013) doi:10.1007/s00220-012-1638-1 [arXiv:1201.0463 [gr-qc]].
  * (31) V. Iyer and R. M. Wald, Phys. Rev. D 50, 846-864 (1994) doi:10.1103/PhysRevD.50.846 [arXiv:gr-qc/9403028 [gr-qc]].