# Probing the scale of grand unification with gravitational waves

Wilfried Buchmuller<EMAIL_ADDRESS>Deutsches Elektronen
Synchrotron DESY, 22607 Hamburg, Germany  Valerie Domcke
<EMAIL_ADDRESS>Deutsches Elektronen Synchrotron DESY, 22607 Hamburg,
Germany  Hitoshi Murayama<EMAIL_ADDRESS>Deutsches Elektronen
Synchrotron DESY, 22607 Hamburg, Germany  Department of Physics, University of
California, Berkeley, CA 94720, USA  Lawrence Berkeley National Laboratory,
Berkeley, CA 94720, USA  Kavli IPMU (WPI), University of Tokyo, Kashiwa
277-8583, Japan  Kai Schmitz<EMAIL_ADDRESS>Theoretical Physics
Department, CERN, 1211 Geneva 23, Switzerland

###### Abstract

The spontaneous breaking of $U(1)_{B-L}$ around the scale of grand unification
can simultaneously account for hybrid inflation, leptogenesis, and neutralino
dark matter, thus resolving three major puzzles of particle physics and
cosmology in a single predictive framework. The $B\\!-\\!L$ phase transition
also results in a network of cosmic strings. If strong and electroweak
interactions are unified in an $SO(10)$ gauge group, containing $U(1)_{B-L}$
as a subgroup, these strings are metastable. In this case, they produce a
stochastic background of gravitational waves that evades current pulsar timing
bounds, but features a flat spectrum with amplitude
$h^{2}\Omega_{\text{GW}}\sim 10^{-8}$ at interferometer frequencies. Ongoing
and future LIGO observations will hence probe the scale of $B\\!-\\!L$
breaking.

††preprint: CERN-TH-2019-215††preprint: DESY 19-210††preprint: IPMU 19-0179

Introduction. The grand unified (GUT) gauge group $SO(10)$ contains
$B\\!-\\!L$, the difference between baryon and lepton number, as a local
symmetry. As shown in Ref. Buchmuller:2012wn , the decay of a false vacuum of
unbroken $B\\!-\\!L$ symmetry is an intriguing and testable mechanism to
generate the initial conditions of the hot early universe (see Ref.
Buchmuller:2013dja  for a review). With $B\\!-\\!L$ broken close to the
unification scale, the false-vacuum phase yields hybrid inflation Linde:1993cn
and ends in tachyonic preheating Felder:2000hj . Decays of the $B\\!-\\!L$
breaking Higgs field and thermal processes produce an abundance of heavy
neutrinos whose decays generate the entropy of the hot early universe, the
baryon asymmetry via leptogenesis Fukugita:1986hr , and dark matter (DM) in
the form of the lightest supersymmetric particle (LSP) Ellis:1983ew .

The cosmological evolution during and after inflation is determined by a small
set of parameters, including the $B\\!-\\!L$ breaking scale and the masses of
the $B\\!-\\!L$ Higgs boson, heavy neutrino, and gravitino. In this Letter, we
revisit the model in Ref. Buchmuller:2012wn  and perform for the first time a
global analysis identifying the viable parameter space that simultaneously
explains the primordial power spectrum of the cosmic microwave background
(CMB), the matter-antimatter asymmetry, and the DM relic density.

In addition, we study the network of cosmic strings (CSs) that is produced
during the $B\\!-\\!L$ phase transition. We embed our model in an $SO(10)$
GUT, which renders the $B\\!-\\!L$ strings metastable due to the
nonperturbative (and hence exponentially suppressed) production of $SO(10)$
monopole–antimonopole pairs. As a consequence, we find that the cosmic-string
network generates a stochastic background of gravitational waves (GWs) that
evades current pulsar timing constraints but that is testable in ongoing and
future GW observations. For a mild hierarchy between the $SO(10)$ and
$U(1)_{B-L}$ breaking scales, LIGO can probe the scale of $B\\!-\\!L$
breaking. The latter result is a direct consequence of the symmetry breaking
pattern and is independent of many details of the model, resulting in a
prediction of the GW spectrum in terms of the two dimensionful parameters of
the cosmic string network. This represents a clear target for current and
future GW experiments.

Supersymmetric $B\\!-\\!L$ model. The example analyzed in Ref.
Buchmuller:2012wn  is based on the gauge group $SU(3)\times SU(2)\times
U(1)_{Y}\times U(1)_{B-L}\equiv G_{\rm SM}\times U(1)_{B-L}$ and
superpotential

$\displaystyle W$ $\displaystyle=W_{\rm
MSSM}+h_{ij}^{\nu}\mathbf{5}_{i}^{*}n_{j}^{c}H_{u}+\frac{1}{\sqrt{2}}\,h_{i}^{n}n_{i}^{c}n_{i}^{c}S_{1}$
$\displaystyle+\lambda\,\Phi\left(\frac{v_{B-L}^{2}}{2}-S_{1}S_{2}\right)+W_{0}\,.$
(1)

Here, $S_{1,2}$ are chiral superfields whose vacuum expectation values break
$B\\!-\\!L$. In unitary gauge, they correspond to the physical $B\\!-\\!L$
Higgs superfield, $S_{1,2}=S/\sqrt{2}$. $\Phi$ is the inflaton superfield,
$n^{c}_{i}$ contain the charge conjugates of the right-handed neutrinos, the
Standard Model (SM) leptons are arranged in the $SU(5)$ multiplets
$\mathbf{5}^{*}=(d^{c},\ell)$ and $\mathbf{10}=(q,u^{c},e^{c})$, and the SM
Higgs fields are contained in $H_{u,d}$. $W_{\text{MSSM}}$ is the MSSM
superpotential,

$\displaystyle
W_{\text{MSSM}}=h^{u}_{ij}\mathbf{10}_{i}\mathbf{10}_{j}H_{u}+h^{d}_{ij}\mathbf{5}^{*}_{i}\mathbf{10}_{j}H_{d}\,,$
(2)

while $W_{0}$ is a constant ensuring the correct zero point of the scalar
potential in supergravity.

The Yukawa couplings are chosen according to the Froggatt–Nielsen (FN) model
Froggatt:1978nt ; Buchmuller:1998zf , which is known to yield a satisfactory
description of quark, charged-lepton, and neutrino masses and mixing angles,
with

$\displaystyle h_{ij}\sim\eta^{Q_{i}+Q_{j}}\ ,\quad\lambda\sim\eta^{Q_{\Phi}}\
,$ (3)

$\eta\simeq 1/\sqrt{300}$ and FN charges listed in Table 1. As usual, the
Yukawa couplings have unspecified $\mathcal{O}(1)$ coefficients that are not
$SU(5)$-symmetric. For simplicity, the FN charges are restricted to
$b=c=d-1=e/2$. This implies the following mass relations for the heavy
Majorana neutrinos $N_{1,2,3}$ and the $B\\!-\\!L$ Higgs field $S$,

$\begin{split}M_{2}&\simeq M_{3}\simeq m_{S}\simeq\eta^{-2}M_{1}\ ,\\\
M_{1}&\sim\eta^{2d}\,v_{B-L}\ ,\quad v_{B-L}\sim\eta^{2a}\frac{v^{2}_{\rm
EW}}{\overline{m}_{\nu}}\ ,\end{split}$ (4)

where $v_{B-L}$ is the $B\\!-\\!L$ breaking scale, $v_{\rm EW}\simeq
174\,\textrm{GeV}$ is the electroweak scale, and
$\overline{m}_{\nu}=\sqrt{m_{2}m_{3}}$. The decays $S\rightarrow
N_{2}N_{2},N_{3}N_{3}$ are thus forbidden, so that leptogenesis is dominated
by the decay chain $S\rightarrow N_{1}N_{1},N_{1}\rightarrow\ell H_{u}$. The
light neutrino masses $m_{1,2,3}$ (with normal hierarchy) are the singular
values of the seesaw mass matrix,

$\displaystyle m_{\nu}=-m_{D}\,M^{-1}m_{D}^{T}\ ,$ (5)

where $m_{D}$ and $M$ denote the Dirac and heavy Majorana neutrino mass
matrices, respectively. As we will discuss in more detail below, successful
inflation requires $v_{B-L}\simeq\textrm{few}\times 10^{15}~{}\text{GeV}$,
corresponding to $a=0$. Hence, top and bottom Yukawa couplings are of the same
order in $\eta$, and $v_{\text{EW}}\simeq\left\langle H_{u}\right\rangle$.
Finally, an important quantity for leptogenesis is the effective light
neutrino mass

$\displaystyle\widetilde{m}_{1}=\frac{(m_{D}^{\dagger}m_{D}^{\vphantom{dagger}})_{11}}{M_{1}}\
.$ (6)

Parametrically, one has $\widetilde{m}_{1}\sim\overline{m}_{\nu}\simeq
0.02~{}\text{eV}$. However, since $\widetilde{m}_{1}$ strongly depends on
unknown $\mathcal{O}(1)$ factors for the Yukawa couplings, $\widetilde{m}_{1}$
is treated as a free parameter in the range
$10^{-5}~{}\text{eV}\leq\widetilde{m}_{1}\leq 10^{-1}~{}\text{eV}$. In
summary, the free parameters of the model can be chosen to be $v_{B-L}$,
$M_{1}$, $\widetilde{m}_{1}$ and $W_{0}$.

Table 1: FN flavor charge assignment (from Ref. Buchmuller:2012wn ). $\psi_{i}$ | $\mathbf{10}_{3}$ | $\mathbf{10}_{2}$ | $\mathbf{10}_{1}$ | $\mathbf{5}^{*}_{3}$ | $\mathbf{5}^{*}_{2}$ | $\mathbf{5}^{*}_{1}$ | $n^{c}_{3}$ | $n^{c}_{2}$ | $n^{c}_{1}$ | $H_{u,d}$ | $S_{1,2}$ | $\Phi$
---|---|---|---|---|---|---|---|---|---|---|---|---
$Q_{i}$ | 0 | 1 | 2 | $a$ | $a$ | $a+1$ | $b$ | $c$ | $d$ | 0 | 0 | $e$

GUT embedding and cosmic strings. The embedding of $G_{\rm SM}\times
U(1)_{B-L}$ into a larger GUT group, spontaneously broken at
$v_{\text{GUT}}>v_{B-L}$, determines the $B\\!-\\!L$ charges of $S_{1,2}$.
Ref. Buchmuller:2012wn  is based on an $SU(5)$ FN model, and the $B\\!-\\!L$
charges are $q_{S}\equiv q_{S_{2}}=-q_{S_{1}}=2$, $q_{n^{c}_{i}}=1$,
$q_{\Phi}=0$. In this case, the final unbroken gauge group is $G_{\rm
SM}\times\mathbb{Z}_{2}$, which results in stable CSs Dror:2019syi . However,
the model can be modified such that $q_{S}=1$. This is achieved by the
replacement in Eq. (Probing the scale of grand unification with gravitational
waves)

$\displaystyle\frac{1}{\sqrt{2}}h_{i}^{n}n_{i}^{c}n_{i}^{c}S_{1}\rightarrow\frac{1}{M_{*}}h_{i}^{n}n_{i}^{c}n_{i}^{c}S_{1}S_{1}\,,$
(7)

where the mass scale $M_{*}$ is larger than the GUT scale. The heavy neutrino
masses are now $M_{i}=h^{n}_{i}v_{B-L}^{2}/M_{*}$. For $M_{*}$ close to
$v_{B-L}$, the mass relations in Eq. (4) can be preserved by changing the
$\mathcal{O}(1)$ factors. In this form, the model can be embedded into an
$SO(10)$ FN model Asaka:2003fp . $S_{1}$ and $S_{2}$ become part of a
$\mathbf{16}^{*}$ and $\mathbf{16}$, respectively, and a further
$\mathbf{16}^{*}+\mathbf{16}$ pair together with a $\mathbf{10}$ is introduced
to achieve an embedding of quarks and leptons that only respects $SU(5)$ but
not $SO(10)$ Nomura:1998gm . In this way, the $SU(5)$ FN charges in Table 1
can be obtained Asaka:2003fp .

With $q_{S}=1$, the final unbroken group is $G_{\rm SM}$, and there can be no
stable strings since the first homotopy group vanishes,
$\Pi_{1}\left(SO(10)/G_{\rm SM}\right)=0$ Dror:2019syi . To understand this,
consider the breaking scheme (see also Leblond:2009fq )

$SO(10)\rightarrow G_{\text{SM}}\times U(1)_{B-L}\rightarrow G_{\text{SM}}\ .$
(8)

The first step produces monopoles, while the second step produces CSs. The
result is a network of CSs and monopoles, which is unstable Vilenkin:1982hm ;
Martin:1996ea ; Martin:1996cp . The more interesting scenario is if cosmic
inflation occurs after GUT but before $B\\!-\\!L$ breaking. In this case, the
initial monopole population is diluted away. The CSs can then only decay via
the Schwinger production of monopole–antimonopole pairs, leading to a
metastable CS network. In this case, the decay rate per string unit length is
Leblond:2009fq ; Monin:2008mp ; Monin:2009ch

$\Gamma_{d}=\frac{\mu}{2\pi}\exp\left(-\pi\kappa\right)\,,$ (9)

with $\kappa=m^{2}/\mu$ denoting the ratio between the monopole mass $m\sim
v_{\text{GUT}}$ and the CS tension $\mu$. In Abelian field theories, one
obtains $\mu=2\pi B(\beta)v^{2}_{B-L}$ with $\beta=\lambda/(8g^{2})$ and
$B(\beta)=2.4/\ln(2/\beta)$ for $\beta<0.01$ Hindmarsh:2011qj . For an
appropriate choice of $v_{B-L}<v_{\rm GUT}$, the CSs are sufficiently long-
lived to give interesting signatures but decay before emitting low-frequency
GWs that are strongly constrained by pulsar timing arrays.

Hybrid inflation, baryogenesis, dark matter. The superpotential term in Eq.
(Probing the scale of grand unification with gravitational waves) causing
spontaneous $B\\!-\\!L$ breaking is precisely the superpotential of F-term
hybrid inflation (FHI) Copeland:1994vg ; Dvali:1994ms . It was widely believed
that FHI could not account for the correct scalar spectral index of the CMB
power spectrum; but the analyses in BasteroGil:2006cm ; Rehman:2009nq ;
Nakayama:2010xf ; Buchmuller:2014epa ; Schmitz:2018nhb  showed that FHI is
viable once the effect of supersymmetry (SUSY) breaking on the inflaton
potential is taken into account. With an (almost) vanishing cosmological
constant, SUSY breaking generates the constant term

$\displaystyle W_{0}=\alpha\,m_{3/2}\,M_{\rm Pl}\,,$ (10)

where $\alpha\sim 1$ encodes the details of SUSY breaking and $m_{3/2}$ is the
gravitino mass. Choosing $\alpha=1$ for definiteness, this adds a term linear
in the inflaton field to the inflaton potential Buchmuller:2000zm ,

$\displaystyle
V_{3/2}(\Phi)=-\lambda\,v_{B-L}^{2}\,m_{3/2}\left(\Phi+\Phi^{*}\right)+\cdots\,.$
(11)

FHI now becomes a two-field model of inflation in the complex $\Phi$-plane
Buchmuller:2014epa . The choice of the inflationary trajectory impacts the CMB
observables and even determines whether a graceful exit from inflation is at
all possible. Trajectories far from the real axis require a significant amount
of tuning. For simplicity, we thus restrict ourselves to trajectories along
the real axis. In this case, the observed values of the amplitude and index of
the scalar power spectrum, $A_{s}^{\rm obs}$ and $n_{s}^{\rm obs}$, eliminate
two of the three parameters $v_{B-L}$, $\lambda$, and $m_{3/2}$, where we
recall that we can exchange $\lambda\leftrightarrow M_{1}$ using the FN
relations (4).

We will now demonstrate that the parameter values required for FHI are
consistent with leptogenesis and neutralino DM. First, we define the reheating
temperature $T_{\text{rh}}$ in terms of the total and radiative energy
densities,

$\displaystyle
T=T_{\text{rh}}\quad\Leftrightarrow\quad\rho_{\text{tot}}=2\,\rho_{\text{rad}}\,.$
(12)

$T_{\text{rh}}$ is determined by the decay widths of $S$ and $N_{1}$, which
can be expressed in terms of the masses defined above,

$\displaystyle\Gamma^{0}_{S}$
$\displaystyle=\frac{m_{S}}{8\pi}\left(\frac{M_{1}}{v_{B-L}}\right)^{2}\left(1-\frac{4M^{2}_{1}}{m^{2}_{S}}\right)^{1/2}\,,$
(13) $\displaystyle\Gamma^{0}_{N_{1}}$
$\displaystyle\simeq\frac{1}{4\pi}\frac{\widetilde{m}_{1}}{\overline{m}_{\nu}}\frac{M_{1}^{2}}{v_{B-L}}\,.$

The neutrinos $N_{1}$ produced in $S$ decays are relativistic, which is
accounted for by the averaged Lorentz factor $\gamma$,

$\displaystyle\Gamma^{S}_{N_{1}}=\gamma^{-1}\,\Gamma^{0}_{N_{1}}\
,\quad\gamma^{-1}=\left\langle\frac{M_{1}}{E_{N_{1}}}\right\rangle\ .$ (14)

Equating the decay rates $\Gamma_{S}\equiv\Gamma^{0}_{S}$ and
$\Gamma_{N_{1}}\equiv\Gamma^{S}_{N_{1}}$ with the Hubble parameter $H$ defines
the $S$ and $N_{1}$ decay temperatures,

$\displaystyle
H\left(T_{X}\right)=\Gamma_{X}:\>\>T_{X}=\left(\frac{90}{\alpha_{X}\pi^{2}g_{*}}\right)^{1/4}\sqrt{\Gamma_{X}M_{\rm
Pl}}\,,$ (15)

where $X=\\{S,N_{1}\\}$ and
$\alpha^{-1}_{X}\equiv\rho_{\text{rad}}/\rho_{\text{tot}}|_{T_{X}}$. For
$T_{\rm rh}$, we then find (see Appendix C in Buchmuller:2013lra  for details)

$\displaystyle T_{\mathrm{rh}}$ $\displaystyle\simeq
0.85\times\mathrm{min}\\{T_{N_{1}},T_{S}\\}\,.$ (16)

Figure 1: Viable parameter space (green) for hybrid inflation, leptogenesis,
neutralino DM, and big bang nucleosynthesis. Hybrid inflation and the dynamics
of reheating correlate the parameters $v_{B-L}$, $T_{\text{rh}}$, $m_{3/2}$
and $\widetilde{m}_{1}$ (black curves). Successful leptogenesis occurs outside
the gray-shaded region. Neutralino DM is viable in the green region,
corresponding to a higgsino (wino) with mass $100\leq m_{\rm
LSP}/\textrm{GeV}\leq 1060$ (2680).

For any given values of $\widetilde{m}_{1}$ and $m_{3/2}$, successful hybrid
inflation selects a point in the $v_{B-L}$ – $T_{\text{rh}}$ plane, such that
the values of all parameters are fixed [up to an $\mathcal{O}(1)$ uncertainty
due to the constant $\alpha$ in Eq. (10)]. This is illustrated in Figure 1.
The turn-over at large gravitino masses and high reheating temperatures
reflects the quartic supergravity coupling in the inflaton potential becoming
important. Even larger gravitino masses entail field excursions of ${\cal
O}(M_{\text{Pl}})$ and thus are sensitive to further, model-dependent
supergravity corrections. For large values of $\widetilde{m}_{1}$, the
reheating temperature is set by $S$ decays, such that it becomes independent
of $\widetilde{m}_{1}$. This is reflected in the merging of different
$\widetilde{m}_{1}$ contour lines at the upper edge of the band depicted in
Fig. 1.

The decay of the heavy neutrino $N_{1}$ is responsible for nonthermal
leptogenesis Lazarides:1991wu ; Asaka:1999yd , which dominates in our model
over thermal leptogenesis for $T_{\text{rh}}\lesssim 10^{10}$ GeV
Buchmuller:2012wn . The resulting baryon asymmetry can be estimated as

$\displaystyle\eta_{B}\simeq\eta_{B}^{\rm nt}\simeq C_{\rm
sph}\,\frac{g_{*,s}^{0}}{g_{*,s}}\frac{g_{*,\rho}}{g_{\gamma}}\frac{\pi^{4}}{30\,\zeta(3)}\,\varepsilon_{1}\frac{T_{\rm
rh}}{\gamma\,M_{1}}\,,$ (17)

with $C_{\text{sph}}=8/23$ denoting the sphaleron conversion factor,
$g_{*,\rho}=g_{*,s}=915/4$, $g_{*,s}^{0}=43/11$, and $g_{\gamma}=2$ counting
(effective) numbers of relativistic degrees of freedom and
$\varepsilon_{1}\lesssim 2\times 10^{-6}M_{1}/\left(10^{10}\text{GeV}\right)$
Hamaguchi:2001gw ; Davidson:2002qv ; Buchmuller:2005eh  parameterizing the
$CP$ asymmetry in $N_{1}$ decays. This agrees with the result obtained solving
the corresponding Boltzmann equations within a factor of two Buchmuller:2012wn
. The gray shading in Fig. 1 indicates the region where leptogenesis falls
short of explaining the observed baryon asymmetry.

Gravitino masses of $\mathcal{O}\left(1\right)\,\textrm{TeV}$ or larger point
to a neutralino LSP, which is produced thermally as well as nonthermally in
gravitino decays Buchmuller:2012bt . Gravitinos are in turn generated in
decays of the $B\\!-\\!L$ Higgs field at a rate

$\Gamma_{S}^{3/2}=\frac{x}{32\pi}\left(\frac{\left<S\right>}{M_{\text{Pl}}}\right)^{2}\frac{m_{S}^{3}}{M_{\text{Pl}}^{2}}\,,$
(18)

as well as from the thermal bath (for a discussion and references, see
Jeong:2012en ). The parameter $x$ in Eq. (18) encodes details of the
unspecified SUSY-breaking sector, and for definiteness, we will assume all
gaugino masses to be significantly lighter than $m_{3/2}$. Taking into account
that gravitinos must decay early enough to preserve big bang nucleosynthesis
(BBN) Kawasaki:2017bqm  and $m_{\text{LSP}}\gtrsim 100$ GeV PhysRevD.98.030001
, a higgsino or wino LSP can account for the observed DM relic density in the
green-shaded region of Fig. 1. It is highly nontrivial that neutralino DM and
leptogenesis can be successfully realized in the same parameter region. In
summary, the viable parameter region of our model is given by $v_{B-L}\simeq
3.0\cdots 5.8\times 10^{15}\,\textrm{GeV}$, $m_{3/2}\simeq
10\,\textrm{TeV}\cdots 10\,\textrm{PeV}$, and $x\lesssim 0.4$.

Gravitational waves. The network of CSs formed during the $B\\!-\\!L$ phase
transition acts as a source of GWs Buchmuller:2013lra . The dominant
contribution is expected to come from long-lived, sub-horizon CS loops that
eventually decay emitting gravitational radiation. Modeling the evolution and
GW emission of a CS network is a challenging task, resulting in several
competing models in the literature. Among others, major questions are the
ansatz for CSs (field theory versus Nambu–Goto), the correct loop number
density and the most efficient GW production channel (see Ref. Auclair:2019wcv
and references therein for a comprehensive review). Moreover, for metastable
strings, the GW production from fast-moving monopoles requires further
investigation Vilenkin:1982hm ; Leblond:2009fq . For concreteness, we will
base our analysis on the BOS model Blanco-Pillado:2013qja , but the arguments
presented here can easily be applied to any other CS model.

Figure 2: GW spectrum for $G\mu=2\times 10^{-7}$. Different values of
$\sqrt{\kappa}$ are indicated in different colors; the blue curve corresponds
to a CS network surviving until today. The dot-dashed lines depict the
analytical estimate (25). The (lighter) gray-shaded areas indicate the
sensitivities of (planned) experiments SKA Smits:2008cf , LISA Audley:2017drz
, LIGO LIGOScientific:2019vic  and ET ET , the crosses within the SKA band
indicate constraints by the IPTA Verbiest:2016vem .

The present-day GW spectrum can be expressed as Auclair:2019wcv

$\displaystyle\Omega_{\text{GW}}(f)=\frac{\partial\rho_{\text{GW}}(f)}{\rho_{c}\partial\ln
f}=\frac{8\pi f(G\mu)^{2}}{3H_{0}^{2}}\sum_{n=1}^{\infty}C_{n}(f)\,P_{n}\,,$
(19)

where $\rho_{\text{GW}}$ denotes the GW energy density, $\rho_{c}$ is the
critical energy density of the universe, $G\mu$ denotes the dimensionless
string tension, $H_{0}=100\,h\,\textrm{km}/\textrm{s}/\textrm{Mpc}$ is today’s
Hubble parameter, $P_{n}\simeq 50/\zeta[4/3]\,n^{-4/3}$ is the power spectrum
of GWs emitted by the $n^{\rm th}$ harmonic of a CS loop, and $C_{n}(f)$
indicates the number of loops emitting GWs that are observed at a given
frequency $f$,

$\displaystyle
C_{n}(f)=\frac{2n}{f^{2}}\int_{z_{\text{min}}}^{z_{\text{max}}}dz\>\frac{\mathcal{N}\left(\ell\left(z\right),\,t\left(z\right)\right)}{H\left(z\right)(1+z)^{6}}\,,$
(20)

which is a function of the number density of CS loops $\mathcal{N}(\ell,t)$,
with $\ell=2n/((1+z)f)$, selecting the loops that contribute to the spectrum
at frequency $f$ today. The loop number density can be estimated analytically
and improved with input from numerical simulations. In particular, for loops
generated and decaying in the radiation-dominated (RD) era, it can be
estimated as Blanco-Pillado:2013qja ; Auclair:2019wcv

$\displaystyle\mathcal{N}_{r}(\ell,t)$
$\displaystyle=\frac{0.18}{t^{3/2}(\ell+\Gamma G\mu t)^{5/2}}\,,$ (21)

where $\Gamma\simeq 50$ parametrizes the CS decay rate into GWs,
$\dot{\ell}=-\Gamma G\mu$. The integration range in Eq. (20) accounts for the
lifetime of the CS network, from the formation at $z_{\text{max}}\simeq
T_{\text{rh}}/(2.7~{}\text{K})$ until their decay when the decay rate of a
string loop with average length $\bar{\ell}\simeq 13\,G\mu/H$ equals the
Hubble rate, $\bar{\ell}\,\Gamma_{d}=H$. The latter condition yields

$z_{\text{min}}=\left(\frac{70}{H_{0}}\right)^{1/2}\left(\Gamma\;\Gamma_{d}\;G\mu\right)^{1/4}\,,$
(22)

which coincides with the estimate found in Leblond:2009fq .

Figure 3: Constraints on the CS parameter space from GW searches. The orange,
red, and purple regions are excluded by the existing bounds from EPTA
Shannon:2015ect , LIGO O2 LIGOScientific:2019vic  and PLANCK Ade:2015xua ,
respectively. The shaded green and blue regions indicate the prospective reach
of LIGO at design sensitivity and of LISA. The hatched region indicates the
viable parameter space from Fig. 1.

For a sufficiently long-lived CS network, the spectrum features a plateau over
many orders of magnitude in frequency, accounting for the GW production from
CS loops during the RD era. The amplitude of this plateau is Auclair:2019wcv

$\Omega_{\text{GW}}^{\text{plateau}}\simeq
8.04\,\Omega_{r}\left(\frac{G\mu}{\Gamma}\right)^{1/2}\,.$ (23)

A finite lifetime of the CS network leads to a suppression of the GW spectrum
at low frequencies. Focusing for simplicity only on the RD era, we note that
the scaling of the integrand of Eq. (20) with $z$ changes at
$z_{\ell}=\Gamma/2\,G\mu\,f\,t_{\text{ref}}\,z_{\text{ref}}^{2}$, with
$t_{\text{ref}}=t(z_{\text{ref}})$ denoting an arbitrary reference time during
RD. One can verify that, for $z_{\text{min}}\ll z_{\ell}$ (corresponding to
large frequencies), $\Omega_{\text{GW}}(f)$ becomes scale-invariant, whereas
for $z_{\text{min}}\gg z_{\ell}$, we find $\Omega_{\text{GW}}\propto f^{3/2}$.
With this, we can estimate the turn-over point $f_{*}$ of the spectrum by
solving $z_{\text{min}}=z_{\ell}$ for $f$,

$\displaystyle f_{*}\simeq 3.0\times
10^{14}\,\text{Hz}\>e^{-\pi\kappa/4}\left(\frac{10^{-7}}{G\mu}\right)^{1/2}\,.$
(24)

In summary, we arrive at a simple estimate for the GW spectrum from a
metastable CS network,

$\displaystyle\Omega_{\text{GW}}(f)=\Omega_{\text{GW}}^{\text{plateau}}\;\text{min}\left[(f/f_{*})^{3/2},1\right]\,.$
(25)

Fig. 2 shows the GW spectrum obtained by numerically evaluating Eq. (19) (see
Ref. Auclair:2019wcv  for details) as well as the analytical estimate Eq.
(25). The shaded regions indicate the power-law-integrated sensitivity curves
of current and planned experiments Thrane:2013oya . We see that, for
$G\mu=2\times 10^{-7}$, the constraint from the European Pulsar Timing Array
(EPTA) Shannon:2015ect  enforces $\sqrt{\kappa}\lesssim 8$. The analytical
estimate accurately explains the plateau value and the $f^{3/2}$ drop-off at
low frequencies, showing overall good agreement with the numerical result.

Comparing the predicted GW spectrum (25) with the observational bounds and
prospects depicted in Fig. 2, we can map out the regions in the $\kappa$ –
$G\mu$ parameter plane that are already excluded or that will be probed in the
near future, see Fig. 3. In particular, cosmological $B\\!-\\!L$ breaking with
stable CSs is excluded, as are GUT monopole masses above $m=5.4\times 10^{16}$
GeV. Assuming a mild hierarchy between the GUT and $B\\!-\\!L$ scales,
$m/v_{B-L}\gtrsim 6$, the entire remaining parameter space is testable with
LIGO.

Conclusions. We have presented a minimal extension of the SM with $U(1)_{B-L}$
symmetry that simultaneously explains inflation, leptogenesis, and neutralino
DM. Gravitinos are unstable and heavier than $10~{}\text{TeV}$. Remarkably,
the viable parameter space of our model automatically implies a large
stochastic signal in GWs that is well within the reach of LIGO’s design
sensitivity. It would be interesting to relax some of our model-building
assumptions, such as the underlying flavor model, in future work. However, the
general philosophy behind our model — inflation ending in a GUT-scale phase
transition in combination with leptogenesis and dark matter in a SUSY
extension of the SM — provides a testable framework for the physics of the
early universe. A characteristic feature of this framework is a stochastic
background of gravitational waves emitted by metastable cosmic strings.

###### Acknowledgements.

Acknowledgements. We thank Ryusuke Jinno for helpful discussions on the GW
production from metastable cosmic strings. This project has received funding
from the Deutsche Forschungsgemeinschaft under Germany’s Excellence Strategy –
EXC 2121 “Quantum Universe” – 390833306 (V. D.) and from the European Union’s
Horizon 2020 Research and Innovation Programme under grant agreement number
796961, “AxiBAU” (K. S.). The work of H. M. was supported by the NSF grant
PHY-1915314, by the U.S. DOE Contract DE-AC02-05CH11231, by the JSPS Grant-in-
Aid for Scientific Research JP17K05409, MEXT Grant-in-Aid for Scientific
Research on Innovative Areas JP15H05887, JP15K21733, by WPI, MEXT, Japan, and
Hamamatsu Photonics.

## References

  * (1) W. Buchmuller, V. Domcke, and K. Schmitz, Spontaneous $B\\!-\\!L$ Breaking as the Origin of the Hot Early Universe, Nucl. Phys. B862 (2012) 587–632, [arXiv:1202.6679].
  * (2) W. Buchmuller, V. Domcke, K. Kamada, and K. Schmitz, A Minimal Supersymmetric Model of Particle Physics and the Early Universe, in Cosmology and Particle Physics beyond Standard Models (L. Ávarez-Gaumé, G. S. Djordjević, and D. Stojković, eds.), CERN-Proceedings-2014-001.  arXiv:1309.7788.
  * (3) A. D. Linde, Hybrid inflation, Phys. Rev. D49 (1994) 748–754, [astro-ph/9307002].
  * (4) G. N. Felder, J. Garcia-Bellido, P. B. Greene, L. Kofman, A. D. Linde, and I. Tkachev, Dynamics of symmetry breaking and tachyonic preheating, Phys. Rev. Lett. 87 (2001) 011601, [hep-ph/0012142].
  * (5) M. Fukugita and T. Yanagida, Baryogenesis Without Grand Unification, Phys. Lett. B174 (1986) 45–47.
  * (6) J. R. Ellis, J. S. Hagelin, D. V. Nanopoulos, K. A. Olive, and M. Srednicki, Supersymmetric Relics from the Big Bang, Nucl. Phys. B238 (1984) 453–476. [,223(1983)].
  * (7) C. D. Froggatt and H. B. Nielsen, Hierarchy of Quark Masses, Cabibbo Angles and CP Violation, Nucl. Phys. B147 (1979) 277–298.
  * (8) W. Buchmuller and T. Yanagida, Quark lepton mass hierarchies and the baryon asymmetry, Phys. Lett. B445 (1999) 399–402, [hep-ph/9810308].
  * (9) J. A. Dror, T. Hiramatsu, K. Kohri, H. Murayama, and G. White, Testing the Seesaw Mechanism and Leptogenesis with Gravitational Waves, Phys. Rev. Lett. 124 (2020) 041804, [arXiv:1908.03227].
  * (10) T. Asaka, Lopsided mass matrices and leptogenesis in SO(10) GUT, Phys. Lett. B562 (2003) 291–298, [hep-ph/0304124].
  * (11) Y. Nomura and T. Yanagida, Bimaximal neutrino mixing in SO(10)(GUT), Phys. Rev. D59 (1999) 017303, [hep-ph/9807325].
  * (12) L. Leblond, B. Shlaer, and X. Siemens, Gravitational Waves from Broken Cosmic Strings: The Bursts and the Beads, Phys. Rev. D79 (2009) 123519, [arXiv:0903.4686].
  * (13) A. Vilenkin, Cosmological evolution of monopoles connected by strings, Nucl. Phys. B196 (1982) 240–258.
  * (14) X. Martin and A. Vilenkin, Gravitational wave background from hybrid topological defects, Phys. Rev. Lett. 77 (1996) 2879–2882, [astro-ph/9606022].
  * (15) X. Martin and A. Vilenkin, Gravitational radiation from monopoles connected by strings, Phys. Rev. D55 (1997) 6054–6060, [gr-qc/9612008].
  * (16) A. Monin and M. B. Voloshin, The Spontaneous breaking of a metastable string, Phys. Rev. D78 (2008) 065048, [arXiv:0808.1693].
  * (17) A. Monin and M. B. Voloshin, Destruction of a metastable string by particle collisions, Phys. Atom. Nucl. 73 (2010) 703–710, [arXiv:0902.0407].
  * (18) M. Hindmarsh, Signals of Inflationary Models with Cosmic Strings, Prog. Theor. Phys. Suppl. 190 (2011) 197–228, [arXiv:1106.0391].
  * (19) E. J. Copeland, A. R. Liddle, D. H. Lyth, E. D. Stewart, and D. Wands, False vacuum inflation with Einstein gravity, Phys. Rev. D49 (1994) 6410–6433, [astro-ph/9401011].
  * (20) G. R. Dvali, Q. Shafi, and R. K. Schaefer, Large scale structure and supersymmetric inflation without fine tuning, Phys. Rev. Lett. 73 (1994) 1886–1889, [hep-ph/9406319].
  * (21) M. Bastero-Gil, S. F. King, and Q. Shafi, Supersymmetric Hybrid Inflation with Non-Minimal Kahler potential, Phys. Lett. B651 (2007) 345–351, [hep-ph/0604198].
  * (22) M. U. Rehman, Q. Shafi, and J. R. Wickman, Supersymmetric Hybrid Inflation Redux, Phys. Lett. B683 (2010) 191–195, [arXiv:0908.3896].
  * (23) K. Nakayama, F. Takahashi, and T. T. Yanagida, Constraint on the gravitino mass in hybrid inflation, JCAP 1012 (2010) 010, [arXiv:1007.5152].
  * (24) W. Buchmuller, V. Domcke, K. Kamada, and K. Schmitz, Hybrid Inflation in the Complex Plane, JCAP 1407 (2014) 054, [arXiv:1404.1832].
  * (25) K. Schmitz and T. T. Yanagida, Axion Isocurvature Perturbations in Low-Scale Models of Hybrid Inflation, Phys. Rev. D98 (2018), no. 7 075003, [arXiv:1806.06056].
  * (26) W. Buchmuller, L. Covi, and D. Delepine, Inflation and supersymmetry breaking, Phys. Lett. B491 (2000) 183–189, [hep-ph/0006168].
  * (27) W. Buchmuller, V. Domcke, K. Kamada, and K. Schmitz, The Gravitational Wave Spectrum from Cosmological $B-L$ Breaking, JCAP 1310 (2013) 003, [arXiv:1305.3392].
  * (28) G. Lazarides and Q. Shafi, Origin of matter in the inflationary cosmology, Phys. Lett. B258 (1991) 305–309.
  * (29) T. Asaka, K. Hamaguchi, M. Kawasaki, and T. Yanagida, Leptogenesis in inflaton decay, Phys. Lett. B464 (1999) 12–18, [hep-ph/9906366].
  * (30) K. Hamaguchi, H. Murayama, and T. Yanagida, Leptogenesis from N dominated early universe, Phys. Rev. D65 (2002) 043512, [hep-ph/0109030].
  * (31) S. Davidson and A. Ibarra, A Lower bound on the right-handed neutrino mass from leptogenesis, Phys. Lett. B535 (2002) 25–32, [hep-ph/0202239].
  * (32) W. Buchmuller, R. D. Peccei, and T. Yanagida, Leptogenesis as the origin of matter, Ann. Rev. Nucl. Part. Sci. 55 (2005) 311–355, [hep-ph/0502169].
  * (33) W. Buchmuller, V. Domcke, and K. Schmitz, WIMP Dark Matter from Gravitino Decays and Leptogenesis, Phys. Lett. B713 (2012) 63–67, [arXiv:1203.0285].
  * (34) K. S. Jeong and F. Takahashi, A Gravitino-rich Universe, JHEP 01 (2013) 173, [arXiv:1210.4077].
  * (35) M. Kawasaki, K. Kohri, T. Moroi, and Y. Takaesu, Revisiting Big-Bang Nucleosynthesis Constraints on Long-Lived Decaying Particles, Phys. Rev. D97 (2018), no. 2 023502, [arXiv:1709.01211].
  * (36) Particle Data Group Collaboration, M. Tanabashi, K. Hagiwara, K. Hikasa, K. Nakamura, Y. Sumino, and e. a. Takahashi, Review of particle physics, Phys. Rev. D 98 (Aug, 2018) 030001.
  * (37) P. Auclair et al., Probing the gravitational wave background from cosmic strings with LISA, JCAP 04 (2020) 034, [arXiv:1909.00819].
  * (38) J. J. Blanco-Pillado, K. D. Olum, and B. Shlaer, The number of cosmic string loops, Phys. Rev. D89 (2014), no. 2 023512, [arXiv:1309.6637].
  * (39) R. Smits, M. Kramer, B. Stappers, D. R. Lorimer, J. Cordes, and A. Faulkner, Pulsar searches and timing with the square kilometre array, Astron. Astrophys. 493 (2009) 1161–1170, [arXiv:0811.0211].
  * (40) LISA Collaboration, P. Amaro-Seoane et al., Laser Interferometer Space Antenna, arXiv:1702.00786.
  * (41) LIGO Scientific, Virgo Collaboration, B. P. Abbott et al., Search for the isotropic stochastic background using data from Advanced LIGO’s second observing run, Phys. Rev. D100 (2019), no. 6 061101, [arXiv:1903.02886].
  * (42) Einstein Telescope Collaboration http://www.et-gw.eu/index.php/etsensitivities.
  * (43) J. P. W. Verbiest et al., The International Pulsar Timing Array: First Data Release, Mon. Not. Roy. Astron. Soc. 458 (2016), no. 2 1267–1288, [arXiv:1602.03640].
  * (44) R. M. Shannon et al., Gravitational waves from binary supermassive black holes missing in pulsar observations, Science 349 (2015), no. 6255 1522–1525, [arXiv:1509.07320].
  * (45) Planck Collaboration, P. A. R. Ade et al., Planck 2015 results. XIII. Cosmological parameters, Astron. Astrophys. 594 (2016) A13, [arXiv:1502.01589].
  * (46) E. Thrane and J. D. Romano, Sensitivity curves for searches for gravitational-wave backgrounds, Phys. Rev. D88 (2013), no. 12 124032, [arXiv:1310.5300].


aainstitutetext: College of Education for the Future, Beijing Normal
University, Zhuhai 519087, Chinabbinstitutetext: Department of Physics,
Beijing Normal University, Beijing, 100875, China

# Gedanken experiments at high-order approximation: Kerr black hole cannot be
overspun

Aofei Sangab and Jie Jiangab111Corresponding author.
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Sorce and Wald proposed a new version of gedanken experiments to examine the
weak cosmic censorship conjecture (WCCC) in Kerr-Newmann black holes. However,
their discussion only includes the second-order approximation of perturbation
and there exists an optimal condition such that the validity of the WCCC is
determined by the higher-order approximations. Therefore, in this paper, we
extended their discussions into the high-order approximations to study the
WCCC in a nearly extremal Kerr black hole. After assuming that the spacetime
satisfies the stability condition and the perturbation matter fields satisfy
the null energy condition, based on the Noether charge method by Iyer and
Wald, we completely calculate the first four order perturbation inequalities
and discuss the corresponding gedanken experiment to overspin the Kerr black
hole. As a result, we find that the nearly extremal Kerr black holes cannot be
destroyed under the fourth-order approximation of perturbation. Then, by using
the mathematical induction, we strictly prove the $n$th order perturbation
inequality when the first $(n-1)$ order perturbation inequalities are
saturated. Using these results, we discuss the first $100$ order approximation
of the gedanken experiments and find that the WCCC in Kerr black hole is valid
under the higher-order approximation of perturbation. Our investigation
implies that the WCCC might be strictly satisfied in Kerr black holes under
the perturbation level.

###### Keywords:

General relativity, Kerr black hole, weak cosmic censorship conjecture

##  1 Introduction

General relativity, which describes gravitational interactions in terms of
space-time curvature, has made great achievements in explaining many kinds of
gravitational phenomena at various scales. Despite its success, however,
general relativity still faces some theoretical challenges. One problem with
general relativity is the existence of spacetime singularities in the solution
of Einstein’s equations. The Singularity theoremPenrose:1969pc  guarantees
that the formation of singularities is a general feature of gravitational
collapse. At the singularity, the physical quantities diverge. To ensure that
there are no naked singularities in spacetime, Penrose proposed the weak
cosmic censorship conjecturePenrose:1969pc ; Wald:1997wa .

Since the WCCC was proposed, it becomes one of the important conjectures in
classical general relativity and still lacks general proof till now. In 1974,
Wald first proposed a gedanken experimentWald:1974wl  to test the validity of
WCCC. They assumed a test particle absorbed by an extremal Kerr-Newman black
hole and showed that the black hole cannot be over-spun or over-charged
through this process under first-order approximation. After that, through
dropping the test particle, HubenyHubeny:1998ga  considered the second-order
approximation and found some possibilities of destruction. Then, violations of
this type are found in many followup worksdeFelice:2001wj ; Jacobson:2010iu ;
Chirco:2010rq ; Saa:2011wq ; Gao:2012ca . However, Hubeny’s method neglected
some second-order effects such as the self-force and self-energy effects. In
2017, Sorce and WaldSorce:2017dst  proposed a new version of the gedanken
experiments. In this version, they considered a fully dynamic process of some
matter field falling into a nearly extremal Kerr black hole. By assuming the
matter fields satisfy the null energy condition, they used the Noether charge
methodIyer:1994ys  to derive the first- and second-order perturbation
inequalities

$\displaystyle\begin{aligned} &\delta M-\Omega_{H}\delta J-\Phi_{H}\delta
Q\geq 0\,,\\\
&\delta^{2}M-\Omega_{H}\delta^{2}J-\Phi_{H}\delta^{2}Q\geq-\frac{\kappa}{8\pi}\delta^{2}A^{\text{KN}}_{B}\end{aligned}$
(1)

to constrain the mass $M(\lambda)$, angular momentum $J(\lambda)$ and electric
charge $Q(\lambda)$ of the geometry after perturbation. Using these results,
they showed that the existence condition of the event horizon
$h(\lambda)=M^{2}(\lambda)-Q^{2}(\lambda)-J^{2}(\lambda)/M^{2}(\lambda)$ can
be reduced to

$\displaystyle\begin{aligned} h(\lambda)\geq\left(\frac{(J^{2}-M^{4})Q\delta
Q-2JM^{2}\delta J}{M(M^{4}+J^{2})}\lambda+M\varepsilon\right)^{2}\geq
0\end{aligned}$ (2)

under the second-order approximation, in which $\varepsilon=r_{h}/M-1$ denotes
the deviation to the extremal black hole. This result indicates that no
violations can occur under the second-order approximation. Since then, this
method has been extended to many other stationary black holes to test the
validity of WCCC.

However, the discussion from Sorce and Wald only includes the first- and
second-order approximation of the perturbation. According to their result (2),
we can see that there exists a second-order optimal condition in which the
first two order perturbation inequalities are saturated and

$\displaystyle\begin{aligned} \delta
Q=\frac{M^{2}[(M^{4}+J^{2})\varepsilon-2J\delta
J\lambda]}{(M^{4}-J^{2})Q\lambda}\,,\end{aligned}$ (3)

such that $h(\lambda)=0$ under the second-order approximation and we cannot
judge the validity of WCCC in this approximation. For strictness, it is
necessary to consider high order approximation of the gedanken experiments. In
2020, one of our authors studied the high-order gedanken experiments for a
nearly extremal Reissner-Nordström black holes under the spherical
perturbationWang:2020vpn , and found that the black hole cannot be overcharged
in this process. However, one of the drawbacks of their method is that it is
only suitable for the spherical symmetric perturbation process as well as the
spherical spacetime. Astronomical observations show that most stars in the
universe carry large angular momentum but small electric charge. Moreover, the
real physical process in our universe is not spherical. Therefore, in this
paper, we would like to extend the Sorce-Wald investigation into the high-
order approximation in a nearly extremal Kerr black hole with any perturbation
matter fields which satisfy the null energy condition and examine the WCCC in
this process.

The remainder of this paper is organized as follows. In the next section, we
review the Noether charge and variational identity in General Relativity. In
section 3, we introduce the geometry of the Kerr black hole under the
perturbation matter fields. In section 4, we perform the gedanken experiments
proposed by Sorce and WaldSorce:2017dst  to examine the WCCC under the high-
order approximation of perturbation when the matter fields satisfy the null
energy condition. Finally, we give a brief conclusion in section 5.

##  2 Variational identity in Einstein gravity

In this section, we review the Noether charge and variational identity in
Einstein gravity with the Lagrangian four-form

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

Considering an one-parameter family with the parameter $\lambda$, the $k$th-
order variation of the field $\phi$ is defined by

$\displaystyle\delta^{k}\phi=\left.\frac{d^{k}\phi(\lambda)}{d\lambda^{k}}\right|_{\lambda=0}\,,$
(5)

i.e., its $k$th derivative evaluated at $\lambda=0$. Taking the derivative of
the Lagrangian, we have

$\displaystyle\frac{d\bm{L}_{\text{grav}}}{d\lambda}=\bm{E}_{g}^{ab}\frac{dg_{ab}}{d\lambda}+d\bm{\Theta}\left(g,\frac{dg}{d\lambda}\right)\,,$
(6)

where

$\displaystyle\begin{aligned}
\bm{E}_{g}^{ab}&=-\frac{\bm{\epsilon}}{16\pi}G^{ab}\,,\\\
\bm{\Theta}_{abc}\left(g,\frac{dg}{d\lambda}\right)&=\frac{1}{16\pi}\bm{\epsilon}_{dabc}g^{de}g^{fg}\left(\nabla_{g}\frac{dg_{ef}}{d\lambda}-\nabla_{e}\frac{dg_{fg}}{d\lambda}\right)\end{aligned}$
(7)

are the on-shell equation of motion and the symplectic potential three-form
separately. For any two-parameter family with parameters $\lambda_{1}$ and
$\lambda_{2}$, we can define a symplectic current three form

$\displaystyle\bm{\omega}_{abc}\left(g;\frac{\partial
g}{\partial\lambda_{1}},\frac{\partial
g}{\partial\lambda_{2}}\right)=\frac{1}{16\pi}\bm{\epsilon}_{dabc}w^{d}\,,$
(8)

in which

$\displaystyle w^{a}=P^{abcdef}\left(\frac{\partial
g_{bc}}{\partial\lambda_{2}}\nabla_{d}\frac{\partial
g_{ef}}{\partial\lambda_{1}}-\frac{\partial
g_{bc}}{\partial\lambda_{1}}\nabla_{d}\frac{\partial
g_{ef}}{\partial\lambda_{2}}\right)$ (9)

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}$
(10)

The Noether current three-form $\bm{J}_{\zeta}$ associated with the vector
field $\zeta^{a}$ is defined by

$\displaystyle\begin{aligned}
\bm{J}_{\zeta}=\bm{\Theta}(g,\mathcal{L}_{\zeta}g)-\zeta\cdot\bm{L}_{\text{grav}}\,.\end{aligned}$
(11)

For another, the Noether current can also be expressed asIyer:1994ys

$\displaystyle\begin{aligned}
\bm{J}_{\zeta}=\bm{C}_{\zeta}+d\bm{Q}_{\zeta}\end{aligned}$ (12)

with the Noether charge two-form

$\displaystyle\begin{aligned}
\left(\bm{Q}_{\zeta}\right)_{ab}=-\frac{1}{16\pi}\bm{\epsilon}_{abcd}\nabla^{c}\zeta^{d}\end{aligned}$
(13)

and the constraint

$\displaystyle\begin{aligned}
(\bm{C}_{\zeta})_{abc}=\frac{1}{8\pi}\bm{\epsilon}_{eabc}\zeta^{d}G_{d}{}^{e}\,.\end{aligned}$
(14)

The derivative of the identity

$\displaystyle\begin{aligned}
\bm{\Theta}(g,\mathcal{L}_{\zeta}g)-\zeta\cdot\bm{L}_{\text{grav}}=\bm{C}_{\zeta}+d\bm{Q}_{\zeta}\end{aligned}$
(15)

gives the variational identity

$\displaystyle\begin{aligned} \quad
d&\left[\frac{d\bm{Q}_{\zeta}}{d\lambda}-\zeta\cdot\bm{\Theta}\left(g,\frac{dg}{d\lambda}\right)\right]=\bm{\omega}\left(g;\frac{dg}{d\lambda},\mathcal{L}_{\zeta}g\right)-\zeta\cdot\bm{E}_{g}^{ab}\frac{dg_{ab}}{d\lambda}-\frac{d\bm{C}_{\zeta}}{d\lambda}\,.\end{aligned}$
(16)

It is necessary to note that the above identity holds for any spacetime
configuration even it does not satisfy the on-shell equation of motion
$G_{ab}=0$.

##  3 Perturbed geometry of Kerr black hole

In this section, we would like to introduce the perturbation geometry of the
Kerr black hole and introduce the process that some perturbation matters fall
into a nearly extremal Kerr black hole. In this case, the full Lagrangian is
given by

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

where we use $\bm{L}_{\text{mt}}$ to denote the Lagrangian of the perturbation
matter fields. From the variation of the Lagrangian, we can find the equation
of motion is

$\displaystyle G_{ab}=T_{ab}\,,$ (18)

where $T_{ab}$ is the stress-energy tensor of the perturbation matter field.
Now, we consider a perturbation to the Kerr black hole, i.e., we have
$T_{ab}=0$ in the background geometry. The metric of the Kerr solution in the
Boyer-Lindquist coordinate system can be written as

$\displaystyle\begin{aligned}
ds^{2}=&-\left[1-\frac{2Mr}{\rho^{2}}\right]dt^{2}+\frac{\rho^{2}}{\Delta}dr^{2}+\rho^{2}d\theta^{2}\\\
&+\left[(r^{2}+a^{2})^{2}-\Delta
a^{2}\sin^{2}\theta\right]\frac{\sin^{2}\theta}{\rho^{2}}d\varphi^{2}-\frac{4aMr\sin^{2}\theta}{\rho^{2}}dtd\varphi\end{aligned}$
(19)

with

$\displaystyle\begin{aligned}
&a=J/M\,,\quad\quad\rho^{2}=r^{2}+a^{2}\cos^{2}\theta\,,\quad\quad\Delta(r)=r^{2}-2Mr+a^{2}\,,\end{aligned}$
(20)

in which $M$ and $J$ are the mass and angular momentum of the spacetime and
they are defined by

$\displaystyle\begin{aligned} \delta
M=\int_{\infty}\delta\bm{Q}_{t}-t\cdot\bm{\Theta}(g,\delta
g)\,,\quad\quad\delta J=\int_{\infty}\delta\bm{Q}_{\varphi}\,.\end{aligned}$
(21)

Here

$\displaystyle t^{a}=\left(\frac{\partial}{\partial
t}\right)^{a}\,,\quad\varphi^{a}=\left(\frac{\partial}{\partial\varphi}\right)^{a}$
(22)

are the killing vectors related to the time evolution and axial symmetry at
asymptotically infinity, individually.

The horizon of the black hole is given by the positive root of
$\Delta(r_{h})=0$, i.e.,

$\displaystyle\begin{aligned}
r_{h}=M+\sqrt{M^{2}-J^{2}/M^{2}}\,.\end{aligned}$ (23)

That is to say, when $M^{2}\geq J^{2}/M^{2}$, the solution describes a black
hole geometry; when $M^{2}<J^{2}/M^{2}$, the solution describes a naked
singularity. Especially, when $M^{4}=J^{2}$, it describes an extremal black
hole.

For the black hole case, the area, angular velocity, and surface gravity of
the event horizon is given by

$\displaystyle\begin{aligned}
A=4\pi\left(r_{h}^{2}+J^{2}/M^{2}\right)\,,\quad\quad\Omega_{H}=\frac{J/M}{r_{h}^{2}+J^{2}/M^{2}}\,,\quad\quad\kappa=\frac{M^{2}-J^{2}/M^{2}}{2Mr_{h}}\,.\end{aligned}$
(24)

To test the WCCC, we consider a perturbation process that some collision
matters fall into the future horizon during a finite timeSorce:2017dst , i.e.,
the perturbation matter fields vanish at sufficiently late times as well as
some early time before perturbation (See Fig.1). Let $\phi(\lambda)$ to be a
one-parameter family in the configuration space, and each element
$\phi(\lambda)$ in this family represents a perturbation process and
$\phi=\phi(0)$ describes the background geometry, which can be described by
the Kerr metric with the mass $M$ and angular momentum $J$. Here we introduce
$\phi$ to denote the collection of the metric $g_{ab}$ and the perturbation
matter fields. We assume that the perturbation vanishes at early times, i.e.,
the spacetime before the perturbation is the same as the background, and we
also assume that the spacetime satisfies the stability conditionSorce:2017dst
; Hollands:2012sf , which means that the late-time geometry can also be
described by the Kerr metric with different mass $M(\lambda)$ and angular
momentum $J(\lambda)$ and therefore we have
$\bm{E}_{g}^{ab}(\lambda)=\bm{C}_{\zeta}(\lambda)=0$ at late times.

Figure 1: Plot showing the perturbed geometry of the Kerr black hole and the
choice of the hypersurface $\Sigma=\mathcal{H}\cup\Sigma_{1}$.

Under the stability condition, examining the WCCC is equivalent to checking
whether the late-time geometry also describes a black hole. Therefore, we can
define a function

$\displaystyle\begin{aligned}
h(\lambda)=M(\lambda)^{4}-J(\lambda)^{2}\end{aligned}$ (25)

If $h(\lambda)\geq 0$, there exists an event horizon after perturbation and
the WCCC is satisfied in this process. If $h(\lambda)<0$, the black hole is
destroyed and the WCCC is violated by the perturbation process.

In this paper, we would like to consider some perturbation matter falling into
a nearly extremal Kerr black hole. Without loss of generality, we set $M=1$
for the mass of the background black hole. For the nearly extremal black
holes, we can define a small quantity $\epsilon$ given by

$\displaystyle\begin{aligned} J=\sqrt{1-\epsilon^{2}}\end{aligned}$ (26)

to describe the deviation from the extremal black hole. To make the black hole
easier to be destroyed, this parameter is chosen to be in the same order with
$\lambda$ Sorce:2017dst . Expanding $h(\lambda)$ at $\lambda=0$, we have

$\displaystyle\begin{aligned} h(\lambda)=&1-J^{2}+\lambda\left(4\delta
M-2J\delta J\right)+\lambda^{2}\left(6\delta M^{2}+2\delta^{2}M-\delta
J^{2}-J\delta^{2}J\right)\\\ &+\frac{\lambda^{3}}{3}(12\delta M^{3}+18\delta
M\delta^{2}M+2\delta^{3}M-3\delta J\delta^{2}J-J\delta^{3}J)\\\
&+\frac{\lambda^{4}}{12}\left(12\delta M^{4}+72\delta
M\delta^{2}M+18\delta^{2}M^{2}+24\delta M\delta^{3}M\right.\\\
&\left.+2\delta^{4}M-3\delta^{2}J^{2}-4\delta
J\delta^{3}J-J\delta^{4}J\right)+\mathcal{O}(\lambda^{4})\,.\end{aligned}$
(27)

Therefore, the key point to check the signature of $h(\lambda)$ under the
perturbation process is to get the constraint of the variational quantities
$\delta^{k}M$ and $\delta^{k}J$ from the null energy condition of the
perturbation matter fields.

##  4 Gedanken experiments to overspin the black holes

Next, we would like to derive the perturbation inequality by assuming the
perturbation matter field satisfies the null energy condition and examining
whether the horizon exists after the perturbation. The null energy condition
of the perturbation matter fields demand

$\displaystyle\begin{aligned} T_{ab}(\lambda)l^{a}l^{b}\geq 0\,,\end{aligned}$
(28)

in which $l^{a}$ is any future-directed null vector in the configuration
$\phi(\lambda)$.

With a same setup as Ref. Sorce:2017dst , we choose a hypersurface
$\Sigma=\mathcal{H}\cup\Sigma_{1}$ as shown in Fig.1 satisfying the following
conditions: (1) $\mathcal{H}$ is a null hypersurface which starts from the
bifurcation surface $B$ before the perturbation and finally ends at a very
late-time cross-section $B_{1}$ where all the perturbation matter fields
vanish. (2) $\Sigma_{1}$ is a spacelike hypersurface that starts from $B_{1}$
and extends to infinity. (4) $\mathcal{H}$ is the Killing horizon after
perturbation and therefore it is the event horizon in the background geometry.
Let $\xi^{a}$ be the null generator of the hypersurface $\mathcal{H}$. In the
background geometry, it becomes the Killing vector of the event horizon, i.e.,

$\displaystyle\begin{aligned}
\xi^{a}=t^{a}+\Omega_{H}\varphi^{a}\end{aligned}$ (29)

On the background, $\xi^{a}$ is the Killing vector field of the background
horizon and the tangent vector of $\mathcal{H}$. Choose a gauge to make the
hypersurface $\Sigma$ and the vector field $\xi^{a}$ fixed under variation.
Replacing $\zeta^{a}$ by $\xi^{a}$, integration of the variational identity
(16) on $\Sigma$ gives

$\displaystyle\begin{aligned}
\int_{\infty}\left[\frac{d\bm{Q}_{\xi}}{d\lambda}-\xi\cdot\bm{\Theta}\left(g,\frac{dg}{d\lambda}\right)\right]&=\int_{\Sigma}\bm{\omega}\left(g;\frac{dg}{d\lambda},\mathcal{L}_{\xi}g\right)-\int_{\mathcal{H}}\frac{d\bm{C}_{\xi}}{d\lambda}\,,\end{aligned}$
(30)

where we have used that the assumptions that the perturbation vanishes on $B$,
i.e., $dg(\lambda)/d\lambda=0$ on $B$, and $\xi^{a}$ is a tangent vector on
$\mathcal{H}$ and the perturbation matter fields vanish at late times, i.e.,
we have $\bm{E}_{g}^{ab}(\lambda)=\bm{C}_{\xi}(\lambda)=0$ on $\Sigma_{1}$. In
the following, we would like to perform this identity to get the constraint
from the null energy condition of the perturbation matter fields.

###  4.1 First-order approximation of perturbation

With a calculation same with Sorce and Wald, after evaluating the identity
(30) at $\lambda=0$ and using the stationary condition of the background
spacetime $\mathcal{L}_{\xi}g_{ab}=0$, we further obtain

$\displaystyle\begin{aligned} \delta M-\Omega_{H}\delta
J=-\int_{\mathcal{H}}\delta\bm{C}_{\xi}=\int_{\mathcal{H}}\tilde{\bm{\epsilon}}\delta
T_{ab}k^{a}\xi^{b}\,,\end{aligned}$ (31)

where $\tilde{\bm{\epsilon}}$ is the induced volume element of $\mathcal{H}$
which is defined by
$\bm{\epsilon}_{ebcd}=-4k_{[e}\tilde{\bm{\epsilon}}_{bcd]}$ and
$k^{a}\propto\xi^{a}$ is the future-directed normal vector of $\mathcal{H}$.
Considering that $T_{ab}=0$ in the background geometry, the null energy
condition under the first-order approximation implies that $\delta
T_{ab}\xi^{a}\xi^{b}\geq 0$. Then, we can get the first-order perturbation
inequality

$\displaystyle\begin{aligned} \delta M-\Omega_{H}\delta J\geq 0\end{aligned}$
(32)

with the angular momentum

$\displaystyle\begin{aligned}
\Omega_{H}=\frac{1}{2}\sqrt{\frac{1-\epsilon}{1+\epsilon}}\,.\end{aligned}$
(33)

Together with Eq. (26) and Eq. (27), we have

$\displaystyle\begin{aligned} h(\lambda)=\lambda(4\delta M-2\delta J)\geq
0\end{aligned}$ (34)

under the first-order approximation and therefore the black hole cannot be
overspun under the first-order approximation of perturbation. However, this
does not mean that the WCCC is valid in the perturbation process since there
exists a first-order optimal condition $\delta M=\Omega_{H}\delta J$ such that
$h(\lambda)=0$ under the first-order approximation and thus the signature of
$h(\lambda)$ is determined by the higher-order approximation. That is to say,
to strictly test the WCCC in the perturbation process, we need to consider the
second- or even higher-order approximations.

Finally, we would like to discuss the first-order optimal condition. From Eq.
(31), we can see that this condition implies

$\displaystyle\begin{aligned} \delta T_{ab}\xi^{a}\xi^{b}=0\,.\end{aligned}$
(35)

Let $u$ be the parameter of $\xi^{a}$, the expansion and shear of the null
hypersurface $\mathcal{H}$ are defined by

$\displaystyle\begin{aligned}
\vartheta=\frac{1}{2}\gamma^{ab}\mathcal{L}_{\xi}\gamma_{ab}\,,\quad\quad\sigma_{ab}=\frac{1}{2}\mathcal{L}_{\xi}\gamma_{ab}-\frac{\vartheta}{2}\gamma_{ab}\,,\end{aligned}$
(36)

where $\gamma_{ab}$ is the induced metric of the cross-section on
$\mathcal{H}$. Considering the Raychaudhuri equation is given by

$\displaystyle\begin{aligned}
\frac{d\vartheta}{du}=-\frac{1}{2}\vartheta^{2}-\sigma_{ab}\sigma^{ab}-\xi^{a}\xi^{b}T_{ab}+\kappa\vartheta\,,\end{aligned}$
(37)

the optimal condition gives

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

where we used the fact that $\vartheta=\sigma_{ab}=0$ in the background
geometry. Combing with the assumption that the perturbation vanishes (i.e.,
$\delta\vartheta=0$) in the region near the bifurcation surface, the above
result also implies that

$\displaystyle\delta\vartheta=0$ (39)

on $\mathcal{H}$.

###  4.2 Second-order approximation of peturbation

Next, we will derive the second-order perturbation inequality under the first-
order optimal condition. Taking a derivative of Eq. (30), a same calculation
as Ref. Sorce:2017dst  gives

$\displaystyle\begin{aligned}
\delta^{2}M-\Omega_{H}\delta^{2}J=\mathcal{E}_{\Sigma}(g,\delta
g)-\int_{\mathcal{H}}\delta^{2}\bm{C}_{\xi}\,,\end{aligned}$ (40)

where we define the canonical energy as

$\displaystyle\begin{aligned} \mathcal{E}_{\Sigma}(g,\delta
g)=\int_{\Sigma}\bm{\omega}(g;\delta g,\mathcal{L}_{\xi}\delta
g)\,.\end{aligned}$ (41)

Using the first-order optimal condition $\delta T_{ab}\xi^{a}\xi^{b}=0$, we
have

$\displaystyle\begin{aligned}
\int_{\mathcal{H}}\delta^{2}\bm{C}_{\xi}=-\int_{\mathcal{H}}\tilde{\bm{\epsilon}}k^{a}\xi^{b}\delta^{2}T_{ab}\leq
0\,.\end{aligned}$ (42)

Thus, we can obtain

$\displaystyle\begin{aligned}
\delta^{2}M-\Omega_{H}\delta^{2}J\geq\mathcal{E}_{\Sigma}(g,\delta
g)\end{aligned}$ (43)

The canonical energy can be divided into two parts:

$\displaystyle\begin{aligned}
\mathcal{E}_{\Sigma}=\mathcal{E}_{\mathcal{H}}+\mathcal{E}_{\Sigma_{1}}=\int_{\mathcal{H}}\bm{\omega}(g;\delta
g,\mathcal{L}_{\xi}\delta g)+\int_{\Sigma_{1}}\bm{\omega}(g;\delta
g,\mathcal{L}_{\xi}\delta g)\end{aligned}$ (44)

Next, we refer to the method by Sorce and Wald to calculate
$\mathcal{E}_{\Sigma_{1}}(\phi,\delta\phi)$. We consider another one-parameter
family $\phi_{(2)}(\alpha)$ in which every element is given by the Kerr
spacetime with the mass and angular momentum

$\displaystyle\begin{aligned} M_{(2)}(\alpha)=M+\alpha\delta M\,,\quad\quad
J_{(2)}(\alpha)=J+\alpha\delta J\,.\end{aligned}$ (45)

Here $\delta M$ and $\delta J$ are chosen to be agree with the quantities
evaluated by the one-parameter family $\phi(\lambda)$. Since there is only the
first-order variation of $\phi(\lambda)$ in
$\mathcal{E}_{\Sigma_{1}}(\phi,\delta\phi)$, we have

$\displaystyle\begin{aligned}
\mathcal{E}_{\Sigma_{1}}(\phi,\delta\phi)=\mathcal{E}_{\Sigma_{1}}(\phi_{(2)},\delta\phi_{(2)})\,.\end{aligned}$
(46)

In this family, considering that $\delta^{2}M_{(2)}=\delta^{2}J_{(2)}=0$, we
can obtain Sorce:2017dst ; Hollands:2012sf

$\displaystyle\begin{aligned} \mathcal{E}_{\Sigma_{1}}(g,\delta
g)=-\frac{\kappa}{8\pi}\delta^{2}A_{(2)}\,,\end{aligned}$ (47)

in which $A_{(2)}(\alpha)$ is the area of the bifurcation surface $B$ in the
configuration $\phi_{(2)}(\alpha)$. Through a straightforward calculation, we
can obtain

$\displaystyle\begin{aligned}
\delta^{2}A_{(2)}=-4\pi\left(1+\frac{1}{\epsilon}\right)\delta
J^{2}\,,\end{aligned}$ (48)

where we used the first-order optimal condition $\delta M-\Omega\delta J=0$ to
simplify.

Then, we would like to evaluate the horizon part of the canonical energy. To
do this, we introduce the Gaussian null coordinate system of the hypersurface
$\mathcal{H}$. Let the null tangent vector $\xi^{a}=(\partial/\partial u)^{a}$
be the first coordinate basis. Another coordinate basis $(\partial/\partial
z)^{a}$ is given by another null vector and satisfies

$\displaystyle\begin{aligned} \left.\left(\frac{\partial}{\partial
u}\right)^{a}\left(\frac{\partial}{\partial
z}\right)_{a}\right|_{\mathcal{H}}=1\,.\end{aligned}$ (49)

Then, the metric near the null hypersurface $\mathcal{H}$ 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}$
(50)

Here $z=0$ gives the location of $\mathcal{H}$. It is necessary to mention
that this coordinate system is compatible with the gauge choice in previous
calculations, i.e., fixing the coordinate in the variation also makes
$\xi^{a}$ and $\mathcal{H}$ fix.

Using the metric, we can find $\delta g_{ab}=\delta\gamma_{ab}$ on
$\mathcal{H}$. Using the definition of the expansion, the first-order optimal
condition $\delta\vartheta=0$ implies

$\displaystyle g^{ab}\mathcal{L}_{\xi}\delta
g_{ab}=\gamma^{ab}\mathcal{L}_{\xi}\delta\gamma_{ab}=0\,,\quad\text{and}\quad\delta\sigma_{ab}=\frac{1}{2}\mathcal{L}_{\xi}\delta\gamma_{ab}$
(51)

on $\mathcal{H}$. Considering the assumption that the perturbation vanishes at
early times, the first condition of the above expressions also implies

$\displaystyle\begin{aligned} g^{ab}\delta
g_{ab}=\gamma^{ab}\delta\gamma_{ab}=0\end{aligned}$ (52)

on $\mathcal{H}$.

According to Eq. (9), $\mathcal{E}_{\mathcal{H}}(g,\delta g)$ can be written
as

$\displaystyle\begin{aligned} \mathcal{E}_{\mathcal{H}}(g,\delta
g)&=-\frac{1}{16\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\xi_{a}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)\\\
&=-\frac{1}{8\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\bar{P}^{zbcdef}\mathcal{L}_{\xi}\delta
g_{bc}\nabla_{d}\delta
g_{ef}+\frac{1}{16\pi}\int_{\mathcal{H}}\mathcal{L}_{\xi}(\bar{\bm{\epsilon}}\bar{P}^{zbcdef}\delta
g_{bc}\nabla_{d}\delta g_{ef})\,,\end{aligned}$ (53)

where we have defined

$\displaystyle\begin{aligned}
\bar{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}\,,\end{aligned}$
(54)

and used the first-order optimal condition $g^{ab}\mathcal{L}_{\xi}\delta
g_{ab}=0$ on $\mathcal{H}$. Here we introduce the notations

$\displaystyle\begin{aligned} \bar{\bm{\epsilon}}=\sqrt{\gamma}du\wedge
d\theta^{1}\wedge
d\theta^{2}\,,\quad\quad\hat{\bm{\epsilon}}=\sqrt{\gamma}d\theta^{1}\wedge
d\theta^{2}\,.\end{aligned}$ (55)

With a straightforward calculation using the Gaussian null coordinate system,
it is easy to find

$\displaystyle\begin{aligned} \bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\delta
g_{ef}=-\frac{1}{2}\delta\gamma_{ij;u}=-\delta\sigma_{ij}+\delta\gamma_{k(i}\Gamma^{k}{}_{j)u}\,.\end{aligned}$
(56)

On background, we have $\Gamma^{i}{}_{ju}=0$. Besides, according to the
assumption by Sorce and WaldSorce:2017dst , the perturbation is physically
stationary on $B_{1}$, which means $\delta\sigma_{ab}=0$ on $B_{1}$. Together
with the first-order optimal condition (51), we have

$\displaystyle\begin{aligned}
\mathcal{E}_{H}(\phi,\delta\phi)&=\frac{1}{4\pi}\int_{H}\bar{\bm{\epsilon}}\delta\sigma_{ab}\delta\sigma^{ab}+\frac{1}{16\pi}\int_{B_{1}}\hat{\bm{\epsilon}}\delta
g^{ab}\delta\sigma_{ab}\\\
&=\frac{1}{4\pi}\int_{H}\bar{\bm{\epsilon}}\delta\sigma_{ab}\delta\sigma^{ab}\geq
0\,.\end{aligned}$ (57)

Summarizing the above results, the second-order perturbation inequality gives

$\displaystyle\begin{aligned}
\delta^{2}M-\Omega_{H}\delta^{2}J\geq-\frac{\kappa}{8\pi}\delta^{2}A_{(2)}=\frac{\delta
J^{2}}{4}\,.\end{aligned}$ (58)

where we used

$\displaystyle\kappa=\frac{\epsilon}{2+2\epsilon}$ (59)

and the first-order optimal condition.

After considering the first-order optimal condition and the second-order
perturbation inequality, we can get

$\displaystyle\begin{aligned} h(\lambda)\geq(\epsilon-\lambda\delta J)^{2}\geq
0\end{aligned}$ (60)

under the second-order approximation of perturbation, i.e., we neglected the
higher-order terms
$\mathcal{O}(\epsilon^{3},\lambda^{3},\epsilon\lambda^{2},\dots)$. This result
implies that the nearly extremal Kerr black hole cannot be overspun in the
second-order approximation of perturbation. Same to the first-order case,
there also exists a second-order optimal condition in which the first two
perturbation inequalities are saturated and

$\displaystyle\delta J=\epsilon/\lambda\,,$ (61)

in which we cannot determine the validity of the WCCC in the second-order
approximation. Therefore, we need to extend the discussion into the third- or
even higher-order approximation.

With the same discussion as the first-order case, it is not hard to verify
that the second-order optimal condition will give the additional conditions

$\displaystyle\begin{aligned}
&\delta^{2}T_{ab}\xi^{a}\xi^{b}=0\,,\quad\quad\delta\sigma_{ab}=0\,.\end{aligned}$
(62)

Together with $\delta\vartheta=0$, these also gives

$\displaystyle\begin{aligned}
\mathcal{L}_{\xi}\delta\gamma_{ab}=\delta\gamma_{ab}=0\end{aligned}$ (63)

on $\mathcal{H}$. Using the above conditions and combing with the first-order
optimal condition, the second variation of the Raychaudhuri equation gives

$\displaystyle\frac{d\delta^{2}\vartheta}{du}=\kappa\delta^{2}\vartheta\,.$
(64)

Considering that $\delta^{2}\vartheta=0$ on $\mathcal{H}$ at early times, this
implies that $\delta^{2}\vartheta=0$ on the whole $\mathcal{H}$.

###  4.3 Third-order approximation of perturbation

In this subsection, we extend the discussion into the third-order
approximation. Taking two variation to Eq. (30), similar calculation can give

$\displaystyle\begin{aligned}
\delta^{3}M-\Omega_{H}\delta^{3}J&=\delta\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta\mathcal{E}_{\Sigma_{1}}(g,\delta
g)+\int_{H}\tilde{\bm{\epsilon}}k^{a}\xi^{b}\delta^{3}T_{ab}\\\
&\geq\delta\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta\mathcal{E}_{\Sigma_{1}}(g,\delta g)\,,\end{aligned}$ (65)

where we used the null energy condition $\delta^{3}T_{ab}\xi^{a}\xi^{b}\geq 0$
of the perturbation matter field under the third-order approximation in the
second-order optimal condition.

First, we evaluate the canonical energy contributed by $\Sigma_{1}$.
Considering that $\delta\mathcal{E}_{\Sigma_{1}}(g,\delta g)$ only depends on
the first two order variation of $g_{ab}$. We can also calculate it by
introducing a one-parameter family $\phi_{(3)}(\alpha)$ in which any element
describes a Kerr solution with the mass and angular momentum

$\displaystyle\begin{aligned} &M_{(3)}(\alpha)=M+\alpha\delta
M+\frac{1}{2}\alpha^{2}\delta^{2}M\,,\\\ &J_{(3)}(\alpha)=J+\alpha\delta
J+\frac{1}{2}\alpha^{2}\delta^{2}J\,,\end{aligned}$ (66)

where $\delta M,\delta^{2}M,\delta J$ and $\delta^{2}J$ are chosen to agree
with the variations from the configuration $\phi(\lambda)$ under the first-two
order approximation. With a similar calculation, we can also get

$\displaystyle\begin{aligned} \delta\mathcal{E}_{\Sigma_{1}}(g,\delta
g)&=-\frac{\kappa}{8\pi}\delta^{3}A_{(3)}=\frac{3\delta
J}{8}\left(2\delta^{2}J-\sqrt{\frac{1-\epsilon}{1+\epsilon}}\delta
J^{2}\right)\,.\end{aligned}$ (67)

Finally, we turn to evaluate the canonical energy from the hypersurface
$\mathcal{H}$. From the second-order optimal condition, we also have

$\displaystyle\begin{aligned}
&g^{ab}\mathcal{L}_{\xi}\delta^{2}g_{ab}=g^{ab}\delta^{2}g_{ab}=0\,,\quad\quad\delta^{2}\sigma_{ab}=\frac{1}{2}\mathcal{L}_{\xi}\delta^{2}\gamma_{ab}\,,\quad\quad\delta
g_{ab}=\delta\gamma_{ab}=\delta\sigma_{ab}=0\end{aligned}$ (68)

on $\mathcal{H}$. In addition, since $\xi^{a}$ is killing vector on the
background, it is easy to check

$\displaystyle\begin{aligned} \delta(\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef})=\nabla_{d}(\mathcal{L}_{\xi}\delta^{2}g_{ef})=\mathcal{L}_{\xi}(\nabla_{d}\delta^{2}g_{ef})\end{aligned}$
(69)

According to Eq. (9), $\delta\mathcal{E}_{\mathcal{H}}(g,\delta g)$ can be
written as

$\displaystyle\begin{aligned} \delta\mathcal{E}_{\mathcal{H}}(g,\delta
g)&=-\frac{1}{8\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\bar{P}^{zbcdef}\mathcal{L}_{\xi}\delta^{2}g_{bc}\nabla_{d}\delta
g_{ef}+\frac{1}{16\pi}\int_{\mathcal{H}}\mathcal{L}_{\xi}(\bar{\bm{\epsilon}}\bar{P}^{zbcdef}\delta
g_{bc}\nabla_{d}\delta^{2}g_{ef})\\\
&=\frac{1}{4\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\delta^{2}\sigma_{ab}\delta\sigma^{ab}+\frac{1}{16\pi}\int_{B_{1}}\hat{\bm{\epsilon}}\delta
g^{ab}\delta^{2}\sigma_{ab}=0\end{aligned}$ (70)

where we used the second-order optimal condition Eq.
(LABEL:optimalsecondorder).

Summarizing the above results, the third-order perturbation inequality can be
shown as

$\displaystyle\begin{aligned}
\delta^{3}M-\Omega_{H}\delta^{3}J\geq-\frac{\kappa}{8\pi}\delta^{3}A_{(3)}\,.\end{aligned}$
(71)

After taking these conditions into account and together with the second-order
optimal condition, it is not hard to show

$\displaystyle\begin{aligned} h(\lambda)\geq 0\end{aligned}$ (72)

under the third-order approximation of perturbation, i.e., we neglect the
higher-order terms
$\mathcal{O}(\epsilon^{4},\lambda^{4},\lambda\epsilon^{3}\dots)$. This result
shows that the WCCC cannot be violated under the fourth-order approximation.
Similarly, there also exists a third-order optimal condition leading us to
consider the higher-order approximations, i.e., the second-order optimal
condition is satisfied and the third-order perturbation inequality is
saturated.

Under the third-order optimal condition, a similar calculation gives the
following additional conditions:

$\displaystyle\begin{aligned}
\delta^{3}T_{ab}\xi^{a}\xi^{b}=0\,,\quad\quad\delta^{3}\vartheta=0\end{aligned}$
(73)

on $\mathcal{H}$

###  4.4 Fourth-order approximation of perturbation

Taking two variations to Eq. (30) and considering the null energy condition
$\delta^{4}T_{ab}\xi^{a}\xi^{b}\geq 0$ of the perturbation matter fields under
the third-order optimal condition, the similar calculation can give

$\displaystyle\begin{aligned}
\delta^{4}M-\Omega_{H}\delta^{4}J&=\delta^{2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta^{2}\mathcal{E}_{\Sigma_{1}}(g,\delta
g)+\int_{\mathcal{H}}\tilde{\bm{\epsilon}}k^{a}\xi^{b}\delta^{4}T_{ab}\\\
&\geq\delta^{2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta^{2}\mathcal{E}_{\Sigma_{1}}(g,\delta g)\,.\end{aligned}$ (74)

For the canonical energy contributed by $\Sigma_{1}$. Similarly, it can be
calculated in the one-parameter family $\phi_{(4)}(\alpha)$ constructed by the
Kerr spacetime with

$\displaystyle\begin{aligned} &M_{(4)}(\alpha)=M+\alpha\delta
M+\frac{1}{2}\alpha^{2}\delta^{2}M+\frac{1}{3!}\alpha^{3}\delta^{3}M\,,\\\
&J_{(4)}(\alpha)=J+\alpha\delta
J+\frac{1}{2}\alpha^{2}\delta^{2}J+\frac{1}{3!}\alpha^{3}\delta^{3}J\,.\end{aligned}$
(75)

Then, we can get

$\displaystyle\begin{aligned} \delta^{2}\mathcal{E}_{\Sigma_{1}}(g,\delta
g)&=-\frac{\kappa}{8\pi}\delta^{4}A_{(4)}\\\ &=\delta
J\delta^{3}J+\frac{3}{4}\delta^{2}J^{2}-\frac{9\sqrt{1-\epsilon}}{4\sqrt{1+\epsilon}}\delta
J^{2}\delta^{2}J+\frac{3(3-5\epsilon)}{16(1+\epsilon)}\delta
J^{4}\,.\end{aligned}$ (76)

The third-order optimal condition gives

$\displaystyle\begin{aligned}
&g^{ab}\mathcal{L}_{\xi}\delta^{3}g_{ab}=g^{ab}\delta^{3}g_{ab}=0\end{aligned}$
(77)

on $\mathcal{H}$, where we have used the second-order optimal condition
$\delta g_{ab}=0$ on $\mathcal{H}$. Then,
$\delta^{2}\mathcal{E}_{\mathcal{H}}(g,\delta g)$ can be shown as

$\displaystyle\begin{aligned} \delta^{2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)=-\frac{1}{8\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\left[\mathcal{L}_{\xi}\delta^{2}\gamma_{bc}\delta(\bar{P}^{zbcdef}\nabla_{d}\delta
g_{ef})-\delta^{2}\gamma_{bc}\delta(\bar{P}^{zbcdef}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef})\right]\,.\quad\end{aligned}$ (78)

We calculate the variation term in the first term and find

$\displaystyle\begin{aligned}
\delta(\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\delta
g_{ef})&=\delta(-\delta\sigma_{ij}+\delta\gamma_{k(i}\Gamma^{k}{}_{j)u})\\\
&=-\delta^{2}\sigma_{ij}+\sum_{m=0}^{1}\delta^{m+1}\gamma_{k(i}\delta^{1-m}\Gamma^{k}{}_{j)u}\\\
&=-\delta^{2}\sigma_{ij}\,,\end{aligned}$ (79)

as well as

$\displaystyle\begin{aligned}
\delta(\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef})&=\bar{P}^{z}{}_{ij}{}^{def}\mathcal{L}_{\xi}\nabla_{d}\delta^{2}g_{ef}=\mathcal{L}_{\xi}(\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\delta^{2}g_{ef})\\\
&=\mathcal{L}_{\xi}(-\delta^{2}\sigma_{ij}+\delta^{2}\gamma_{k(i}\Gamma^{k}{}_{j)u})=-\mathcal{L}_{\xi}\delta^{2}\sigma_{ij}\,,\end{aligned}$
(80)

where we have used the third-optimal condition, $\delta
g_{ab}=\delta\gamma_{ab}=0$ and $\Gamma^{i}{}_{iu}=0$ on $\mathcal{H}$. Using
the above results and considering the stability condition
$\delta^{2}\sigma_{ab}=0$ on $B_{1}$, we can further obtain

$\displaystyle\begin{aligned} \delta^{2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)&=\frac{1}{4\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\mathcal{L}_{\xi}\delta^{2}\gamma_{bc}\delta^{2}\sigma^{bc}+\frac{1}{8\pi}\int_{B_{1}}\hat{\bm{\epsilon}}\delta^{2}\gamma_{bc}\delta^{2}\sigma^{bc}\\\
&=\frac{1}{2\pi}\int_{H}\bar{\bm{\epsilon}}\delta^{2}\sigma_{bc}\delta^{2}\sigma^{bc}\geq
0\,,\end{aligned}$ (81)

Combing these results, the fourth-order perturbation inequality can be given
by

$\displaystyle\begin{aligned}
\delta^{4}M-\Omega_{H}\delta^{4}J\geq-\frac{\kappa}{8\pi}\delta^{4}A_{(4)}\,.\end{aligned}$
(82)

After a straightforward calculation, we find

$\displaystyle\begin{aligned}
h(\lambda)\geq\frac{1}{4}(\epsilon^{2}-\lambda^{2}\delta^{2}J)^{2}\geq
0\end{aligned}$ (83)

under the fourth-order approximation of $\lambda$ and $\epsilon$. Therefore
the WCCC cannot be violated under the fourth-order approximation. Similarly,
there exits a fourth-order optimal condition in which
$\delta^{2}J=\epsilon^{2}/\lambda^{2}$, fourth-order perturbation inequality
is saturated, and the third-order optimal condition is satisfied, such that
the fourth-order approximation cannot judge the signature of $h(\lambda)$.

Similar discussions can show that the fourth-order optimal condition can also
give

$\displaystyle\begin{aligned}
&\delta^{4}T_{ab}\xi^{a}\xi^{b}=0\,,\quad\quad\delta^{2}\sigma_{ab}=0\,.\end{aligned}$
(84)

These imply

$\displaystyle\begin{aligned}
\delta^{4}\vartheta=0\,,\quad\quad\delta^{2}\gamma_{ab}=\delta^{2}g_{ab}=0\,.\end{aligned}$
(85)

###  4.5 $k$th-order approximation of perturbation

####  4.5.1 $k$th-order perturbation inequality

In this subsection, we would like to discuss the gedanken experiment under the
higher-order perturbation inequality. By summarizing the previous results, we
would like to derive the following $n$th-order perturbation inequality

$\displaystyle\begin{aligned}
\delta^{n}M-\Omega_{H}\delta^{n}J\geq-\frac{\kappa}{8\pi}\delta^{n}A_{(n)}\end{aligned}$
(86)

when the first $(n-1)$th order perturbation inequalities are saturated, and
the saturation of the first $n$th order perturbation inequalities give

$\displaystyle\begin{aligned}
&\delta^{i}T_{ab}\xi^{a}\xi^{b}=0\,,\quad\forall\,i\leq n\,,\\\
&\delta^{i}\vartheta=\frac{1}{2}\gamma^{ab}\mathcal{L}_{\xi}\delta^{i}\gamma_{ab}=0\,,\quad\forall\,i\leq
n\,,\\\ &\delta^{j}g_{ab}=\delta^{j}\gamma_{ab}=0\,,\quad\forall\,1\leq
j\leq[n/2]\,,\end{aligned}$ (87)

on the hypersurface $\mathcal{H}$, where $[n/2]$ denote the integer part of
$n/2$. Here $A_{(n)}(\alpha)$ is the area of bifurcation surface $B$ in the
one-parameter family $\phi_{(n)}(\alpha)$ constructed by the Kerr solution
with the parameters

$\displaystyle\begin{aligned}
M_{(n)}(\alpha)=\sum_{i=0}^{n-1}\frac{\alpha^{i}\delta^{i}M}{i!}\,,\quad\quad
J_{(n)}(\alpha)=\sum_{i=0}^{n-1}\frac{\alpha^{i}\delta^{i}J}{i!}\,.\end{aligned}$
(88)

Next, we would like to perform the mathematical induction to prove the above
result, i.e., we first assume that they are satisfied at $n$th order and then
prove that they hold for $(n+1)$th order. However, we can note that the
condition of Eq. (LABEL:optimaln) is different for the odd and even orders.
Therefore, we need to prove the above results for $n=2k-1$ and $n=2k$
separately.

(1) Case of $n=2k$

We assume that the Eqs. (86) and (LABEL:optimaln) are satisfied when $n=2k-1$
for $k\geq 1$. Taking $(2k-1)$ variation to Eq. (30) and considering the null
energy condition $\delta^{2k}T_{ab}\xi^{a}\xi^{b}\geq 0$ of the perturbation
matter fields when the first $(2k-1)$ order perturbation inequalities are
saturated, we can further obtain

$\displaystyle\begin{aligned}
\delta^{2k}M-\Omega_{H}\delta^{2k}J&=\delta^{2k-2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta^{2k-2}\mathcal{E}_{\Sigma_{1}}(g,\delta
g)+\int_{H}\tilde{\bm{\epsilon}}k^{a}\xi^{b}\delta^{2k}T_{ab}\\\
&\geq\delta^{2k-2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta^{2k-2}\mathcal{E}_{\Sigma_{1}}(g,\delta g)\,.\end{aligned}$ (89)

For the second term in the right side of the above expression, since
$\delta^{2k}\mathcal{E}_{\Sigma_{1}}(g,\delta g)$ only contains the first $2k$
order variational of $g_{ab}$, it can be evaluated by the one-parameter family
$\phi_{(2k)}(\alpha)$ and finally we have

$\displaystyle\begin{aligned} \delta^{2k-2}\mathcal{E}_{\Sigma_{1}}(g,\delta
g)&=-\frac{\kappa}{8\pi}\delta^{2k}A_{(2k)}\,.\end{aligned}$ (90)

From Eq. (9), using the Gaussian null coordinate system (50),
$\delta^{2k-2}\mathcal{E}_{\mathcal{H}}(g,\delta g)$ can be written as

$\displaystyle\begin{aligned} \delta^{2k-2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)&=-\frac{1}{16\pi}\delta^{2k-2}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\xi_{a}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)\\\
&=-\frac{1}{16\pi}\delta^{2k-2}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\left(\mathcal{L}_{\xi}\delta\gamma_{bc}P^{zbcdef}\nabla_{d}\delta
g_{ef}-\delta\gamma_{bc}P^{zbcdef}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef}\right)\,.\end{aligned}$ (91)

Using the conditions (LABEL:optimaln) with $n=2k-1$, especially

$\displaystyle\begin{aligned}
&\delta^{j}g_{ab}=\delta^{j}\gamma_{ab}=0\,,\quad\quad\forall 1\leq j\leq
k-1\,,\\\
&g^{ab}\mathcal{L}_{\xi}\delta^{i}\gamma_{ab}=\gamma^{ab}\mathcal{L}_{\xi}\delta^{i}\gamma_{ab}=0\,,\quad\forall\,i\leq
2k-1\end{aligned}$ (92)

on $\mathcal{H}$, we can obtain

$\displaystyle\begin{aligned}
\delta^{m}[g^{bc}\mathcal{L}_{\xi}\delta\gamma_{bc}]=\sum_{i=0}^{m}C_{m}^{i}\delta^{i}(g^{bc})\mathcal{L}_{\xi}\delta^{2k-i+1}\gamma_{bc}=g^{bc}\mathcal{L}_{\xi}\delta^{m+1}\gamma_{bc}=0\end{aligned}$
(93)

on $\mathcal{H}$ for any $m\leq 2k-2$, where $C_{m}^{i}$ is the binomial
coefficient. Considering the assumption that the perturbation vanishes at
early time, the optimal condition
$g^{ab}\mathcal{L}_{\xi}\delta^{i}\gamma_{ab}=0$ implies that

$\displaystyle\begin{aligned}
g^{ab}\delta^{i}\gamma_{ab}=0\,,\quad\quad\forall i\leq 2k-1\,.\end{aligned}$
(94)

Then, similar calculation can give $\delta^{2k-2}[g^{bc}\delta\gamma_{bc}]=0$.
Therefore, we can replace $P^{zbcdef}$ by $\bar{P}^{zbcdef}$. Again using
these conditions, from the variation of Eq. (56), for any integer $0\leq m\leq
k-2$, it is not hard to verify

$\displaystyle\begin{aligned}
\delta^{m}[\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\delta
g_{ef}]&=\delta^{m}[-\delta\sigma_{ij}+\delta\gamma_{l(i}\Gamma^{l}{}_{j)u}]\\\
&=-\delta^{m+1}\sigma_{ij}+\sum_{i=0}^{m}C_{m}^{i}\delta^{i+1}\gamma_{l(i}\delta^{m-i}[\Gamma^{l}{}_{j)u}]=0\,.\end{aligned}$
(95)

Similar features are also satisfied for
$\delta^{m}[\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef}]$. Then, (91) can be simplified to

$\displaystyle\begin{aligned}
&\delta^{2k-2}\mathcal{E}_{\mathcal{H}}=-\frac{C_{2k-2}^{k-1}}{16\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\left[\mathcal{L}_{\xi}\delta^{k}\gamma_{bc}\delta^{k-1}(\bar{P}^{zbcdef}\nabla_{d}\delta
g_{ef})-\delta^{k}\gamma_{bc}\delta^{k-1}(\bar{P}^{zbcdef}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef})\right]\,.\end{aligned}$

Same calculation as Eq. (95), we can further obtain

$\displaystyle\begin{aligned}
\delta^{k-1}[\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\delta
g_{ef}]=-\delta^{k}\sigma_{ij}\,,\quad\quad\delta^{k-1}[\bar{P}^{z}{}_{ij}{}^{def}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef}]=-\mathcal{L}_{\xi}\delta^{k}\sigma_{ij}\,.\end{aligned}$ (97)

Then, straightforward calculation gives

$\displaystyle\begin{aligned} \delta^{2k-2}\mathcal{E}_{\mathcal{H}}(g,\delta
g)&=\frac{C_{2k-2}^{k-1}}{16\pi}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\left[\mathcal{L}_{\xi}\delta^{k}\gamma_{bc}\delta^{k}\sigma^{bc}-\delta^{k}\gamma_{bc}\mathcal{L}_{\xi}(\delta^{k}\sigma^{bc})\right]\\\
&=\frac{C_{2k-2}^{k-1}}{8\pi}\int_{H}\bar{\bm{\epsilon}}\mathcal{L}_{\xi}\delta^{k}\gamma_{bc}\delta^{k}\sigma^{bc}+\frac{C_{2k-2}^{k-1}}{16\pi}\int_{B_{1}}\tilde{\bm{\epsilon}}\delta^{k}\gamma_{bc}\delta^{k}\sigma^{bc}\\\
&=\frac{C_{2k-2}^{k-1}}{8\pi}\int_{H}\bar{\bm{\epsilon}}\mathcal{L}_{\xi}\delta^{k}\gamma_{bc}\delta^{k}\sigma^{bc}\geq
0\,,\end{aligned}$ (98)

where we used the stability condition $\delta^{k}\sigma_{ab}=0$ on $B_{1}$.

So far, we have proved $\delta^{2k-2}\mathcal{E}_{\mathcal{H}}\geq 0$ under
the condition (LABEL:optimaln) with $n=2k-1$. The $(2k)$th-order perturbation
inequality becomes

$\displaystyle\begin{aligned}
\delta^{2k}M-\Omega_{H}\delta^{2k}J\geq-\frac{\kappa}{8\pi}\delta^{2k}A_{(2k)}\,.\end{aligned}$
(99)

Next, we need to show that the saturation of the $(2k)$th order perturbation
inequality can give the results in Eq. (LABEL:optimaln) with $n=2k$. From the
above calculation, it is easy to see that saturation of the $(2k)$th
inequality needs the following additional conditions

$\displaystyle\begin{aligned}
\delta^{2k}T_{ab}\xi^{a}\xi^{b}=0\,,\quad\quad\delta^{k}\sigma_{ab}=0\end{aligned}$
(100)

on $\mathcal{H}$. From the Raychaudhuri equation (37), we have

$\displaystyle\begin{aligned}
\frac{d\delta^{2k}\vartheta}{du}-\kappa\delta^{2k}\vartheta=-\frac{1}{2}\sum_{i=0}^{2k}C_{2k}^{i}\delta^{i}\vartheta\delta^{2k-i}\vartheta-\sum_{i=0}^{2k}C_{2k}^{i}\delta^{i}\sigma_{ab}\delta^{2k-i}\sigma^{ab}-\xi^{a}\xi^{b}\delta^{2k}T_{ab}=0\,,\end{aligned}$
(101)

where we have used the condition (LABEL:optimaln) with $n=2k-1$ and
$\delta^{2k}T_{ab}\xi^{a}\xi^{b}=0$ on $\mathcal{H}$. Considering the fact
that the perturbation vanishes at early times, this implies
$\delta^{2k+2}\vartheta=0$ on $\mathcal{H}$. Together with the condition
(LABEL:optimaln) of $n=2k-1$, we can easily to see that saturation of the
first $(2k)$ order perturbation inequalities is equivalent to the condition
(LABEL:optimaln) with $n=2k$. We have completed the proof of the case with
$n=2k$.

(2) Case of $n=2k+1$

Next, we will prove Eq. (86) and (LABEL:optimaln) with $n=2k+1$ when the
condition (LABEL:optimaln) with $n=2k$ is satisfied. Same calculation can give

$\displaystyle\begin{aligned}
\delta^{2k+1}M-\Omega_{H}\delta^{2k+1}J&=\delta^{2k-1}\mathcal{E}_{\mathcal{H}}(g,\delta
g)+\delta^{2k-1}\mathcal{E}_{\Sigma_{1}}(g,\delta
g)+\int_{H}\tilde{\bm{\epsilon}}k^{a}\xi^{b}\delta^{2k+1}T_{ab}\\\
&\geq-\frac{\kappa}{8\pi}\delta^{2k+1}A_{(2k+1)}+\delta^{2k-1}\mathcal{E}_{\mathcal{H}}(g,\delta
g)\,,\end{aligned}$ (102)

where we have used the null energy condition of the perturbation matter field
$\delta^{2k+1}T_{ab}\xi^{a}\xi^{b}\geq 0$. For the last term, we have

$\displaystyle\begin{aligned}
\delta^{2k-1}\mathcal{E}_{\mathcal{H}}=-\frac{1}{16\pi}\delta^{2k-1}\int_{\mathcal{H}}\bar{\bm{\epsilon}}\left(\mathcal{L}_{\xi}\delta\gamma_{bc}P^{zbcdef}\nabla_{d}\delta
g_{ef}-\delta\gamma_{bc}P^{zbcdef}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef}\right)\,.\end{aligned}$ (103)

Using the condition (LABEL:optimaln) with $n=2k$, especially

$\displaystyle\begin{aligned}
&\delta^{j}g_{ab}=\delta^{j}\gamma_{ab}=0\,,\quad\quad\forall 1\leq j\leq
k\,,\\\
&g^{ab}\mathcal{L}_{\xi}\delta^{i}\gamma_{ab}=\gamma^{ab}\mathcal{L}_{\xi}\delta^{i}\gamma_{ab}=0\,,\quad\forall\,i\leq
2k\end{aligned}$ (104)

on $\mathcal{H}$, Using these conditions as well as a same calculation as Eqs.
(93) and (95), we can further obtain

$\displaystyle\begin{aligned}
\delta^{m}[g^{bc}\mathcal{L}_{\xi}\delta\gamma_{bc}]=\delta^{m}[g^{bc}\delta\gamma_{bc}]=0\,,\quad\quad\forall
m\leq 2k-1\end{aligned}$ (105)

on $\mathcal{H}$, and therefore

$\displaystyle\begin{aligned} \delta^{m}[P^{z}{}_{ij}{}^{def}\nabla_{d}\delta
g_{ef}]=\delta^{m}[P^{z}{}_{ij}{}^{def}\nabla_{d}\mathcal{L}_{\xi}\delta
g_{ef}]=0\,,\quad\quad\forall m\leq k-1\end{aligned}$ (106)

on $\mathcal{H}$. Using the above results, from Eq. (103), it is not hard to
see

$\displaystyle\begin{aligned} \delta^{2k-1}\mathcal{E}_{\mathcal{H}}(g,\delta
g)=0\,.\end{aligned}$ (107)

Finally, the $(2k+1)$th-order perturbation inequality can be written as

$\displaystyle\begin{aligned}
\delta^{2k+1}M-\Omega_{H}\delta^{2k+1}J\geq-\frac{\kappa}{8\pi}\delta^{2k+1}A_{(2k+1)}\,.\end{aligned}$
(108)

Saturations this inequality demands the additional condition

$\displaystyle\delta^{(2k+1)}T_{ab}\xi^{a}\xi^{b}=0$ (109)

on $\mathcal{H}$. Using the condition (LABEL:optimaln) with $n=2k$ and the
Raychaudhuri, straight forward calculation can show that the above condition
also gives $\delta^{(2k+1)}\vartheta=0$ and finally proved that the saturation
of the first $(2k+1)$ order perturbation inequality is equivalent to the
condition (LABEL:optimaln) with $n=2k+1$, i.e., we have completed the proof.
$\Box$

Finally, we turn to calculate the quantity $\kappa\delta^{n}\bar{A}_{(n)}$,
where we define $\bar{A}=A/8\pi$. From the expression (24) of the area $A$, we
have the identity

$\displaystyle\begin{aligned}
\bar{A}_{(n)}^{2}(\alpha)-2\bar{A}_{(n)}(\alpha)M_{(n)}^{2}(\alpha)-J_{(n)}^{2}(\alpha)=0\,.\end{aligned}$
(110)

in which these quantity is evaluated in the configuration
$\phi_{(n)}(\alpha)$. In this case, we have the setup such that
$\delta^{i}M_{(n)}=\delta^{i}M$ and $\delta^{i}J_{(n)}=\delta^{i}J\quad\forall
i\leq n-1$, and thus we have

$\displaystyle\begin{aligned}
\delta^{i}\bar{A}=\delta^{i}\bar{A}_{(n)}\,,\quad\quad\forall i\leq
n-1\,.\end{aligned}$ (111)

The saturation condition of the first $(n-1)$ order perturbation inequalities
$\delta^{i}\vartheta=0$ for any $i\leq n-1$ implies that
$\delta^{i}\bar{A}=0\,,\forall i\leq n-1$, i.e., we also have
$\delta^{i}\bar{A}_{(n)}=0\,,\forall i\leq n-1$. Using the above results, we
can further obtain

$\displaystyle\begin{aligned}
\frac{\kappa}{8\pi}\delta^{n}A_{(n)}=\frac{1}{2}\sum_{i=1}^{n-1}C_{n}^{i}\left(\delta^{i}M\delta^{n-i}M-\frac{1}{2(1+\epsilon)}\delta^{i}J\delta^{n-2}J\right)\,.\end{aligned}$
(112)

####  4.5.2 Gedanken experiments to overspun the black hole at higher-order
approximation

Figure 2: Plot showing the value of $A_{2m}$ with different integer $m$

After obtaining the above results, it is not hard to extend the discussion
into any arbitrary high-order approximations. With a straightforward
extension, we evaluated $h(\lambda)$ under the first $100$ order approximation
of $\lambda$ and $\epsilon$. After considering the perturbation inequality and
optimal conditions(i.e., the lower approximation of $h(\lambda)$ vanishes), we
find that the $n$th order approximation of $h(\lambda)$ can always be
expressed as

$\displaystyle\begin{aligned} &h(\lambda)\geq 0\,,\quad\forall\,n=2k+1\\\
&h(\lambda)\geq
B_{n/2}(\lambda^{n/2}\delta^{n/2}J-A_{n/2}\epsilon^{n/2})^{2}\geq
0\quad\forall\,n=2k\,,\end{aligned}$ (113)

with the non-negative integer $k$. Here $B_{n/2}$ is a positive quantity and
$A_{2m+1}=0$ for any integer $m\geq 1$. In Fig.2, we show the value of the
parameter $A_{2m}$ of different integer $m$. These results show that the
nearly extremal Kerr black hole cannot be overspun by the perturbation matter
fields with the null energy condition under the higher-order approximation of
perturbation. It indicates that the WCCC in Kerr black holes might be
straightly valid under the perturbation process.

##  5 Conclusion

In this paper, we extended the Sorce-Wald gedanken experiments into the
higher-order approximations to study the WCCC of a Kerr black hole under the
perturbation process. First of all, under the assumptions that the spacetime
is stable under the perturbation and the matter fields satisfy the null energy
condition, we derived the first four order perturbation inequalities under the
corresponding lower optimal conditions. After considering these constraints,
we showed that the nearly Kerr black hole cannot be overspun under the fourth-
order approximation of perturbation. Then, based on the mathematical
induction, we proved that the $n$th-order perturbation inequality can be
generally expressed as

$\displaystyle\begin{aligned}
\delta^{n}M-\Omega_{H}\delta^{n}J\geq-\frac{\kappa}{8\pi}\delta^{n}A_{(n)}\end{aligned}$
(114)

when the first $(n-1)$th order perturbation inequalities are saturated. Using
the perturbation inequalities, we extended the Sorce-Wald discussion to the
first 100 order approximation. Our result showed that the WCCC is always valid
in the perturbation process as long as the matter fields satisfy the null
energy condition, which implies that the WCCC might be strictly true under the
non-approximation level of the perturbation.

## References

  * (1) R. Penrose, Gravitational collapse: The role of general relativity, Riv. Nuovo Cim. 1 , 252-276 (1969).
  * (2) R. M. Wald, “Gravitational collapse and cosmic censorship,” 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, “Overcharging a black hole and cosmic censorship,” Phys. Rev. D 59, 064013 (1999) doi:10.1103/PhysRevD.59.064013 [arXiv:gr-qc/9808043 [gr-qc]].
  * (5) F. de Felice and Y. Q. Yu, “Turning a black hole into a naked singularity,” Class. Quant. Grav. 18, 1235-1244 (2001) doi:10.1088/0264-9381/18/7/307.
  * (6) T. Jacobson and T. P. Sotiriou, “Destroying black holes with test bodies,” J. Phys. Conf. Ser. 222, 012041 (2010) doi:10.1088/1742-6596/222/1/012041 [arXiv:1006.1764 [gr-qc]].
  * (7) G. Chirco, S. Liberati and T. P. Sotiriou, “Gedanken experiments on nearly extremal black holes and the Third Law,” Phys. Rev. D 82, 104015 (2010) doi:10.1103/PhysRevD.82.104015 [arXiv:1006.3655 [gr-qc]].
  * (8) A. Saa and R. Santarelli, “Destroying a near-extremal Kerr-Newman black hole,” Phys. Rev. D 84, 027501 (2011) doi:10.1103/PhysRevD.84.027501 [arXiv:1105.3950 [gr-qc]].
  * (9) 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]].
  * (10) J. Sorce and R. M. Wald, “Gedanken experiments to destroy a black hole. II. Kerr-Newman black holes cannot be overcharged or overspun,” Phys. Rev. D 96, no.10, 104014 (2017) doi:10.1103/PhysRevD.96.104014 [arXiv:1707.05862 [gr-qc]].
  * (11) V. Iyer and R. M. Wald, “Some properties of Noether charge and a proposal for dynamical black hole entropy,” Phys. Rev. D 50, 846-864 (1994) doi:10.1103/PhysRevD.50.846 [arXiv:gr-qc/9403028 [gr-qc]].
  * (12) X. Y. Wang and J. Jiang, “Gedanken experiments at high-order approximation: nearly extremal Reissner-Nordström black holes cannot be overcharged,” JHEP 05, 161 (2020) doi:10.1007/JHEP05(2020)161 [arXiv:2004.12120 [hep-th]].
  * (13) S. Hollands and R. M. Wald, “Stability of Black Holes and Black Branes,” Commun. Math. Phys. 321, 629-680 (2013). doi:10.1007/s00220-012-1638-1 [arXiv:1201.0463 [gr-qc]].