††institutetext: Science Institute, University of Iceland
Dunhaga 3, 107 Reykjavík, Iceland.

# Krylov Complexity of Open Quantum Systems: From Hard Spheres to Black Holes

Vyshnav Mohan<EMAIL_ADDRESS>
###### Abstract

We examine the complexity of quasi-static chaotic open quantum systems. As a
prototypical example, we analytically compute the Krylov complexity of a
slowly leaking hard-sphere gas using Berry’s conjecture. We then connect it to
the holographic complexity of a $d+1$-dimensional evaporating black hole using
the Complexity=Volume proposal. We model the black hole spacetime by stitching
together a sequence of static Schwarzschild patches across incoming negative
energy null shock waves. Under certain identification of parameters, we find
the late time complexity growth rate during each quasi-static equilibrium to
be the same in both systems.

##  1 Introduction

In recent years, the complexity of quantum systems has shown much promise as a
tool in diagnosing quantum chaos Parker:2018yvk ; Rabinovici:2022beu ;
Espanol:2022cqr ; Hashimoto:2023swv ; Camargo:2023eev . Complexity tracks the
degrees of freedom of the system, especially at very long time scales where
other quantum information theoretic measures like entanglement entropy have
saturated Susskind:2014moa . This feature has been particularly fruitful in
theories involving black holes, where the volume of spacelike extremal
codimension-1 surfaces has been conjectured to measure the complexity of the
corresponding state in the dual boundary theory Susskind:2014rva ;
Stanford:2014jda . The proposal, dubbed the Complexity=Volume (CV)
prescription, adds a non-trivial entry to the holographic dictionary.

In a quantum mechanical system, we can quantify complexity using Krylov
Complexity (or K-complexity) Parker:2018yvk ; Barbon:2019wsy , which measures
the growth of operators by treating its evolution as a particle hopping on a
semi-infinite chain111For applications of Krylov complexity to various
systems, refer to Dymarsky:2021bjq ; Avdoshkin:2022xuw ; Camargo:2022rnt ;
Caputa:2021ori ; Adhikari:2022whf ; Kundu:2023hbk ; Patramanis:2023cwz ;
Adhikari:2022oxr ; Bhattacharjee:2022vlt ; Vasli:2023syq ; Rabinovici:2023yex
; Jian:2020qpp ; Erdmenger:2023wjg ; Hornedal:2022pkc .. Under Hamiltonian
evolution, a wavefunction centered around the first site will spread deeper
into the chain. Krylov complexity is defined as the average position of the
particle on the chain as a function of time, thereby capturing the “spread” of
the operator in the operator space.

In section 2, we will use this formulation of operator complexity to study the
growth of operators in a slowly leaking hard sphere gas. We will assume that
the gas is leaking from a small box into a bigger box, as in Krishnan:2021faa
. Moreover, we will work in the semiclassical, low density limit of the
system. A box of hard sphere gas is classically chaotic and satisfies Berry’s
conjecture, which states that the high-lying energy eigenstates behave as if
they have been picked from a Gaussian ensemble berry1983semiclassical ;
berry1991chaos ; Srednicki:1994mfb . Berry’s conjecture played a crucial role
in arriving at Eigenstate Thermalization Hypothesis (ETH) Srednicki:1994mfb
(see also PhysRevA.43.2046 ) and makes the system analytically tractable.

However, since we are dealing with an open quantum system, the canonical
Krylov complexity calculations would not work here. This is because the
evolution is non-unitary, and the Lanczos algorithm one usually employs to
calculate the K-complexity leads to unsatisfactory results
Bhattacharya:2022gbz ; Liu:2022god ; Bhattacharjee:2022lzy ;
Bhattacharya:2023zqt . We will sidestep this difficulty by working with a
slowly leaking gas. This allows us to focus on a time period, which we will
refer to as an epoch, during which the boxes equilibrate separately, and there
is no overall exchange of particles. During an epoch, the Lanczos algorithm
gives meaningful results, and we can compute the Krylov complexity using
standard techniques. Moreover, we will also show that we can patch together
adjacent epochs to form a continuous curve.

Using ETH and very mild assumptions on the off-diagonal matrix elements of the
operators, we will argue how the hard sphere calculation can be generalized to
any chaotic open quantum system. This prompts us to look at our slowly leaking
hard sphere gas model as an excellent prototype where analytic calculations
can be carried out explicitly.

In section 3, we will use the CV prescription to calculate the holographic
complexity of an evaporating black hole. We will model an evaporating black
hole by stitching together a sequence of static Schwarzschild spacetimes
across incoming negative energy null shock waves. Each Schwarzschild patch is
characterized by a constant mass that decreases as we go across the shock
waves. These patches correspond to periods where the black hole is effectively
not evaporating. Therefore, we will refer to them as epochs, akin to our
slowly leaking gas calculation. We will calculate the volume of boundary-
anchored extremal codimension-1 surfaces in this background. Under an
identification of parameters, we show that the late time rate of growth of
complexity during an epoch matches the slowly leaking gas calculation.

##  2 Krylov Complexity of a Slowly Leaking Gas

Consider an operator $O$ that acts on the states of a quantum mechanical
system. If the Hamiltonian of the system is $H$, then the time evolution of
this operator is given by

$\displaystyle O(t)=e^{iHt}Oe^{-iHt}=e^{i\mathcal{L}t}O$ (1)

where $\mathcal{L}(O)=[H,O]$ is the Liouvillian superoperator. Krylov
complexity measures the spread of $O(t)$ in the Krylov subspace, the Hilbert
space spanned by $\mathcal{L}^{n}O$. We will not provide a pedagogical review
of Krylov complexity here as detailed reviews can be found elsewhere
Barbon:2019wsy ; Parker:2018yvk ; Caputa:2021sib .

The microcanonical refinement of Krylov complexity, as introduced in
Kar:2021nbm , will be our focus for analysis. Additionally, we will use the
moment method to compute Krylov complexity. Let us quickly review this
construction. Our starting point is the thermal two-point function of the
operator:

$\displaystyle
G(t)=\sum_{i,j}e^{-\frac{\beta}{2}\left(E_{i}+E_{j}\right)}e^{it\left(E_{i}-E_{j}\right)}\left|\left\langle
E_{i}\left|O\right|E_{j}\right\rangle\right|^{2}$ (2)

Here $\beta$ is the inverse temperature, and $E_{i,j}$ are the energy
eigenvalues of the system. Approximating the sum by an integral over the
density of eigenstates $\rho(E)$, we get

$\displaystyle G(t)=\int_{0}^{\infty}dEe^{-\beta
E}\int_{-2E}^{2E}d\omega\rho\left(E+\frac{\omega}{2}\right)\rho\left(E-\frac{\omega}{2}\right)\left|\left\langle
E+\frac{\omega}{2}\left|O\right|E-\frac{\omega}{2}\right\rangle\right|^{2}e^{i\omega
t}$ (3)

where we have defined the average energy and the energy difference as follows:

$\displaystyle E=\frac{E_{i}+E_{j}}{2}\quad\text{and}\quad\omega=E_{i}-E_{j}.$
(4)

The Liouvillian is sensitive only to the energy differences $\omega$, and it
does not mix different average energy $E$ sectors Kar:2021nbm . Therefore, we
can work with a fixed $E$ and then average over all the other energy sectors
at the end of the calculation. The fixed energy two point function can be
obtained by taking its inverse Laplace transform:

$\displaystyle
G_{E}(t)=\int_{-2E}^{2E}d\omega\rho\left(E+\frac{\omega}{2}\right)\rho\left(E-\frac{\omega}{2}\right)\left|\left\langle
E+\frac{\omega}{2}\left|P_{1}\right|E-\frac{\omega}{2}\right\rangle\right|^{2}e^{i\omega
t}.$ (5)

The key elements in our analysis are the moments of these two point functions.
They are given by

$\displaystyle\mu_{n}^{E}=\frac{\left.\left(-i\frac{d}{dt}\right)^{n}G_{E}(t)\right|_{t=0}}{G_{E}(0)}$
(6)

Using the Hankel transformation matrix $M_{ij}=\mu_{i+j}^{E}$, we can
immediately calculate the Lanczos coefficients $b_{n}^{E}$:

$\displaystyle\left(b_{1}^{E}\right)^{2n}\left(b_{2}^{E}\right)^{2n-2}\ldots\left(b_{n}^{E}\right)^{2}=\operatorname{det}\left[M_{ij}\right]_{0\leq
i,j\leq n}$ (7)

The Lanczos coefficients are handy objects because they contain all the
information about the dynamics of the operator $O$. Moreover, they completely
determine the Krylov complexity of the operator. To see this, let us note that
the fixed energy Krylov complexity is given in terms of the $K$-wavefunctions
$\phi_{E,n}$ as follows Kar:2021nbm

$\displaystyle K^{E}(t)=\sum_{n=0}^{D_{O}-1}n\left|\phi_{E,n}(t)\right|^{2},$
(8)

where $D_{O}$ is the dimensionality of the Krylov subspace. The
$K$-wavefunctions are, in turn, related to the Lanczos coefficients through
the Schrödinger equation:

$\displaystyle\dot{\phi}_{E,n}(t)=b_{n+1}^{E}\phi_{E,n+1}(t)-b_{n}^{E}\phi_{E,n-1}(t).$
(9)

Using the initial condition $\phi_{E,n}(0)=\delta_{n0}$, we can solve for the
$K$-wavefunctions and compute $K^{E}$.

Krylov Complexity $K^{E}$ measures the complexity growth of operators within a
specific energy sector. Operationally, this definition of complexity
corresponds to the choice of a “microcanonical” inner product on the Krylov
subspace:

$\displaystyle(A|B)_{E}=\int_{-2E}^{2E}d\omega\rho\left(E+\frac{\omega}{2}\right)\rho\left(E-\frac{\omega}{2}\right)\left\langle
E+\frac{\omega}{2}\left|A^{*}\right|E-\frac{\omega}{2}\right\rangle\left\langle
E+\frac{\omega}{2}\left|B\right|E-\frac{\omega}{2}\right\rangle.$

If we had instead used a thermal inner product that involved integrating over
various energy sectors, there would have been a risk of the signatures of
chaos being washed out (See Kar:2021nbm  for details). Therefore, the natural
thing to do is to compute the growth in a fixed energy sector and then take an
average over various values of $E$. Since different energy sectors are
thermally populated, we can define the thermal Krylov complexity by taking a
Laplace transform of its fixed energy counterpart Kar:2021nbm :

$\displaystyle K_{th}(t)=\frac{\int_{0}^{\infty}dEe^{-\beta
E}\mathcal{C}(E)K^{E}(t)}{\int_{0}^{\infty}dEe^{-\beta E}\mathcal{C}(E)}$ (11)

where the normalization constant $\mathcal{C}(E)$ is given by

$\displaystyle\mathcal{C}(E)=\int_{-2E}^{2E}d\omega\rho\left(E+\frac{\omega}{2}\right)\rho\left(E-\frac{\omega}{2}\right)\left|\left\langle
E+\frac{\omega}{2}\left|O\right|E-\frac{\omega}{2}\right\rangle\right|^{2}.$
(12)

Our primary focus is on the thermal Krylov complexity. In the following
subsections, we will use the moments of the thermal two point functions to
calculate $K_{th}(t)$ of a slowly leaking hard sphere gas.

###  2.1 Warm-up: Single Box

Before we look at the slowly leaking gas, it is instructive to look at a
single box of hard sphere gas. Consider a cubic box of edge length $L+2a$
enclosing $N$ hard spheres. We will assume that the hard spheres are identical
and have radius $a$. The classical Hamiltonian of the system is given by

$\displaystyle
H=\sum_{i=1}^{N}\frac{\mathbf{p}_{i}^{2}}{2m}+\sum_{i<j}V\left(\left|\mathbf{x}_{i}-\mathbf{x}_{j}\right|\right)$
(13)

where

$\displaystyle V(r)=\begin{cases}+\infty&\text{ for }r<2a\\\ 0&\text{ for
}r>2a\end{cases}$ (14)

This system is classically chaotic and shows eigenstate thermalization when
treated quantum mechanically Srednicki:1994mfb . Let us denote the energy
eigenfunctions of the system by $\psi(X)$, where $X$ is the $3N$-dimensinal
position vector. The wavefunctions are defined on the domain

$\displaystyle
D=\left\\{\mathrm{x}_{1},\ldots,\mathrm{x}_{N}\big{|}x_{i1,2,3}\in\left[-\frac{1}{2}L,\frac{1}{2}L\right];\left|\mathbf{x}_{i}-\mathbf{x}_{j}\right|\geq
2a\right\\}$ (15)

We will impose the boundary condition that $\psi(X)$ vanishes on $\partial D$.
The wavefunctions with energy $E_{\alpha}$ can be chosen to be of the
following form Srednicki:1994mfb :

$\displaystyle\psi_{\alpha}(\mathbf{X})=\mathcal{N}_{\alpha}\int
d^{3N}PA_{\alpha}(\mathbf{P})\delta\left(\mathbf{P}^{2}-2mE_{\alpha}\right)\exp(i\mathbf{P}\cdot\mathbf{X}/\hbar)$
(16)

where $\mathcal{N}_{\alpha}$ is the normalization constant. We can choose the
wavefunction to be everywhere real by imposing
$A_{\alpha}^{*}(\mathbf{P})=A_{\alpha}(-\mathbf{P})$.

Let us define the thermal wavelength $\lambda$ of the system as follows:

$\displaystyle\lambda=\sqrt{\frac{2\pi\hbar^{2}}{mkT_{\alpha}}}$ (17)

where the termperature $T_{\alpha}$ is related to the energy through the
relation $E_{\alpha}=\frac{3}{2}NkT_{\alpha}$. When the thermal wavelength
$\lambda\lesssim a$, then the high-lying energy eigenstates are expected to
satisfy the Berry’s conjecture, which states that $A_{\alpha}(\mathbf{P})$ can
be regarded as a Gaussian variable with the two-point function

$\displaystyle\left\langle
A_{\alpha}(\mathbf{P})A_{\beta}\left(\mathbf{P}^{\prime}\right)\right\rangle_{\mathrm{EE}}=\delta_{\alpha\beta}\frac{\delta^{3N}\left(\mathbf{P}+\mathbf{P}^{\prime}\right)}{\delta\left(\mathbf{P}^{2}-\mathbf{P}^{\prime
2}\right)}.$ (18)

The subscript $\mathrm{EE}$ stands for eigenstate ensemble, the fictitious
Gaussian ensemble the high-lying wavefunctions can be thought to belong to.
The wavefunction in the momentum space is given by:

$\displaystyle\widetilde{\psi}_{\alpha}(\mathbf{P})\equiv
h^{-3N/2}\int_{D}d^{3N}X\psi_{\alpha}(\mathbf{X})\exp(-i\mathbf{P}\cdot\mathbf{X}/\hbar)$
(19)

We will assume that we always work in the regime where Berry’s conjecture
holds ($\lambda\lesssim a$). Moreover, we will also assume that the density of
the gas is low, that is, $Na\ll L^{3}$. These assumptions give us enormous
analytic control over the system. In particular, it is easy to see that the
averaged two point functions are given by Srednicki:1994mfb

$\displaystyle\left\langle\widetilde{\psi}_{\alpha}^{*}(\mathbf{P})\widetilde{\psi}_{\beta}\left(\mathbf{P}^{\prime}\right)\right\rangle_{\mathrm{EE}}=\delta_{\alpha\beta}\mathcal{N}_{\alpha}^{2}h^{3N}\delta\left(\mathbf{P}^{2}-2mE_{\alpha}\right)\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right)$
(20)

where

$\displaystyle\delta_{D}^{3N}(\mathbf{K})\equiv
h^{-3N}\int_{D}d^{3N}X\exp(i\mathbf{K}\cdot\mathbf{X}/\hbar).$ (21)

###  2.2 Krylov Complexity from Moments

Consider the operator $P_{1}$ that measures the momentum of one of the
particles. The matrix elements of this operator in the energy eigenbasis are
given by

$\displaystyle\langle E_{m}|P_{1}|E_{n}\rangle=\int d^{3}p_{1}d^{3}p_{2}\cdots
d^{3}p_{N}\ |p_{1}|\
\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}_{n}(\textbf{P}),$ (22)

where $|p_{1}|$ is the magnitude of the momenta of one of the particles. Let
us calculate the Krylov complexity of this operator by using the moments of
the thermal two point function. Using (6), we get

$\displaystyle\mu_{2n}^{E}=\frac{1}{\mathcal{C}(E)}\int_{-2E}^{2E}d\omega\
\rho\left(E+\frac{\omega}{2}\right)\rho\left(E-\frac{\omega}{2}\right)\left|\left\langle
E+\frac{\omega}{2}\left|P_{1}\right|E-\frac{\omega}{2}\right\rangle\right|^{2}\omega^{2n}$
(23)

$\mathcal{C}(E)$ is the normalization constant we defined in (12). The density
of eigenstates $\rho$ is given by Srednicki:1994mfb

$\displaystyle\rho(E)=\frac{1}{\Gamma(3N/2)E}\left(\frac{mL^{2}E}{2\pi\hbar^{2}}\right)^{\frac{3N}{2}}.$
(24)

Now let us calculate the average of the moments in the eigenstate ensemble:

$\displaystyle\left\langle\mu_{2n}^{E}\right\rangle_{\text{EE}}=\frac{1}{\mathcal{C}(E)}\int_{-2E}^{2E}d\omega\
\rho_{0}(E,\omega)\left|\left\langle
E+\frac{\omega}{2}\left|P_{1}\right|E-\frac{\omega}{2}\right\rangle\right|^{2}_{\text{EE}}\omega^{2n}$
(25)

Here we have used the shorthand notation $\rho_{0}(E,\omega)$ for the product
of the density of states. Using (22), we see that:

$\displaystyle\left|\left\langle
E_{m}\left|P_{1}\right|E_{\ell}\right\rangle\right|^{2}_{\text{EE}}=\int
d^{3N}Pd^{3N}P^{\prime}\
|p_{1}||p^{\prime}_{1}|\left\langle\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}_{\ell}(\textbf{P})\widetilde{\psi}^{*}_{\ell}(\textbf{P}^{\prime})\widetilde{\psi}_{m}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}$
(26)

The four-point function can be broken down into two point functions using Wick
contractions:

$\displaystyle\begin{aligned}
\left\langle\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}_{\ell}(\textbf{P})\widetilde{\psi}^{*}_{\ell}(\textbf{P}^{\prime})\widetilde{\psi}_{m}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}&=\left\langle\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}_{\ell}(\textbf{P})\right\rangle_{\text{EE}}\left\langle\widetilde{\psi}^{*}_{\ell}(\textbf{P}^{\prime})\widetilde{\psi}_{m}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}\\\
&+\left\langle\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}^{*}_{\ell}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}\left\langle\widetilde{\psi}_{\ell}(\textbf{P})\widetilde{\psi}_{m}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}\\\
&+\left\langle\widetilde{\psi}_{\ell}(\textbf{P})\widetilde{\psi}^{*}_{\ell}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}\left\langle\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}_{m}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}\end{aligned}$
(27)

From (20), we can see that contracting two eigenfunctions will produce a delta
function in its indices. The only terms in (27) which would contribute to the
moment calculation are the ones without any $\delta_{m\ell}$ factor. This is
because (25) contains a factor $(E_{\ell}-E_{m})$, multpyling the four-point
functions. Therefore, only the last term in (27) would contribute. Let us look
at the integral of this term separately:

$\displaystyle\Phi_{m\ell}=\int d^{3N}Pd^{3N}P^{\prime}\
p_{1}p^{\prime}_{1}\left\langle\widetilde{\psi}_{\ell}(\textbf{P})\widetilde{\psi}^{*}_{\ell}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}\left\langle\widetilde{\psi}^{*}_{m}(\textbf{P})\widetilde{\psi}_{m}(\textbf{P}^{\prime})\right\rangle_{\text{EE}}$
(28)

Using (20), we see that

$\displaystyle\Phi_{m\ell}=\mathcal{N}_{m}^{2}\mathcal{N}_{\ell}^{2}(Lh)^{3N}\int
d^{3N}Pd^{3N}P^{\prime}\
p_{1}p^{\prime}_{1}\delta\left(\mathbf{P}^{2}-2mE_{\ell}\right)\delta\left(\mathbf{P}^{\prime}{}^{2}-2mE_{m}\right)\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right)$
(29)

where we have used Srednicki:1994mfb

$\displaystyle(\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right))^{2}=(L/h)^{3N}\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right).$
(30)

Now let us look at the $m=\ell$ component. This corresponds to the case where
$\omega=0$ and $E_{m}=E_{\ell}=E$. In the low density limit, we can
essentially replace
$\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right)$ with a dirac
delta. This gives us

$\displaystyle\Phi_{mm}=\mathcal{N}^{4}(Lh)^{3N}\int d^{3N}P\
p_{1}^{2}\delta\left(\mathbf{P}^{2}-2mE\right)$ (31)

Choosing the normalization constant as in Srednicki:1994mfb

$\displaystyle\mathcal{N}^{-2}\equiv\mathcal{N}_{i}^{-2}=L^{3N}\frac{(2\pi
mE)^{\frac{3N}{2}},}{\Gamma(3N/2)E}$ (32)

we get

$\displaystyle\begin{aligned} \Phi_{mm}&=\mathcal{N}^{2}h^{3N}\int d^{3}p_{1}\
p_{1}^{2}\left(2\pi mkT\right)^{-3/2}e^{-\mathbf{p}_{1}^{2}/2mkT}\\\
&=4\pi\mathcal{N}^{2}h^{3N}\int dp_{1}\ p_{1}^{4}\left(2\pi
mkT\right)^{-3/2}e^{-\mathbf{p}_{1}^{2}/2mkT}\\\
&=3mkT\mathcal{N}^{2}(h)^{3N}\\\
&=(h/L)^{3N}3mkT\frac{\Gamma(3N/2)(2mE)}{\left(2m\pi E\right)^{3N/2}}\\\
&\equiv\Phi_{E}\end{aligned}$ (33)

Now let us return to the non-diagonal elements of $\Phi_{m\ell}$. From our
definition (21), it is easy to see that
$\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right)$ is a sharply
peaked function that is zero almost everywhere. Let us choose this function to
be a Gaussian distribution as in Srednicki:1994mfb :

$\displaystyle\delta_{D}^{3N}\left(\mathbf{P}-\mathbf{P}^{\prime}\right)\simeq(L/h)^{3N}\exp\left[-\left(\mathbf{P}-\mathbf{P}^{\prime}\right)^{2}L^{2}/4\pi\hbar^{2}\right]$
(34)

This gives us Srednicki:1994mfb

$\displaystyle\Phi_{ij}\simeq\Phi_{ii}\exp\left[-m\left(E_{i}-E_{j}\right)^{2}L^{2}/8\pi\hbar^{2}E_{i}\right]$
(35)

In the notation $E_{i}=E+\frac{\omega}{2}$ and $E_{j}=E-\frac{\omega}{2}$, the
expression simplifies to

$\displaystyle\Phi_{ij}\simeq\Phi_{E}\exp\left[-\frac{m\omega^{2}L^{2}}{8\pi\hbar^{2}E}\right]$
(36)

Plugging in the expressions, we find that the following integral gives the
moments:

$\displaystyle\left\langle\mu_{2n}^{E}\right\rangle_{\text{EE}}=\frac{1}{\mathcal{C}(E)}\int_{-2E}^{2E}d\omega\
\rho_{0}\left(E,\omega\right)\omega^{2n}\Phi_{E}\exp\left[-\frac{m\omega^{2}L^{2}}{8\pi\hbar^{2}E}\right]$
(37)

####  2.2.1 A Tale of Two Saddles

We will compute the moment integral (37) using a saddle point approximation.
We claim that for small moments $n$, the saddle point $\omega^{*}$ is located
at

$\displaystyle\omega^{*}\ll 2E.$ (38)

To see this, let us work out the saddle point equation up to the leading order
in $\omega/2E$. This gives us

$\displaystyle\frac{2n}{\omega_{*}}-\frac{m\omega^{*}L^{2}}{4\pi\hbar^{2}E}=0$
(39)

Solving the above equation, we find that the saddle point is located at

$\displaystyle\omega_{*}=\sqrt{\frac{8\lambda^{2}E^{2}n}{3L^{2}N}}$ (40)

where $\lambda$ is the thermal wavelength at the energy $E$. When
$n\ll\frac{3NL^{2}}{2\lambda^{2}}$, $\omega^{*}$ is much smaller than $2E$, as
advertised. Substituting the saddle point in (37), we get:

$\displaystyle\left\langle\mu_{2n}^{E}\right\rangle_{\text{EE}}\simeq\frac{1}{\mathcal{C}(E)}\rho_{0}\left(E,\omega_{*}\right)\
\omega_{*}^{2n}\ \Phi_{E}\
\exp\left[-\frac{m\omega_{*}^{2}L^{2}}{8\pi\hbar^{2}E}\right]$ (41)

When $n$ is sufficiently large, we can calculate the Lanczos coefficients
using the relation Barbon:2019wsy ; Parker:2018yvk ; Kar:2021nbm

$\displaystyle\mu_{2n}\sim\left(b_{n}\right)^{2n}e^{o(n)}.$ (42)

We find that

$\displaystyle\left\langle
b_{n}^{E}\right\rangle_{\text{EE}}\sim\sqrt{\frac{8\lambda^{2}E^{2}n}{3L^{2}N}}.$
(43)

This behavior is termed “Lanczos ascent” in literature Kar:2021nbm . When
$n\gg\frac{3NL^{2}}{2\lambda^{2}}$, $\omega_{*}\gg 2E$. Since this violates
our assumption (38), (40) ceases to be a saddle point of the moment integral.
For large $n$, (40) gets replaced by a new saddle point located at
$\omega_{*}\simeq 2E$. We can obtain this saddle point mechanically by noting
that if222Note that the $\frac{\lambda}{L}$ is a tiny number.

$\displaystyle\left|\frac{\omega}{2}-E\right|\ll\frac{2\lambda^{2}E}{L^{2}},$
(44)

the saddle point equation is given by

$\displaystyle\frac{d}{d\omega}\left(\rho_{0}\left(E,\omega\right)\omega^{2n}\right)=0.$
(45)

Solving the equation, we find that the new saddle point is located at

$\displaystyle\omega_{*}=2E\sqrt{\frac{n}{\frac{3N-2}{2}+n}}.$ (46)

Figure 1: The figure on the left shows the growth of the Lanczos coefficients
w.r.t $n$. There is an initial scrambling phase (shaded in green) where
$b_{n}^{E}$ grows as $\sqrt{n}$, and then it saturates to a constant value.
The figure on the right shows the corresponding Krylov Complexity growth.
During the scrambling phase, the K-complexity grows quadratically. Then it
transitions into a linear growth proportional to the average energy of the
box.

When $n\gg\frac{3NL^{2}}{2\lambda^{2}}$,

$\displaystyle\omega_{*}\simeq 2E.$ (47)

It is easy to verify that (46) satisfies (44) when
$n\gg\frac{3NL^{2}}{2\lambda^{2}}$. Therefore, the saddle point of the
integral changes when we look at higher moments. This has interesting
consequences for the behavior of Lanczos coefficients. In particular, we can
use (42) to see that:

$\displaystyle\left\langle b_{n}^{E}\right\rangle_{\text{EE}}\sim 2E.$ (48)

The saturation of Lanczos coefficients is referred to as “Lanczos plateau”.
This gives us

$\displaystyle\left\langle
b_{n}^{E}\right\rangle_{\text{EE}}\sim\begin{cases}\sqrt{\frac{8\lambda^{2}E^{2}n}{3L^{2}N}},&n<\frac{3NL^{2}}{2\lambda^{2}}\\\
E,&n>\frac{3NL^{2}}{2\lambda^{2}}\end{cases}$ (49)

From Barbon:2019wsy ; Fan:2022xaa ; noh2021operator , we can see that these
Lanczos coefficients result in an initial scrambling phase where the Krylov
complexity grows quadratically. Following the scrambling phase, K-complexity
switches to linear growth. Working out the growth rates and reinstating
factors of $\hbar$, we find that the Krylov complexity is given by

$\displaystyle\left\langle
K^{E}(t)\right\rangle_{\text{EE}}\sim\begin{cases}\frac{8\lambda^{2}E^{2}}{3L^{2}N\hbar^{2}}t^{2},&t<t_{S}\\\
\frac{E}{\hbar}\ t,&t>t_{S}\end{cases}$ (50)

where $t_{S}$ is the scrambling time. We can determine $t_{S}$ by noting that
the dynamics of the operator can be thought of as a particle moving on a one-
dimensional semi-infinite chain. The sites on the chain are labeled by $n$,
and Krylov complexity (8) measures the average position $\langle n\rangle$.
The saddle points change when $n\sim\frac{3NL^{2}}{2\lambda^{2}}\equiv n_{S}$.
Therefore, scrambling time is the amount of time Krylov complexity takes to
reach $n_{S}$. This gives us

$\displaystyle t_{S}\simeq\frac{\beta L^{2}\hbar}{2\lambda^{2}}\ \text{where}\
\beta=\frac{1}{k_{B}T}$ (51)

Using (11), we can compute the thermal K-complexity in the eigenstate
ensemble:

$\displaystyle\left\langle
K_{th}\right\rangle_{\text{EE}}=\frac{\int_{0}^{\infty}dEe^{-\beta
E}\mathcal{C}(E)\left\langle
K^{E}(t)\right\rangle_{\text{EE}}}{\int_{0}^{\infty}dEe^{-\beta
E}\mathcal{C}(E)}$ (52)

It is easy to see that we get

$\displaystyle\left\langle
K_{th}\right\rangle_{\text{EE}}\sim\begin{cases}\frac{8\lambda_{*}^{2}E_{*}^{2}}{3L^{2}N\hbar^{2}}t^{2},&t<t_{S}\\\
\frac{E_{*}}{\hbar}\ t,&t>t_{S}\end{cases}$ (53)

where $E_{*}$ and $\lambda_{*}$ are the average energy and thermal wavelength
of the box. Plotting these functions, we obtain figure 1.

###  2.3 Slowly Leaking Gas

Figure 2: Consider two cubic boxes in contact with each other. We fill the box
on the left with $N$ hard spheres. Let us assume that we are in the semi-
classical limit where we can localize particles. At time $t=0$, we poke a
small hole on their common wall so that the gas leaks slowly into the right
box. Let us denote the instantaneous number of particles in the boxes by
$N_{L}$ and $N_{R}$.

Now let us look at the slowly leaking gas model used in Krishnan:2021faa .
Consider two cubic boxes sharing a common side as in figure 2. Let us assume
that the left box has $N$ hard spheres while the right box is empty. At time
$t=0$, we make a small hole in their shared wall so that the gas leaks slowly
into the right box. We assume that the gas is leaking so slowly that there
exists a time scale over which both the boxes have separately equilibrated. We
will refer to this period as an epoch. During each epoch, the number of
particles in each box, denoted by $N_{L,R}$, remains roughly constant. Since
we are working in the semiclassical limit, we can always use either $N_{L}$ or
$N_{R}$ to characterize each epoch.

####  2.3.1 Krylov Complexity during an epoch

Within an epoch, there is no net exchange of particles between the boxes.
Consequently, when examining the full Hamiltonian, it becomes evident that
during each epoch, the Hamiltonian takes on the following factorized
structure:

$\displaystyle H\simeq H_{L}\otimes\mathds{1}_{R}+\mathds{1}_{L}\otimes H_{R}$
(54)

Consider an operator $P_{1}$, which measures the momentum of one of the
particles in the left box. The operator has the following form

$\displaystyle P_{1}=P_{1,L}\otimes\mathds{1}_{R}.$ (55)

Here $P_{1,L}$ is an operator that acts only on the left box. Now let us look
at the Krylov complexity of this operator. If we are computing the Krylov
complexity w.r.t to the entire system, we will get the results of the previous
section - a scrambling phase followed by linear growth. This is because the
left and right boxes form a closed classically chaotic system, and we
anticipate a version of Berry’s conjecture to be applicable to the eigenstates
of the entire system. However, we are after the Krylov complexity of the left
box alone, which is an open quantum system.

The evolution of an operator $O_{L}$ acting on the left box is given by the
master equation:

$\displaystyle\dot{O_{L}}(t)=\frac{i}{\hbar}\left(\mathcal{L}_{L}+\mathcal{L_{D}}\right)$
(56)

Here, $\mathcal{L}_{L}$ denotes the Liouvillian of the left box, and
$\mathcal{L}_{D}$ represents a dissipative term arising from the interaction
between the two boxes. The dissipative term makes the evolution non-Hermitian,
and there is no consensus on extending Krylov complexity calculations to open
quantum systems Bhattacharjee:2022lzy ; Bhattacharya:2022gbz ; Liu:2022god .
However, the “effective” Hamiltonian has no interaction term during an epoch.
Therefore, the evolution of an operator $O_{L}$ would be controlled only by
$\mathcal{L}_{L}$, or equivalently, the Hamiltonian of the left box. This
allows us to carry over our definition of moments (23) to the two-box system:

$\displaystyle\left(\mu_{2n}^{E}\right)_{L}=\frac{1}{\mathcal{C}(E_{L})}\int_{-2E_{L}}^{2E_{L}}d\omega\
\rho_{0}(E_{L},\omega_{L})\left|\left\langle
E_{L}+\frac{\omega_{L}}{2}\left|O_{L}\right|E_{L}-\frac{\omega_{L}}{2}\right\rangle\right|^{2}\omega_{L}^{2n}$
(57)

where $E_{L}\pm\frac{\omega_{L}}{2}$ are the energy eigenvalues of the left
box. Now let us look at the operator $P_{1}$. The operator acting on the left
box is given by tracing out the degrees of freedom of the right box. This
gives us $\hbox{\rm Tr}_{R}{\left(P_{1}\right)}$. The moments of this operator
are given by

$\displaystyle\left(\mu_{2n}^{E}\right)_{L}=\frac{1}{\mathcal{C}(E_{L})}\int_{-2E_{L}}^{2E_{L}}d\omega\
\rho_{0}(E_{L},\omega_{L})\left|\left\langle
E_{L}+\frac{\omega_{L}}{2}\left|\hbox{\rm
Tr}_{R}{\left(P_{1}\right)}\right|E_{L}-\frac{\omega_{L}}{2}\right\rangle\right|^{2}\omega_{L}^{2n}$
(58)

Using (55), we get

$\displaystyle\left(\mu_{2n}^{E}\right)_{L}=\frac{1}{\mathcal{C}(E_{L})}\int_{-2E_{L}}^{2E_{L}}d\omega\
\rho_{0}(E_{L},\omega_{L})\left|\left\langle
E_{L}+\frac{\omega_{L}}{2}\left|P_{1,L}\right|E_{L}-\frac{\omega_{L}}{2}\right\rangle\right|^{2}\omega_{L}^{2n}$
(59)

which is precisely the same expression we had in (23). Therefore, the
calculations in the previous section will go through - The K-complexity will
grow linearly, following a quadratic scrambling phase.

Figure 3: The figure shows the thermal Krylov complexity of a slowly leaking
gas as a function of time. The complexity continues to rise, but eventually
levels off. The inset figures show how K-complexity grows during two different
epochs. As we move into future epochs, we observe a decrease in the late time
growth rate.

####  2.3.2 Stitching Together Epochs

From (53), we can see that the late time linear growth rate of Krylov
complexity during an epoch is given by the average energy of the left box,
which we will denote by $E_{L}$. Using $E=\frac{3}{2}Nk_{B}T$, we get

$\displaystyle\frac{dK_{th}(t)}{dt}\sim\frac{E_{L}}{\hbar}=\frac{3}{2\hbar}N_{L}k_{B}T$
(60)

Now, we will compute the time dependence of $N_{L}$. Suppose the right box is
sufficiently large for the gas to leak out completely. If the area of the hole
is given by $A$, the leakage rate is given by Schroederbook

$\displaystyle\frac{dN}{dt}=-\frac{A}{2L^{3}}\sqrt{\frac{kT}{m}}N.$ (61)

We can immediately integrate the above equation if we assume that the
temperature of the left box remains constant. This gives us

$\displaystyle\begin{aligned}
N(t)&=Ne^{-\frac{A}{2L^{3}}\sqrt{\frac{kT}{m}}t}\\\ &\equiv
Ne^{-\frac{t}{t_{L}}}\end{aligned}$ (62)

where $t_{L}$ is the leakage time of the system. Treating $N(t)$ as an
instantaneous value, we can integrate (60) to obtain

$\displaystyle
K_{th}(t)\sim\frac{N\sqrt{mkT}L^{3}}{A\hbar}\left(1-e^{-\frac{t}{t_{L}}}\right)$
(63)

Plotting the above function, we get figure 3. Krylov complexity keeps
increasing and eventually levels off when $t=t_{L}$.

##  3 Holographic Complexity of an Evaporating Black Hole

In this section, we will study the holographic complexity of a slowly
evaporating black hole using the Complexity=Volume prescription
Susskind:2014rva ; Stanford:2014jda . To make the calculations tractable, let
us model the black hole by patching together a sequence of $k$ static
Schwarzschild spacetimes across negative energy null shock waves. The metric
of the $d+1$ dimensional black hole is then given by

$\displaystyle ds^{2}=-F(r,v)dv^{2}+2dvdr+r^{2}d\Omega_{d-1}$ (64)

where

$\displaystyle F(r,v)=1-\frac{f(v)}{r^{d-2}}.$ (65)

We will choose $f$ to have the following profile

$\displaystyle f(v)=\begin{cases}\omega_{1}^{d-2},&v<v_{1}\\\ \cdots\\\
\omega_{i}^{d-2},&v_{i-1}<v<v_{i}\\\ \cdots\\\ 0,&v_{k-1}<v<v_{k}\end{cases}$
(66)

The mass of each patch is given by

$\displaystyle M_{i}=\frac{(d-1)\Omega_{d-1}}{16\pi G_{N}}\omega_{i}^{d-2},$
(67)

where

$\displaystyle M_{1}>M_{2}>\cdots>0.$ (68)

For the patched-up spacetime to be a good approximation to an evaporating
black hole, we will assume the width of each patch to be much smaller than the
time scales at which the black hole mass changes considerably. Moreover, we
will also assume that the width is larger than the scrambling time of the
black hole. Each patch corresponds to a period where the black hole is
effectively not evaporating. Therefore, we will adopt the terminology from our
previous section and refer to each patch as an epoch. In particular, the patch
between $v=v_{i-1}$ and $v=v_{i}$ shock waves will be labeled as the $i$-th
epoch. During an epoch, the black hole has a constant mass, and we can rewrite
$F(r,v)$ as $F_{i}(r)$.

###  3.1 Penrose Diagram

Figure 4: Penrose diagram of the black hole spacetime when there are three
epochs. We will assume the mass of each epoch to satisfy the relation
$M_{1}>M_{2}>M_{3}\neq 0$. The horizon is indicated by the 45-degree line. The
shock waves are marked by the red lines. A spacelike surface, indicated by the
teal color, will be disconnected in the diagram. The end points of these
disconnected surfaces follow the ordering (73).

To understand the structure of the spacetime, it is instructive to draw its
Penrose diagram. To simplify the discussion, we will restrict ourselves to
$3+1$-dimensions in this subsection. During an epoch, the tortoise coordinate
is given by

$\displaystyle
r^{*}_{i}(r)=\int^{r}\frac{dr^{\prime}}{F_{i}(r^{\prime})}=r+2G_{N}M_{i}\log\left(\left|\frac{r}{2G_{N}M_{i}}-1\right|\right)$
(69)

The corresponding outgoing Eddington–Finkelstein coordinates $u_{i}$ are
defined as follows

$\displaystyle u_{i}\equiv v-2r^{*}_{i}(r)$ (70)

It is easy to see that $u$ is discontinuous across the boundaries of the
epochs. Therefore, if we go along a continuous curve, the coordinate $u$ has a
“jump” in its value as soon as we cross a shock wave. Consequently, employing
identical $u$ and $v$ coordinates throughout all epochs would render
continuous curves discontinuous in the corresponding Penrose diagram. This
discontinuity is the cost we have to pay to keep the Penrose undeformed
Chapman:2018dem ; Chapman:2018lsv .

During an epoch, the outgoing null Kruskal coordinate $U_{i}$ can be defined
as follows Misner:1973prb :

$\displaystyle U_{i}=\begin{cases}-e^{-\frac{u_{i}}{4G_{N}M}},&\text{Outside
the horizon}\\\ e^{-\frac{u_{i}}{4G_{N}M}},&\text{Inside the
horizon}\end{cases}$ (71)

Since $v$ is globally defined, we can define the outgoing null Kruskal
coordinate $V$ everywhere as

$\displaystyle V=e^{\frac{v}{4G_{N}M}}.$ (72)

Using the definition of Kruskal coordinates, we can see that the horizon of
the epochs will always be at $U_{i}=0$. Now let us figure out what happens to
a continuous curve as it crosses a shock wave. We will denote the radial
coordinate of the surface at the location of the shock wave by $r_{S}$. Using
(70) and (71), it is easy to verify that when $r_{S}$ is either in the
interior or sufficiently far away from the horizon, we have the relation

$\displaystyle U_{i}(r_{S})>U_{i+1}(r_{S})$ (73)

Now let us draw the Penrose diagram of the spacetime. We will use the same $U$
and $V$ coordinates across all the epochs. This gives us figure 4. A connected
surface will be discontinuous in this diagram. We can locate the end points of
these disconnected pieces by using (73).

###  3.2 Complexity=Volume

We can study the growth of complexity using the Complexity=Volume conjecture
Stanford:2014jda ; Susskind:2014rva . Let us assume that there are only three
epochs, with masses satisfying the relation $M_{1}>M_{2}>M_{3}\neq 0$. We will
see that extending our results to an arbitrary number of epochs is
straightforward. The Penrose diagram of this spacetime is given in figure 4.

Consider spherically symmetric spacelike codimension-1 surfaces anchored onto
a cutoff surface in the asymptotic region. We will assume that the cutoff
surface is at $r=r_{\infty}$ and the boundary anchoring time is denoted by
$t$. The volume of this surface is given by

$\displaystyle\mathcal{V}=\Omega_{d-1}\int d\lambda
r^{d-1}\sqrt{-F(r,v)\dot{v}^{2}+2\dot{v}\dot{r}}$ (74)

Here $\Omega_{d-1}$ is the dimensionless area of the $(d-1)$-dimensional unit
sphere. Let us rewrite the volume integral as the following summation

$\displaystyle\mathcal{V}=\sum_{i}\mathcal{V}_{i}$ (75)

where $\mathcal{V}_{i}$ is the volume of the portion of the surface located
within the $i$-th epoch. We can then express each of these terms as follows:

$\displaystyle\mathcal{V}_{i}=\Omega_{d-1}\int d\lambda
r^{d-1}\sqrt{-F_{i}(r)\dot{v}^{2}+2\dot{v}\dot{r}}\equiv\Omega_{d-1}\int
d\lambda\mathcal{L}_{i}$ (76)

We can see that $\mathcal{L}_{i}$ is independent of $v$. Therefore, there is a
conserved quantity associated with the integral, which we will denote by
$E_{i}$:

$\displaystyle
E_{i}=\frac{\partial\mathcal{L}_{i}}{\partial\dot{v}}=\frac{r^{d-1}(\dot{r}-F_{i}\dot{v})}{\sqrt{-F_{i}\dot{v}^{2}+2\dot{v}\dot{r}}}$
(77)

The volume integral is reparametrization invariant, allowing us to choose

$\displaystyle r^{d-1}\sqrt{-F_{i}\dot{v}^{2}+2\dot{v}\dot{r}}=1$ (78)

The maximal volumes are given by extremizing the action in (74). This gives us
the following equations of motion:

$\displaystyle\begin{aligned} E_{i}&=r^{2(d-1)}(\dot{r}-F_{i}(r)\dot{v})\\\
r^{2(d-1)}\dot{r}^{2}&=F_{i}(r)+r^{-2(d-1)}E_{i}^{2}\end{aligned}$ (79)

The part of the surface within an epoch can have a turning point where
$\dot{r}$ vanishes. We will denote this point by $r_{\text{i,min}}$. We will
see later that it is not necessary for the surface to have a turning point in
an epoch. However, these turning points will be crucial in calculating late
time growth rates during an epoch.

To characterize various features of the extremal volume surfaces, it is
convenient to define an effective potential as follows:

$\displaystyle V_{i}(r)=F_{i}(r)r^{2(d-1)}+E_{i}^{2}.$ (80)

In particular, we can obtain the turning point $r_{\text{i,min}}$ from the
zero of the effective potential

$\displaystyle V_{i}(r_{\text{i,min}})=0\quad\implies\quad
F_{i}(r_{\text{i,min}})r_{\text{i,min}}^{2(d-1)}+E_{i}^{2}=0$ (81)

Using the equations of motion, we can rewrite the volume integral as follows:

$\displaystyle\mathcal{V}=\Omega_{d-1}\sum_{i}\int
dr\frac{r^{2(d-1)}}{\sqrt{F_{i}(r)r^{2(d-1)}+E_{i}^{2}}}$ (82)

The Complexity=Volume proposal suggests that the complexity of the black hole
at time $t$ is given by the volume of these extremal surfaces Stanford:2014jda
:

$\displaystyle\mathcal{C}_{i}=\frac{\mathcal{V}(t)}{G_{N}\omega_{i}}$ (83)

where $\omega_{i}$ is the horizon radius during the $i$-th epoch.

First Epoch: Now let us look at the case where the boundary anchoring point is
in the first epoch (Refer figure 5). This case reduces to the calculation of
extremal volumes in a single-sided Schwarzschild black hole. Since we are in
the first epoch, there will only be one term in the sum (82). The surface will
always have a turning point. At late times, $r_{\text{1,min}}$ approaches a
constant radial surface in the interior of the black hole, located at the
critical point of effective potential Stanford:2014jda ; Carmi:2017jqz . We
will refer to this surface as the accumulation surface. The location of this
surface is given by

Figure 5: The left (right) figure shows the early (left) time behavior of the
extremal volume surface when the boundary anchoring points are in the first
epoch. We have not included the cutoff surface in the diagrams to avoid
cluttering. The 45-degree line depicts the horizon. The grey dashed lines in
the interior of the black hole correspond to the accumulation surfaces (Ref
(84)) of the epochs.

$\displaystyle V^{\prime}_{1}(R_{\text{1,min}})=0\implies\
R_{\text{1,min}}=\omega_{1}\left(\frac{d}{2d-2}\right)^{\frac{1}{d-2}}.$ (84)

Since the width of each epoch is much larger than the scrambling time, the
turning points will necessarily approach the accumulation surface at late
times. Now let us calculate the growth rate of (74) as a function of the
boundary anchoring time. From (79), we can see that

$\displaystyle
t+r_{\infty}^{*}-r^{*}\left(r_{1,\min}\right)=\int_{v_{1,\min}}^{v_{\infty}}dv=\int_{r_{1,\min}}^{r_{\infty}}dr\left[\frac{-E_{1}}{F_{1}(r)\sqrt{F_{1}(r)r^{2(d-1)}+E_{1}^{2}}}+\frac{1}{F_{1}(r)}\right]$

Using (3.2), we can rewrite the volume integral (82) as follows:

$\displaystyle\frac{\mathcal{V}}{\Omega_{d-1}}=\int_{r_{1,\min}}^{r_{\max}}dr\left[\frac{\sqrt{F_{1}(r)r^{2(d-1)}+E_{1}^{2}}}{F_{1}(r)}-\frac{E_{1}}{F_{1}(r)}\right]+E_{1}\left(t+r_{\infty}^{*}-r^{*}\left(r_{1,\min}\right)\right)$
(86)

Taking a derivative w.r.t $t$ and using Leibniz integral rule, we find the
simple relation

$\displaystyle\begin{aligned} \frac{d\mathcal{V}}{dt}&=\Omega_{d-1}E_{1}.\\\
&=\Omega_{d-1}\sqrt{-F_{1}(r_{1,\min})}r_{1,\min}^{d-1}\end{aligned}$ (87)

where we have used (81) to rewrite the energy in terms of the metric
components. At late times, $r_{1,\min}$ approaches $R_{1,\min}$, a constant.
Using (84), we find that

$\displaystyle\frac{d\mathcal{V}}{dt}=c_{d}\Omega_{d-1}\omega_{1}^{d-1}$ (88)

where

$\displaystyle
c_{d}=\sqrt{\frac{d-2}{d}}\left(\frac{d}{2d-2}\right)^{\frac{d-1}{d-2}}$ (89)

Therefore, the late time complexity growth during the first epoch is given by

$\displaystyle\frac{d\mathcal{C}_{1}}{dt}=\frac{1}{G_{N}\omega_{1}}\frac{d\mathcal{V}}{dt}=\frac{16\pi
c_{d}}{d-1}M_{1}$ (90)

Second epoch: Now let us look at the case where the boundary anchoring point
is in the second epoch (see figure 6).

Figure 6: The left (right) figure shows the early (left) time behavior of the
extremal volume surface when the boundary anchoring points are in the second
epoch. The extremal surface is disconnected, and the endpoints of the surface
are located using the relation (73).

Following the discussion in section 3.1, the extremal volume surfaces will be
disconnected in the Penrose diagram. The endpoints satisfy the ordering
relation (73). We will label the radial coordinate of the point where the
extremal surface intersects the epoch’s boundary by $r_{S}$. The volume
functional (82) is the sum of two terms:

$\displaystyle\mathcal{V}=\Omega_{d-1}\int_{r_{1,\min}}^{r_{S}}dr\frac{r^{2(d-1)}}{\sqrt{F_{1}(r)r^{2(d-1)}+E_{1}^{2}}}+\Omega_{d-1}\int_{r_{S}}^{r_{\infty}}dr\frac{r^{2(d-1)}}{\sqrt{F_{2}(r)r^{2(d-1)}+E_{2}^{2}}}$
(91)

The surface will always have a turning point in the first epoch. Moreover, the
turning point will be at the accumulation surface of the first epoch
$R_{\text{1,min}}$. Three possibilities arise when we consider the turning
point of the second epoch. The surface will not have a turning point during
the very early stages. However, as boundary anchoring time increases, the
surface will develop a turning point in the interior of the black hole. At
late times, the turning point $r_{2,\min}$ will approach the accumulation
surface of the epoch, given by

$\displaystyle V^{\prime}_{2}(R_{\text{2,min}})=0\implies\
R_{\text{2,min}}=\omega_{2}\left(\frac{d}{2d-2}\right)^{\frac{1}{d-2}}.$ (92)

Now let us calculate the growth rate of these volumes. As in the previous
section, we have

$\displaystyle
v_{1}-r^{*}\left(r_{1,\min}\right)=\int_{v_{1,\min}}^{v_{1}}dv=\int_{r_{1,\min}}^{r_{S}}dr\left[\frac{-E_{1}}{F_{1}(r)\sqrt{F_{1}(r)r^{2(d-1)}+E_{1}^{2}}}+\frac{1}{F_{1}(r)}\right]$
(93)

and

$\displaystyle
t+r_{\infty}^{*}-v_{1}=\int_{v_{1}}^{v_{\infty}}dv=\int_{r_{S}}^{r_{\infty}}dr\left[\frac{-E_{2}}{F_{2}(r)\sqrt{F_{2}(r)r^{2(d-1)}+E_{2}^{2}}}+\frac{1}{F_{2}(r)}\right]$
(94)

This gives us

$\displaystyle\begin{aligned}
\frac{\mathcal{V}}{\Omega_{d-1}}&=\int_{r_{1,\min}}^{r_{S}}dr\left[\frac{\sqrt{F_{1}(r)r^{2(d-1)}+E_{1}^{2}}}{F_{1}(r)}-\frac{E_{1}}{F_{1}(r)}\right]+E_{1}\left(v_{1}-r^{*}\left(r_{1,\min}\right)\right)\\\
&+\int_{r_{S}}^{r_{\max}}dr\left[\frac{\sqrt{F_{2}(r)r^{2(d-1)}+E_{2}^{2}}}{F_{2}(r)}-\frac{E_{2}}{F_{2}(r)}\right]+E_{2}\left(t+r_{\infty}^{*}-v_{1}\right)\end{aligned}$
(95)

Taking a time derivative w.r.t $t$, we get

$\displaystyle\begin{aligned}
\frac{1}{\Omega_{d-1}}\frac{d\mathcal{V}}{dt}&=E_{2}\\\
&+\frac{dr_{S}}{dt}\left[\frac{\sqrt{F_{1}(r)r^{2(d-1)}+E_{1}^{2}}}{F_{1}(r)}-\frac{E_{1}}{F_{1}(r)}-\frac{\sqrt{F_{2}(r)r^{2(d-1)}+E_{2}^{2}}}{F_{2}(r)}-\frac{E_{2}}{F_{2}(r)}\right]_{r=r_{S}}\end{aligned}$
(96)

From (79), it is easy to see that the term in the second line is proportional
to $\dot{v}_{1}(r_{S})-\dot{v}_{2}(r_{S})$. Since $v$ is continuous across
each epoch, this function vanishes, and we get the simple expression:

$\displaystyle\begin{aligned}
\frac{d\mathcal{V}}{dt}&=\Omega_{d-1}E_{2}=\Omega_{d-1}\sqrt{-F_{2}(r_{2,\min})}r_{2,\min}^{d-1}\end{aligned}$
(97)

At late times, $r_{2,\min}$ approaches the accumulation surface (92).
Substituting (92) in the above equation, we find the late time complexity
growth during the second epoch to be given by

$\displaystyle\begin{aligned}
\frac{d\mathcal{C}_{2}}{dt}&=\frac{1}{G_{N}\omega_{2}}\frac{d\mathcal{V}}{dt}\\\
&\simeq\frac{\Omega_{d-1}}{G_{N}\omega_{2}}\sqrt{-F_{2}(R_{2,\min})}R_{2,\min}^{d-1}\\\
&=\frac{\Omega_{d-1}}{G_{N}}\sqrt{\frac{d-2}{d}}\left(\frac{d}{2d-2}\right)^{\frac{d-1}{d-2}}\omega_{2}^{d-2}\\\
&=\frac{16\pi c_{d}}{d-1}M_{2},\end{aligned}$ (98)

where $c_{d}$ is defined in (89).

i-th epoch: It is easy to extend the results to the $i$-th epoch. Complexity
will undergo an initial transitional phase, after which it will settle to a
linear growth characterized by the growth rate:

$\displaystyle\frac{d\mathcal{C}_{i}}{dt}\simeq\frac{16\pi c_{d}}{d-1}M_{i}$
(99)

Now, let us look at the overall trend of complexity growth. We will do this by
taking a “continuum” limit where the width of each epoch goes to zero. Note
that we cannot formally take this limit since the scale separation between the
scrambling time and the width of an epoch was crucial for our calculation to
be valid. One should think of the limit as a means to smooth out the otherwise
jagged curve that would result from combining complexity across various
epochs. Operationally, we proceed by replacing $M(i)$ by the instantaneous
mass of the black hole $M(t)$. Let us assume that we are $3+1$-dimensions. If
the black hole is a perfect black body that satisfies the Stefan-Boltzmann
law, then we have Santi:2019lvi

$\displaystyle M(t)=M_{0}\left(1-\frac{t}{t_{E}}\right)^{1/3}$ (100)

where $t_{E}=\frac{5120\pi G_{N}^{2}M_{0}^{3}}{\hbar c^{4}}$ is the lifetime
of the black hole.

Figure 7: Holographic complexity of an evaporating black hole obtained by
taking the continuum limit of equation (99). Here $t_{E}$ is the lifetime of
the black hole.

Integrating the expression, we get

$\displaystyle\boxed{C(t)=2\sqrt{3}\pi
M_{0}t_{E}\left(1-\left(1-\frac{t}{t_{E}}\right)^{4/3}\right)}$ (101)

Plotting this function, we get figure 7.

##  4 Discussion

Let us briefly review the calculations in section 2. In conjunction with
Berry’s conjecture, the slow leakage assumption played a crucial role in
rendering the Krylov complexity calculation analytically tractable. At first
glance, the former assumption might appear limiting. However, slow leakage is
required if one wants to use intensive quantities like temperature at every
instant of the process. The appearance of an ensemble-averaged semiclassical
description in the context of black holes can be attributed to the existence
of such quasi-static equilibriums Krishnan:2021faa .

Working out details, we saw that the Krylov complexity had the following
behavior:

  * •

During an epoch, complexity goes through a scrambling phase and then
transitions to linear growth.

  * •

Complexity keeps increasing even as we cross the boundary of each epoch.
However, the late time linear growth rate decreases with each successive
epoch.

  * •

When the gas has completely leaked out of the box, complexity levels off.

We claim that these results carry over to any chaotic quasi-static open
quantum system if the operator under consideration satisfies ETH. Let us
briefly outline this calculation by examining a system interacting with its
environment. Consider an operator $O$ which satisfies the ETH ansatz. The off-
diagonal elements of this operator are given by

$\displaystyle\left\langle E_{i}\left|O\right|E_{j}\right\rangle\approx
F\left(E,\omega\right)R_{ij}$ (102)

where $E_{i,j}$ are the energy eigenstates of the system. $R_{ij}$ is a zero
mean, unit variance random matrix whereas $E$ and $\omega$ are given by (4).
During an epoch, the Krylov complexity can be calculated using the Liouvillian
of the system, following the assumptions in 2.3.1. If we assume $F$ to decay
as $\omega\to\infty$, then we can use the arguments in section 3.1 of
Barbon:2019wsy  to see that the Krylov complexity undergoes a scrambling
phase, followed by linear growth. As we have observed in the case of the
slowly leaking gas, the linear growth rate will be proportional to the degrees
of freedom of the system. Therefore, complexity will increase as we go from
one epoch to the other, and it will eventually level off if all the degrees of
freedom leak out of the system. This reproduces the advertised behavior.

The holographic complexity, computed using a gravity calculation, displays the
same behavior we described earlier. We can push this parallel further by
comparing the late time growth rates (60) and (99). From the laws of black
hole thermodynamics, we can see that the mass of the black hole plays the role
of energy Wald:1999vt . Therefore, both calculations result in the same late
time linear growth, provided we correctly identify the thermodynamic
quantities on the black hole side. This provides further evidence to the claim
that black holes can be described by a chaotic open quantum mechanical system
with finite degrees of freedom when observed from the outside Susskind:1993if
; Almheiri:2020cfm . Another manifestation of this proposal can be found in
Krishnan:2021faa , where the entanglement entropy of the slowly leaking gas
model matched the gravitational path integral calculation in Penington:2019kki
. We can put these statements on a firmer footing by thinking of it as a
consequence of black hole complementarity Susskind:1993if , which posits that
the interior of a black hole can be thought to be described by a finite number
of quantum mechanical degrees of freedom living on the stretched horizon of
the black hole. The hard spheres in the left box assume the role of these
degrees of freedom, while the particles in the right box model the outgoing
Hawking radiation, allowing us to make a one-to-one map between an evaporating
black hole and our two-box system.

Returning to the growth rate (99), we find that

$\displaystyle\frac{dC}{dt}\sim M\propto ST$ (103)

where $S$ and $T$ are the entropy and temperature of the black hole during
that epoch. Our results are in tandem with Susskind:2018pmk  and the 2d
gravity calculation performed in Schneiderbauer:2019anh .

One can also calculate complexity using the Complexity=Action (CA)
prescription Brown:2015bva ; Brown:2015lvg ; Schneiderbauer:2020isp . However,
to obtain meaningful results, the inclusion of a counterterm for null
boundaries, similar to the one used in Chapman:2018dem ; Chapman:2018lsv ;
Akhavan:2019zax ; Omidi:2020oit , might be necessary.

###### Acknowledgements.

I thank Friðrik Freyr Gautason, Chethan Krishnan, Watse Sybesma, and Lárus
Thorlacius for their insights and helpful discussions. This work was supported
by the Icelandic Research Fund grant 228952-052.

## References

  * (1) D.E. Parker, X. Cao, A. Avdoshkin, T. Scaffidi and E. Altman, _A Universal Operator Growth Hypothesis_ , _Phys. Rev. X_ 9 (2019) 041017 [1812.08657].
  * (2) E. Rabinovici, A. Sánchez-Garrido, R. Shir and J. Sonner, _Krylov complexity from integrability to chaos_ , _JHEP_ 07 (2022) 151 [2207.07701].
  * (3) B.L. Español and D.A. Wisniacki, _Assessing the saturation of Krylov complexity as a measure of chaos_ , _Phys. Rev. E_ 107 (2023) 024217 [2212.06619].
  * (4) K. Hashimoto, K. Murata, N. Tanahashi and R. Watanabe, _Krylov complexity and chaos in quantum mechanics_ , 2305.16669.
  * (5) H.A. Camargo, V. Jahnke, H.-S. Jeong, K.-Y. Kim and M. Nishida, _Spectral and Krylov Complexity in Billiard Systems_ , 2306.11632.
  * (6) L. Susskind, _Entanglement is not enough_ , _Fortsch. Phys._ 64 (2016) 49 [1411.0690].
  * (7) L. Susskind, _Computational Complexity and Black Hole Horizons_ , _Fortsch. Phys._ 64 (2016) 24 [1403.5695].
  * (8) D. Stanford and L. Susskind, _Complexity and Shock Wave Geometries_ , _Phys. Rev. D_ 90 (2014) 126007 [1406.2678].
  * (9) J.L.F. Barbón, E. Rabinovici, R. Shir and R. Sinha, _On The Evolution Of Operator Complexity Beyond Scrambling_ , _JHEP_ 10 (2019) 264 [1907.05393].
  * (10) A. Dymarsky and M. Smolkin, _Krylov complexity in conformal field theory_ , _Phys. Rev. D_ 104 (2021) L081702 [2104.09514].
  * (11) A. Avdoshkin, A. Dymarsky and M. Smolkin, _Krylov complexity in quantum field theory, and beyond_ , 2212.14429.
  * (12) H.A. Camargo, V. Jahnke, K.-Y. Kim and M. Nishida, _Krylov complexity in free and interacting scalar field theories with bounded power spectrum_ , _JHEP_ 05 (2023) 226 [2212.14702].
  * (13) P. Caputa and S. Datta, _Operator growth in 2d CFT_ , _JHEP_ 12 (2021) 188 [2110.10519].
  * (14) K. Adhikari, S. Choudhury and A. Roy, _Krylov Complexity in Quantum Field Theory_ , _Nucl. Phys. B_ 993 (2023) 116263 [2204.02250].
  * (15) A. Kundu, V. Malvimat and R. Sinha, _State Dependence of Krylov Complexity in $2d$ CFTs_, 2303.03426.
  * (16) D. Patramanis and W. Sybesma, _Krylov complexity in a natural basis for the Schrödinger algebra_ , 2306.03133.
  * (17) K. Adhikari and S. Choudhury, _Cosmological Krylov Complexity_ , _Fortsch. Phys._ 70 (2022) 2200126 [2203.14330].
  * (18) B. Bhattacharjee, X. Cao, P. Nandy and T. Pathak, _Krylov complexity in saddle-dominated scrambling_ , _JHEP_ 05 (2022) 174 [2203.03534].
  * (19) M.J. Vasli, K. Babaei Velni, M.R. Mohammadi Mozaffar, A. Mollabashi and M. Alishahiha, _Krylov Complexity in Lifshitz-type Scalar Field Theories_ , 2307.08307.
  * (20) E. Rabinovici, A. Sánchez-Garrido, R. Shir and J. Sonner, _A bulk manifestation of Krylov complexity_ , _JHEP_ 08 (2023) 213 [2305.04355].
  * (21) S.-K. Jian, B. Swingle and Z.-Y. Xian, _Complexity growth of operators in the SYK model and in JT gravity_ , _JHEP_ 03 (2021) 014 [2008.12274].
  * (22) J. Erdmenger, S.-K. Jian and Z.-Y. Xian, _Universal chaotic dynamics from Krylov space_ , _JHEP_ 08 (2023) 176 [2303.12151].
  * (23) N. Hörnedal, N. Carabba, A.S. Matsoukas-Roubeas and A. del Campo, _Ultimate Speed Limits to the Growth of Operator Complexity_ , _Commun. Phys._ 5 (2022) 207 [2202.05006].
  * (24) C. Krishnan and V. Mohan, _Hints of gravitational ergodicity: Berry’s ensemble and the universality of the semi-classical Page curve_ , _JHEP_ 05 (2021) 126 [2102.07703].
  * (25) M.V. Berry, _Semiclassical mechanics of regular and irregular motion_ , _Les Houches lecture series_ 36 (1983) 171.
  * (26) M. Berry, _Chaos and quantum physics_ , _Les Houches LII/Ed. M. Giannoni, A. Voros, J. Zinn-Justin. Amsterdam: North-Holland_ (1991) .
  * (27) M. Srednicki, _Chaos and Quantum Thermalization_ , cond-mat/9403051.
  * (28) J.M. Deutsch, _Quantum statistical mechanics in a closed system_ , _Phys. Rev. A_ 43 (1991) 2046.
  * (29) A. Bhattacharya, P. Nandy, P.P. Nath and H. Sahu, _Operator growth and Krylov construction in dissipative open quantum systems_ , _JHEP_ 12 (2022) 081 [2207.05347].
  * (30) C. Liu, H. Tang and H. Zhai, _Krylov Complexity in Open Quantum Systems_ , 2207.13603.
  * (31) B. Bhattacharjee, X. Cao, P. Nandy and T. Pathak, _Operator growth in open quantum systems: lessons from the dissipative SYK_ , _JHEP_ 03 (2023) 054 [2212.06180].
  * (32) A. Bhattacharya, P. Nandy, P.P. Nath and H. Sahu, _On Krylov complexity in open systems: an approach via bi-Lanczos algorithm_ , 2303.04175.
  * (33) P. Caputa, J.M. Magan and D. Patramanis, _Geometry of Krylov complexity_ , _Phys. Rev. Res._ 4 (2022) 013041 [2109.03824].
  * (34) A. Kar, L. Lamprou, M. Rozali and J. Sully, _Random matrix theory for complexity growth and black hole interiors_ , _JHEP_ 01 (2022) 016 [2106.02046].
  * (35) Z.-Y. Fan, _Universal relation for operator complexity_ , _Phys. Rev. A_ 105 (2022) 062210 [2202.07220].
  * (36) J.D. Noh, _Operator growth in the transverse-field ising spin chain with integrability-breaking longitudinal field_ , _Physical Review E_ 104 (2021) 034112.
  * (37) D.V. Schroeder, _An introduction to thermal physics_ , Addison Wesley, us ed ed. (1999).
  * (38) S. Chapman, H. Marrochio and R.C. Myers, _Holographic complexity in Vaidya spacetimes. Part I_ , _JHEP_ 06 (2018) 046 [1804.07410].
  * (39) S. Chapman, H. Marrochio and R.C. Myers, _Holographic complexity in Vaidya spacetimes. Part II_ , _JHEP_ 06 (2018) 114 [1805.07262].
  * (40) C.W. Misner, K.S. Thorne and J.A. Wheeler, _Gravitation_ , W. H. Freeman, San Francisco (1973).
  * (41) D. Carmi, S. Chapman, H. Marrochio, R.C. Myers and S. Sugishita, _On the Time Dependence of Holographic Complexity_ , _JHEP_ 11 (2017) 188 [1709.10184].
  * (42) N.S.M. Santi and R. Santarelli, _Mass evolution of Schwarzschild black holes_ , _Braz. J. Phys._ 49 (2019) 897 [1906.07088].
  * (43) R.M. Wald, _The thermodynamics of black holes_ , _Living Rev. Rel._ 4 (2001) 6 [gr-qc/9912119].
  * (44) L. Susskind, L. Thorlacius and J. Uglum, _The Stretched horizon and black hole complementarity_ , _Phys. Rev. D_ 48 (1993) 3743 [hep-th/9306069].
  * (45) A. Almheiri, T. Hartman, J. Maldacena, E. Shaghoulian and A. Tajdini, _The entropy of Hawking radiation_ , _Rev. Mod. Phys._ 93 (2021) 035002 [2006.06872].
  * (46) G. Penington, S.H. Shenker, D. Stanford and Z. Yang, _Replica wormholes and the black hole interior_ , _JHEP_ 03 (2022) 205 [1911.11977].
  * (47) L. Susskind, _Three Lectures on Complexity and Black Holes_ , SpringerBriefs in Physics, Springer, 10, 2018, DOI [1810.11563].
  * (48) L. Schneiderbauer, W. Sybesma and L. Thorlacius, _Holographic Complexity: Stretching the Horizon of an Evaporating Black Hole_ , _JHEP_ 03 (2020) 069 [1911.06800].
  * (49) A.R. Brown, D.A. Roberts, L. Susskind, B. Swingle and Y. Zhao, _Holographic Complexity Equals Bulk Action?_ , _Phys. Rev. Lett._ 116 (2016) 191301 [1509.07876].
  * (50) A.R. Brown, D.A. Roberts, L. Susskind, B. Swingle and Y. Zhao, _Complexity, action, and black holes_ , _Phys. Rev. D_ 93 (2016) 086006 [1512.04993].
  * (51) L. Schneiderbauer, W. Sybesma and L. Thorlacius, _Action Complexity for Semi-Classical Black Holes_ , _JHEP_ 07 (2020) 173 [2001.06453].
  * (52) A. Akhavan and F. Omidi, _On the Role of Counterterms in Holographic Complexity_ , _JHEP_ 11 (2019) 054 [1906.09561].
  * (53) F. Omidi, _Regularizations of Action-Complexity for a Pure BTZ Black Hole Microstate_ , _JHEP_ 07 (2020) 020 [2004.11628].