aainstitutetext: The Laboratory for Quantum Gravity & Strings, Department of
Mathematics & Applied Mathematics, University of Cape Town, Cape Town, South
Africabbinstitutetext: The National Institute for Theoretical and
Computational Sciences, Private Bag X1, Matieland, South Africa

# Krylov Complexity for Jacobi Coherent States

S. Shajidul Haque a,b Jeff Murugan a Mpho Tladi a,b and Hendrik J.R. Van Zyl

###### Abstract

We develop computational tools necessary to extend the application of Krylov
complexity beyond the simple Hamiltonian systems considered thus far in the
literature. As a first step toward this broader goal, we show how the Lanczos
algorithm that iteratively generates the Krylov basis can be augmented to
treat coherent states associated with the Jacobi group, the semi-direct
product of the 3-dimensional real Heisenberg-Weyl group $H_{1}$, and the
symplectic group, $Sp(2,\mathbb{R})\simeq SU(1,1)$. Such coherent states are
physically realized as squeezed states in, for example, quantum optics. With
the Krylov basis for both the $SU(1,1)$ and Heisenberg-Weyl groups being well
understood, their semi-direct product is also partially analytically
tractable. We exploit this to benchmark a scheme to numerically compute the
Lanczos coefficients which, in principle, generalizes to the more general
Jacobi group $H_{n}\rtimes Sp(2n,\mathbb{R})$.

##  1 Introduction

Quantum complexity nielsen00 ; Watrous_2008  aims to quantify how difficult it
is to prepare a desired final state of a quantum system, given some initial
configuration. This requires information about the initial state, the set of
unitary operations at one’s disposal and a measure of the difficulty in
performing each of the unitary operations. As such, the definition of quantum
complexity is not unique and there are, by now, several candidates that are
the subject of much active research (see e.g. reviews Chapman:2021jbh ;
Chen:2021lnq  and references therein).

One such measure, the so-called Krylov- or spread-complexity Parker:2018yvk ,
takes as a starting point the time-evolution of a “simple” reference state in
the Hilbert space. This is intuitively appealing since the time-evolution will
naturally involve the spectrum of the system under consideration which encodes
many of the properties of the system in an invariant way. For example, in
chaotic systems, the spectral statistics of the Hamiltonian are conjectured to
asymptote to those of a random matrix ensemble at late times Bohigas:1983er  .
As such, it is reasonable to expect that Krylov complexity may be used as a
diagnostic of quantum chaotic systems, to supplement other diagnostics
Ali:2019zcj ; Bhattacharyya:2019txx ; Bhattacharyya:2020art ;
Bhattacharyya:2020iic . Furthermore, the space of states explored by
integrable and chaotic systems differ dramatically Balasubramanian:2021mxo ,
with the latter expected to explore a Hilbert space subspace in an ergodic
way. To quantify the complexity of the time-evolved state one may imagine
decomposing it into a linear combination of states from some basis, where the
basis vectors are ordered in increasing complexity. The Krylov basis minimises
the complexity of the time-evolved state over these choices of basis
Caputa:2021sib ; Balasubramanian:2022tpr .

Spread complexity has already been shown to be sensitive to topological phase
transitions Caputa:2022eye ; Caputa:2022yju  as exemplified in the Su-
Shriefer-Heeger (SSH) model 111Previously, the circuit complexity for this
model was studied in Ali:2018aon . of polyacetylene Su:1979ua  polymer chains
as well as the 1-dimensional Kitaev chain Kitaev_2001  p-wave superconductor,
a feature that is robust against a particular choice of reference state
Caputa:2022yju . Key to the computation of the spread complexity in these
quantum systems was the fact that in each case, the Hamiltonian can be
expressed in terms of $SU(2)$ coherent states, for which exact analytic
results relating to spread complexity are known222 See also Muck:2022xfc
where analytic results are obtained from an input of simple Lanczos
coefficients. Caputa:2021sib . Still, as remarkable as both these models are,
they are a set of measure zero in the space of novel quantum systems. Ideally,
we would like to understand the nature and behaviour of Krylov complexity in
more general systems, with the Hamiltonian selected from a larger symmetry
group, or even where perturbative treatments may be unavoidable. This article
aims to initiate development in just this 333 Note that modifications of the
Lanczos algorithm has been developed in the context of diagonalizing sparse
matrices (see e.g. Parlett1979TheLA ; Simon1984TheLA ; Gambolati ). Here we
are interested in a modification aimed at the computations of spread
complexity of quantum states. . Specifically, we will focus on the Jacobi
algebra, a semi-direct sum of $su(1,1)$ and the Heisenberg algebra
Berceanu:2007zz .

Krylov complexity for the Jacobi group is a conservative, yet interesting
extension of the $SU(2)$ case. This is in equal measure because the coherent
states attached to the Jacobi group have an actual physical realization as
squeezed states in quantum optics Berceanu2009TheJG  and, more importantly,
since exact analytic expressions may be obtained for both $SU(1,1)$ and
Heisenberg Krylov complexity Caputa:2021sib , there is hope that similar
expressions for their semi-direct product exist also. Indeed tensor products
of exactly solvable algebras have also been studied Kressner  and, in fact, do
yield analytic expressions (see e.g. Adhikari2022CosmologicalKC ;
Caputa:2022eye ; Caputa:2022yju ). Crucially however, in the case of the
Jacobi algebra, the non-vanishing commutation relations between the elements
of the $su(1,1)$ and Heisenberg algebras give rise to involved combinatorics
in determining the Krylov basis. In this sense, the Jacobi algebra is a simple
prototype for more general Lie algebras.

Conventionally, the Krylov basis is generated by the Lanczos algortihm
Lanczos1950AnIM  (see also the review viswanath2008recursion ), an iterative
algorithm that generalizes the Gram-Schmidt procedure for constructing an
orthonormal basis. However, a naive application of the Lanczos algorithm to
the Jacobi group fails to efficiently compute the Lanczos coefficients and
must be modified. The purpose of this article is to propose such a
modification. Our proposed algorithm generates a family of vectors each step
in the iterative scheme, converging on the Krylov basis and in which at the
$k{th}$ step contains the first $k$ Krylov vectors. Explicit computations in
this modified recursion method may be performed using differential operators
acting on Jacobi coherent states Perelomov . Moreover, by terminating the
procedure at any finite order one obtains an approximate Krylov basis. We
demonstrate that this approximation scheme is accurate for a range of
parameter choices.

This article is structured as follows: In section 2 we review the idea of
Krylov complexity in the context of symmetry groups, paying particular
attention to the cases of $SU(1,1)$ and the Heisenberg group $H_{1}$. This is
followed by a detailed discussion of our proposed modification of the Lanczos
algorithm in section 3 and an explicit computation of the Lanczos coefficients
for the Jacobi group as well as the associate K-complexity in section 4. We
conclude with some discussion and speculation of future directions in section
5.

##  2 Krylov Complexity

We begin with a brief overview of idea of K-complexity for quantum states. For
the computation of state $K$-complexity we require a reference state (chosen
from a Hilbert space with an appropriate inner product) as well as a
Hamiltonian acting on this reference state. The starting point for the
analysis is the time-evolved reference state

$|t\rangle\equiv e^{itH}|\phi_{0}\rangle\,,$ (1)

By expanding this in $t$, one obtains a basis for the space of accessible
target states

$|O_{n})\equiv H^{n}|\phi_{0}\rangle\,,$ (2)

However, this basis is not, in general, orthogonal with respect to to the
Hilbert space inner product. The Krylov basis is obtained by applying a Gram-
Schmidt orthogonalisation procedure on the basis $\\{|O_{n})\\}$. One popular
method by which to achieve this is the Lanczos algorithm, or recursion method
by which we define,

$\displaystyle|K_{0}\rangle$ $\displaystyle\equiv$
$\displaystyle|\phi_{0}\rangle\,,$ $\displaystyle|A_{n+1})$ $\displaystyle=$
$\displaystyle(H-a_{n})|K_{n}\rangle-b_{n}|K_{n-1}\rangle\,,$
$\displaystyle|K_{n}\rangle$ $\displaystyle\equiv$
$\displaystyle(A_{n}|A_{n})^{-\frac{1}{2}}|A_{n})\,,$

and where

$\displaystyle a_{n}$ $\displaystyle=$ $\displaystyle\langle
K_{n}|H|K_{n}\rangle\,,$ $\displaystyle b_{n}$ $\displaystyle=$
$\displaystyle(A_{n}|A_{n})^{\frac{1}{2}}\,.$ (3)

A notable feature of this construction is that it tri-diagonalizes the
Hamiltonian, in the sense that in the Krylov basis,

$H|K_{n}\rangle=a_{n}|K_{n}\rangle+b_{n}|K_{n-1}\rangle+b_{n+1}|K_{n+1}\rangle\,.$
(4)

There are two ways of interpreting this: i) knowledge of the Krylov basis
determines the Lanczos coefficients as Hamiltonian matrix elements and ii) the
Lanczos coefficients give the Krylov basis by means of the recursive Lanczos
algorithm.

Finally, given the Krylov basis, the spread complexity (or state K-complexity)
of a target state $|\psi\rangle$ is given by

$K=\sum_{n=0}^{\infty}n\langle\psi|K_{n}\rangle\langle K_{n}|\psi\rangle\,.$
(5)

This algorithm may be applied, in principle, to any Hamiltonian, however when
the Hamiltonian itself is an element of some symmetry algebra, the computation
simplifies considerably, in large part due to the utility of coherent states
attached to the symmetry Perelomov .

###  2.1 $SU(1,1)$ and Heisenberg K-complexity

Our primary interest in this article is in the study of Krylov complexity when
the Hamiltonian is expressed in terms of elements of the Jacobi algebra.
Physically, this is the most general choice of Hamiltonian up to quadratic
order in the position and momentum operators. As the semi-direct sum of
$su(1,1)$ and the Heisenberg algebra, the Jacobi algebra is generated by the
operators $L_{\pm},L_{0},l_{\pm}$ and $l_{0}$ satisfying the commutation
relations

$\displaystyle\left[L_{0},L_{\pm}\right]=\pm L_{\pm}\ \ \ $ $\displaystyle\ \
\ \left[L_{-},L_{+}\right]=2L_{0}$
$\displaystyle\left[l_{-},l_{+}\right]=l_{0}\ \ \ $ $\displaystyle\ \ \
\left[l_{0},l_{\pm}\right]=0$ (6)

respectively, as well as the “cross-commutators”

$\displaystyle\left[L_{0},l_{\pm}\right]$ $\displaystyle=$
$\displaystyle\pm\frac{1}{2}l_{\pm}$ $\displaystyle\left[L_{+},l_{-}\right]$
$\displaystyle=$ $\displaystyle-l_{+}$ (7)
$\displaystyle\left[L_{-},l_{+}\right]$ $\displaystyle=$ $\displaystyle l_{-}$

with all other commutators vanishing. It is this final set of commutators, of
course, that enlarges the Jacobi algebra from the direct sum of $su(1,1)$ and
Heisenberg algebra.

Since exact analytical results for the Lanczos coefficients and Krylov
complexity for Hamiltonians selected from $su(1,1)$ or the Heisenberg algebra
may be readily obtained (see for example Caputa:2021sib  and references
therein), one might be forgiven for expecting that the Krylov basis (with say,
the lowest weight reference state) for a Jacobi algebra Hamiltonian may also
be computed exactly, and in closed form. However, as we will demonstrate, the
additional commutators in (2.1) result in rather involved combinatorics
associated with the Gram-Schmidt procedure generating the Krylov basis.

Before proceeding to the specifics of the Jacobi K-complexity, let’s first
summarise the simpler cases of $SU(1,1)$ and Heisenberg K-complexity in a form
that will allow for necessary modifications. To this end, it will be
particularly useful to phrase the computation in the language of coherent
states and differential operators, which will generalise nicely when we
consider the Jacobi group K-complexity.

####  2.1.1 $SU(1,1)$

Let’s begin with the more conceptually familiar $su(1,1)$ case. The most
general $su(1,1)$ Hamiltonian can be written as

$H=\alpha L_{+}+\alpha^{*}L_{-}+\gamma L_{0}\,,$ (8)

with (complex) coefficients $\alpha$ and $\gamma$. As reference state we
choose the lowest weight, $|h\rangle$, defined by

$\displaystyle L_{0}|h\rangle=h|h\rangle\,,\quad L_{-}|h\rangle=0\,.$ (9)

The space of states spanned by $\left\\{H^{m}|h\rangle\ ,\ m=0,1,\cdots
n\right\\}$ is then equivalent the space spanned by
$\left\\{L_{+}^{m}|h\rangle\ ,\ m=0,1,\cdots n\right\\}$ for each value of
$n$. Indeed, provided $\alpha\neq 0$, the Krylov basis is given by

$|K_{n}\rangle=\left(\langle
h|L_{-}^{n}L_{+}^{n}|h\rangle\right)^{-\frac{1}{2}}L_{+}^{n}|h\rangle\,.$ (10)

Note that this is independent of the coefficients parametrising the
Hamiltonian. Consequently, the (unnormalised) $SU(1,1)$ Perelemov coherent
states

$|z)=e^{zL_{+}}|h\rangle=\sum_{n=0}^{\infty}\frac{z^{n}}{n!}L_{+}^{n}|h\rangle\,,$
(11)

provide a natural basis in which quantities associated with Krylov complexity
may be computed. For what follows, a particularly convenient representation
for the $su(1,1)$ generators is as differential operators, which act on the
associated coherent states as

$\displaystyle L_{+}|z)$ $\displaystyle=$ $\displaystyle\partial_{z}|z)\,,$
$\displaystyle L_{0}|z)$ $\displaystyle=$
$\displaystyle(h+z\partial_{z})|z)\,,$ (12) $\displaystyle L_{-}|z)$
$\displaystyle=$ $\displaystyle(2zh+z^{2}\partial_{z})|z)\,.$

Two further ingredients for the Lanczos algorithm are required; the overlap
between coherent states,

$(\bar{z}|z)=\sum_{n}\frac{z^{n}\bar{z}^{n}}{(n!)^{2}}\langle
h|L_{-}^{n}L_{+}^{n}|h\rangle=(1-z\bar{z})^{-2h}$ (13)

and the coherent-state expectation value of the Hamiltonian,

$\displaystyle(\bar{z}|H|z)$ $\displaystyle=$ $\displaystyle
h\left(2z\alpha^{*}+2\bar{z}\alpha+(1+z\bar{z})\gamma\right)(1-z\bar{z})^{-(2h+1)}$
(14) $\displaystyle=$
$\displaystyle\sum_{n=1}^{\infty}\frac{\bar{z}^{n}z^{n-1}}{(n-1)!n!}\sqrt{\langle
h|L_{-}^{n-1}L_{+}^{n-1}|h\rangle\langle h|L_{-}^{n}L_{+}^{n}|h\rangle}\
b_{n}$ $\displaystyle+$
$\displaystyle\sum_{n=0}^{\infty}\frac{\bar{z}^{n}z^{n+1}}{(n+1)!n!}\sqrt{\langle
h|L_{-}^{n+1}L_{+}^{n+1}|h\rangle\langle h|L_{-}^{n}L_{+}^{n}|h\rangle}\
b_{n+1}$ $\displaystyle+$
$\displaystyle\sum_{n=0}^{\infty}\frac{\bar{z}^{n}z^{n}}{(n!)^{2}}\langle
h|L_{-}^{n}L_{+}^{n}|h\rangle a_{n}\,.$

Putting this together allows us to extract the Lanczos $a_{n}$ and $b_{n}$
coefficients as

$\displaystyle b_{n}$ $\displaystyle=$
$\displaystyle\left.\hat{O}_{n}^{\dagger}\hat{O}_{n-1}\left(\bar{z}|H|z\right)\right|_{z,\bar{z}\rightarrow
0}\,,$ $\displaystyle a_{n}$ $\displaystyle=$
$\displaystyle\left.\hat{O}_{n}^{\dagger}\hat{O}_{n}(\bar{z}|H|z)\right|_{z,\bar{z}\rightarrow
0}\,,$ (15)

where the family of differential operators,

$\hat{O}_{n}=\frac{1}{\sqrt{\langle
h|L_{-}^{n}L_{+}^{n}|h\rangle}}\partial_{z}^{n}\,.$ (16)

####  2.1.2 Heisenberg

Moving on, similar results may be obtained for the Heisenberg algebra. Here
the most general hermitian Hamiltonian operator we may write down (dropping
the contribution proportional to the identity) is

$H=\beta l_{+}+\beta^{*}l_{-}\,.$ (17)

Again, a natural choice of reference state is the lowest weight state,
$|0\rangle$ satisfying

$l_{-}|0\rangle=0\,.$ (18)

The space of states spanned by $\left\\{H^{m}|0\rangle\ ,\ m=0,1,\cdots
n\right\\}$ is equivalent to the space spanned by $\left\\{l_{+}^{m}|0\rangle\
,\ m=0,1,\cdots n\right\\}$. For any non-zero value of $\beta$, the recursion
method described above yeilds the Krylov basis

$|K_{n}\rangle=\left(\langle
0|l_{-}^{n}l_{+}^{n}|0\rangle\right)^{-\frac{1}{2}}l_{+}^{n}|0\rangle\,.$ (19)

Coherent states in this case are the so-called Glauber coherent states
Glauber1963CoherentAI ,

$|z)=e^{zl_{+}}|0\rangle\,.$ (20)

As before, in order to extract the Lanczos coefficients, we will require the
expectation value of the Hamiltonian in the coherent state basis. This is
simply computed as,

$(\bar{z}|H|z)=(\alpha\bar{z}+z\alpha^{*})e^{z\bar{z}}\,,$ (21)

and along with the family of differential operators

$\hat{O}_{n}=\frac{1}{\sqrt{\langle
0|l_{-}^{n}l_{+}^{n}|0\rangle}}\partial_{z}^{n}\,,$ (22)

provides a simple and exact expression for the Lanczos coefficients

$b_{n}=\hat{O}_{n}^{\dagger}\hat{O}_{n-1}(\bar{z}|H|z)\,.$ (23)

Given the simplicity of these expressions, as well as the (relative) ease with
which they can be derived, one may be justified in hoping that similar
expressions can be derived for Hamiltonian systems that combine these two
algebras in some way.

###  2.2 Jacobi group

We now turn our attention to the Jacobi group and describe the details of our
preferredrepresentation. As a first observation, note that the operators
$L_{0}$, $l_{0}$ and the Casimir
$C=(L_{0})^{2}-\frac{1}{2}(L_{+}L_{-}+L_{-}L_{+})$ may be simultaneously
diagonalized. The simultaneous eigenstate has eigenvalues

$\displaystyle L_{0}|h,n;M\rangle$ $\displaystyle=$
$\displaystyle(h+n)|h,n;M\rangle\,,$ (24) $\displaystyle l_{0}|h,n;M\rangle$
$\displaystyle=$ $\displaystyle M|h,n;M\rangle\,.$

The ladder operators act on this state as

$\displaystyle L_{+}|h,n;M\rangle$ $\displaystyle=$
$\displaystyle\sqrt{(n+1)(2h+n)}|h,n+1;M\rangle\,,$ $\displaystyle
l_{+}|h,n;M\rangle$ $\displaystyle=$
$\displaystyle\sqrt{M}\sqrt{2n+1}|h,n+\frac{1}{2};M\rangle\,,$ $\displaystyle
L_{-}|h,n;M\rangle$ $\displaystyle=$
$\displaystyle\sqrt{n(2h+n-1)}|h,n-1;M\rangle\,,$ $\displaystyle
l_{-}|h,n;M\rangle$ $\displaystyle=$
$\displaystyle\sqrt{M}\sqrt{2n}|h,n-\frac{1}{2};M\rangle\,.$ (25)

In addition, there is one more consistency condition that needs to be imposed.
It isbest seen by studying the overlap

$\langle h,0;M|L_{-}l_{+}l_{+}|h,0;M\rangle=\sqrt{4h}M=M\ \ \ \Rightarrow\ \ \
h=\frac{1}{4}\ \ \ \textnormal{or}\ \ \ M=0\,.$ (26)

We will be interested in the case $M\neq 0$ and will therefore take
$h=\frac{1}{4}$.

The Perelemov coherent states Perelomov  is given by

$|z_{1},z_{2})=e^{z_{1}l_{+}}e^{z_{2}L_{+}}|h,0,M\rangle$ (27)

with the normalised vectors given by
$|z_{1},z_{2}\rangle=(\bar{z}_{1},\bar{z}_{2}|z_{1},z_{2})^{-\frac{1}{2}}|z_{1},z_{2})$.
A useful quantity for our purposes is the coherent state overlap, given by

$\displaystyle\langle
h,0;M|e^{\bar{z}_{2}L_{-}}e^{\bar{z}_{1}l_{-}}e^{z_{1}l_{+}}e^{z_{2}L_{+}}|h,0;M\rangle$
(28) $\displaystyle=$ $\displaystyle\langle
h,0;M|e^{\frac{z_{2}}{1-z_{2}\bar{z}_{2}}L_{+}}e^{\frac{z_{1}+\bar{z}_{1}z_{2}}{1-z_{2}\bar{z}_{2}}l_{+}}e^{-2\log\left(1-z_{1}\bar{z}_{1}\right)L_{0}}e^{\frac{1}{2}\frac{2z_{1}\bar{z}_{1}+\bar{z}_{1}^{2}z_{2}+(z_{1})^{2}\bar{z}_{2}}{1-z_{2}\bar{z}_{2}}l_{0}}e^{\frac{\bar{z}_{1}+z_{1}\bar{z}_{2}}{1-z_{2}\bar{z}_{2}}l_{-}}e^{\frac{\bar{z}_{2}}{1-z_{2}\bar{z}_{2}}L_{-}}|h,0;M\rangle$
$\displaystyle=$
$\displaystyle(1-z_{2}\bar{z}_{2})^{-2h}e^{\frac{M}{2}\frac{2z_{1}\bar{z}_{1}+\bar{z}_{1}^{2}z_{2}+(z_{1})^{2}\bar{z}_{2}}{1-z_{2}\bar{z}_{2}}}\,.$

In arriving at the final expression in (28) we have made extensive use of the
BCH formula. The operator actions (25) may be verified using the overlap by
taking appropriate derivatives. Finally, note that the factor of $M$ may be
absorbed by rescaling $z_{1}\rightarrow\frac{z_{1}}{\sqrt{M}}$. As such, and
especially when obtaining numeric results, we will set $M=1$.

##  3 Jacobi Group Krylov basis

For both the $SU(1,1)$ and Heisenberg symmetry groups a key feature allowing
for exact solutions is that the span of the vectors
$\left\\{H^{m}|\phi_{0}\rangle,\ m=0,1,\cdots n\right\\}$ matches up exactly
with the span of $\left\\{L_{+}^{m}|\phi_{0}\rangle\ ,\ m=0,1,\cdots
n\right\\}$ and $\left\\{l_{+}^{m}|\phi_{0}\rangle\ ,\ m=0,1,\cdots
n\right\\}$ respectively. This plays a central role in readily identifying the
Krylov basis. The Jacobi group however, possesses two sets of ladder
operators, one of which increases the eigenvalue of $L_{0}$ by $1$ while the
other by $\frac{1}{2}$, as in (25). After $n$ applications of a general Jacobi
Hamiltonian one explores a space consisting of specific linear combinations of
the first $n$ half-integer labelled states and the first $(n+1)$ integer
labelled states. Because of this we cannot extract the Krylov basis as easily
as in the $SU(1,1)$ or Heisenberg cases.

As is, by now, well appreciated, the Lanczos algorithm provides an algorithmic
way to generate the Krylov basis vectors. Less well-known though is the fact
that when studying problems involving an infinite-dimensional Krylov subspace
(especially at late times) the algorithm is far less efficient since it
generates the Krylov basis vectors one at a time and a very large number of
basis vectors are required in order to accurately describe late-time behavior.

It is possible to circumvent this problem by modifying the Lanczos algorithm
with the aid of some intuition for the Krylov basis vectors. In the next
subsections we will illustrate this modified iterative method, but first let’s
point out a few features worth highlighting: After $m$ steps in our iteration
the first $m$ Krylov basis vectors are obtained exactly, precisely as one
would in the Lanczos algorithm. Unlike the Lanczos algorithm however, one
obtains not a single vector but a basis of vectors simultaneously. By
terminating the scheme at finite order we obtain an approximate Krylov basis.
The Hamiltonian, when expressed in this approximate basis, is tri-diagonal so
that the scheme provides a systematic way to approximate the Lanczos
coefficients. Alhough our focus in this paper is a Hamiltonian selected from
the Jacobi algebra the scheme is not dependent on this fact and should
generalise.

###  3.1 An iterative scheme

We now turn our attention to an explicit Hamiltonian built from generators of
the Jacobi algebra,

$H=\alpha L_{+}+\alpha^{*}L_{-}+\beta l_{+}+\beta^{*}l_{-}+\gamma L_{0}$ (29)

We will also occasionally refer to the phases of these complex parameters,
using $\alpha=|\alpha|e^{ia}$, $\beta=|\beta|e^{ib}$. As reference state we
choose $|h,0,M\rangle$ which provides a large stability subgroup444Note,
however, that we can easily extend the results we obtain to any choice of
reference state $U|h,0,M\rangle$ with $U\in SU(1,1)\rtimes H_{1}$. Such a
choice of reference state is equivalent to studying the $|h,0,M\rangle$
reference state with $UHU^{\dagger}$ as Hamiltonian. This unitary
transformation may be absorbed as a coordinate transformation of the
constants, $\alpha,\beta$ and $\gamma$. In appendix A we highlight a special
choice of reference state for which the system is exactly solvable., since it
is a lowest weight state of $L_{0}$, annihilated by both $L_{-}$ and $l_{-}$.
Our intuition is guided by the up-ladder part of the Hamiltonian which, after
$n$ applications on the reference state reads

$(\alpha L_{+}+\beta l_{+})^{n}|h,0,M\rangle$ (30)

While these vectors are not orthogonal to the space spanned by
$\left\\{H^{m}|h,0,M\rangle\ ,\ m=1,2,\cdots,n-1\right\\}$, it is
straightforward to identity a component that is; namely

$\left|A_{n}^{(0)}\right\rangle=\frac{\left(\alpha L_{+}^{n}+n\beta
L_{+}^{n-1}l_{+}\right)|h,0,M\rangle}{\sqrt{|\alpha|^{2}\left\langle
L_{-}^{n}L_{+}^{n}\right\rangle+n^{2}|\beta|^{2}\left\langle
L_{-}^{n-1}l_{-}l_{+}L_{+}^{n-1}\right\rangle}}$ (31)

This observation is the starting point for our approximation scheme and
represents the “$0^{\mathrm{t}h}$” order approximation for the Krylov basis.
This family of vectors is what one may have guessed as the Krylov basis based
on the intuition for the $su(1,1)$ and Heisenberg cases. However, these
vectors represent a particular linear combination of $|h,n,M\rangle$ and
$|h,n-\frac{1}{2},M\rangle$. Subsequent actions of the Hamiltonian will mix
with the orthogonal component of this vector and we need to determine the
correct linear combination describing the evolution.

To proceed from this point we thus have to find a family of vectors orthogonal
to $|A_{n}^{(0)}\rangle$ and then find the appropriate linear combination of
these vectors. The first step is easy to perform in this example as we can
determine it uniquely

$\left|B_{n}^{(1)}\right\rangle=\frac{(\beta^{*}L_{+}^{n}-\alpha^{*}L_{+}^{n-1}l_{+})|h,0,M\rangle}{\sqrt{|\beta|^{2}\langle
L_{-}^{n}L_{+}^{n}\rangle+|\alpha|^{2}\langle
L_{-}^{n-1}l_{-}l_{+}L_{+}^{n-1}\rangle}}$ (32)

The second requires finding the appropriate linear combination

$\left|A_{n}^{(1)}\right\rangle=\frac{1}{\sqrt{1+|N_{1,n}|^{2}}}\left(|A_{n}^{(0)}\rangle+N_{1,n}|B_{n-1}^{(1)}\rangle\right)$
(33)

The coefficients $N_{1,n}$ can be determined as

$N_{1,n}=\frac{\left\langle
B_{n-1}^{(1)}|H^{n}|h,0,M\right\rangle}{\left\langle
A_{n}^{(0)}|H^{n}|h,0,M\right\rangle}$ (34)

The formulae (33), (34) constitute the first order approximation to the Krylov
basis in this scheme. From this point on we proceed by iteration. Suppose that
the $i^{\mathrm{th}}$ order approximation to the Krylov basis,
$|A_{n}^{(i)}\rangle$ has been computed. We first find the family of vectors,
$|B_{n}^{(i)}\rangle$, satisfying the orthonormality constraint

$\langle A_{n}^{(i)}|B_{m}^{(i+1)}\rangle=0\ \ \ ;\ \ \ \forall\ n,m$ (35)

and then the appropriate linear combination of these vectors representing the
next order approximation

$\displaystyle\left|A_{n}^{(i+1)}\right\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{1+|N_{i+1,n}|^{2}}}\left(\left|A_{n}^{(i)}\right\rangle+N_{i+1,n}\left|B_{n-\lceil\frac{i}{2}\rceil}^{(i+1)}\right\rangle\right)$
$\displaystyle N_{i+1,n}$ $\displaystyle=$ $\displaystyle\frac{\left\langle
B_{n-\lceil\frac{i}{2}\rceil}^{(i+1)}|H^{n}|h,0,M\right\rangle}{\left\langle
A_{n}^{(i)}|H^{n}|h,0,M\right\rangle}$ (36)

Note that the vector $|B_{n}^{(i+1)}\rangle$ is out of step with (i.e. has a
different $n$-dependence to) the vector $|A_{n}^{(i)}\rangle$. The size of the
step, $j$ is determined by finding the the smallest integer for which $\langle
h,0,M|H^{n}|B_{n-j}^{(i+1)}\rangle\neq 0$. As with the Lanczos algorithm, this
iterative procedure is aimed at obtaining an orthonormal Krylov basis.
However, each step in the iteration process produces not just a single vector,
but a family of vectors with increasing precision as an approximation to the
true Krylov basis. The Hamiltonian expressed in this basis is of a tri-
diagonal form at each step in the iterative process. Furthermore, note that
the set $\left\\{|A_{n}^{(m)}\rangle\ ,\ n=0,1,2,\cdots\right\\}$ contains the
exact Krylov vector $|K_{0}\rangle,|K_{1}\rangle,\cdots|K_{m}\rangle$.

###  3.2 Computation using coherent states

As was the case in the context of the $su(1,1)$ and Heisenberg Hamiltonians,
coherent states provide an efficient computational tool for the implementation
of the iterative procedure outlined above. Here again, it is because the
vectors are replaced by differential operators acting on holomorphic states,
and overlaps can be computed by acting with these differential operators on
scalar functions.

As a starting point, note that with judicious usage of the BCH formulas, the
time-evolution of the reference state under the Hamiltonian (29) is given by,

$\displaystyle e^{itH}|h,0,M\rangle$ (37) $\displaystyle=$
$\displaystyle\left(\cosh(\sqrt{\zeta}t)-i\frac{i\gamma\sinh(t\sqrt{\zeta})}{2\sqrt{\zeta}}\right)^{-2h}\exp\left\\{M\frac{|2\alpha\beta^{*}-\beta\gamma|^{2}}{4\zeta^{2}\cosh(\sqrt{\zeta}t)-2i\zeta^{\frac{3}{2}}\gamma\sinh(\sqrt{\zeta}t)}\right\\}$
$\displaystyle\exp\left\\{M\left(-it\frac{\alpha^{*}\beta^{2}+\alpha(\beta^{*})^{2}-|\beta|^{2}\gamma}{2\zeta}+\frac{\gamma(\alpha^{*}\beta^{2}+\alpha(\beta^{*})^{2})-4|\alpha|^{2}|\beta|^{2}}{4|\alpha|^{2}\zeta}\right)\right\\}$
$\displaystyle\exp\left\\{M\frac{((4|\alpha|^{2}+\gamma^{2})(\alpha^{*}\beta^{2}+\alpha(\beta^{*})^{2})-8|\alpha|^{2}|\beta|^{2}\gamma)(4\sqrt{\zeta}\gamma+i|\alpha|^{2}\sinh(2\sqrt{\zeta}t))}{2\zeta^{\frac{3}{2}}(16\zeta^{2}-\gamma^{4}+|\alpha|^{4}\cosh(2\sqrt{\zeta}t))}\right\\}$
$\displaystyle\exp\left\\{\frac{-i(2\alpha\beta^{*}(\cosh(\sqrt{\zeta}t)-1))-2\sqrt{\zeta}\beta\sinh(\sqrt{\zeta}t)}{2i\zeta\cosh(\sqrt{\zeta}t)+\sqrt{\zeta}\gamma\sinh(\sqrt{\zeta}t)}l_{+}\right\\}\times$
$\displaystyle\exp\left\\{-\frac{2\alpha}{\gamma+2i\sqrt{\zeta}\coth(\sqrt{\zeta}t)}L_{+}\right\\}|h,0,M\rangle\,,$

with $\zeta\equiv|\alpha|^{2}-\frac{\gamma^{2}}{4}$. This in turn may be
combined with (28) to give a bulky though explicit expression for the overlap

$F(\bar{z}_{1},\bar{z}_{2},t)=\langle
h,0,M|e^{\bar{z_{2}}L_{+}}e^{\bar{z}_{1}l_{+}}e^{itH}|h,0,M\rangle$ (38)

This function serves as a generating functional for the overlaps $\langle
A_{n}^{(i)}|H^{n}|h,0,M\rangle$ and $\langle
B_{n-\lceil\frac{i}{2}\rceil}^{(i+1)}|H^{n}|h,0,M\rangle$ needed for (36). In
addition, the coherent state overlap (28) and the expectation value of the
Hamiltonian

$\displaystyle(\bar{z}_{1},\bar{z}_{2}|H|z_{1},z_{2})$ (39) $\displaystyle=$
$\displaystyle
e^{M\frac{\bar{z}_{1}^{2}z2+\bar{z}_{2}z_{1}^{2}+2\bar{z}_{1}z_{2}}{2(1-\bar{z}_{2}z_{2})}}(1-\bar{z}_{2}z_{2})^{-2(h+1)}\left(\alpha\left(\frac{M}{2}(\bar{z}_{1}+z_{1}\bar{z}_{2})^{2}+2h\bar{z}_{2}(1-\bar{z}_{2}z_{2})\right)\right.$
$\displaystyle+$ $\displaystyle\beta
M(\bar{z}_{1}+z_{1}\bar{z}_{2})(1-\bar{z}_{2}z_{2})+\gamma\left(\frac{M}{2}(z_{1}+\bar{z}_{1}z_{2})(\bar{z}_{1}+z_{1}\bar{z}_{2})+h(1-\bar{z}_{2}^{2}z_{2}^{2})\right)$
$\displaystyle+$
$\displaystyle\left.\alpha^{*}\left(\frac{M}{2}(z_{1}+\bar{z}_{1}z_{2})^{2}+2hz_{2}(1-\bar{z}_{2}z_{2})\right)+\beta^{*}M(z_{1}+\bar{z}_{1}z_{2})(1-\bar{z}_{2}z_{2})\right)\,,$

will be important. In the language of Jacobi coherent states the Krylov basis
vectors may be represented by differential operators acting on the state
$|z_{1},z_{2})$. This is analogous to the analysis performed for the $SU(1,1)$
and Heisenberg groups in subsection 2.1. Here the "0th" order Krylov vectors
correspond to the differential operator

$\displaystyle\hat{A}_{n}^{(0)}$ $\displaystyle=$
$\displaystyle\frac{e^{-ia}}{\sqrt{4h|\alpha|^{2}+2nM|\beta|^{2}}\sqrt{n!(2h+1)_{(n-1)}}}\left(\alpha\partial_{z_{2}}^{n}+\beta
n\partial_{z_{2}}^{n-1}\partial_{z_{1}}\right)\,,$ (40)

with $|A_{n}^{(0)}\rangle=\hat{A}_{n}^{(0)}|z_{1},z_{2})|_{z_{1}\rightarrow
0,z_{2}\rightarrow 0}$. Written in the coherent state language, the algorithm
now proceeds as follows: Given $\hat{A}_{n}^{(i)}$ such that

$\left.(\hat{A}_{m}^{(i)})^{\dagger}\hat{A}_{n}^{(i)}(\bar{z}_{1},\bar{z}_{2}|z_{1},z_{2})\right|_{z_{1},\bar{z}_{1}\rightarrow
0,z_{2},\bar{z}_{2}\rightarrow 0}=\delta_{m,n}\,,$ (41)

we need to find the family of differential operators $\hat{B}_{n}^{(i+1)}$
such that

$\left.(\hat{B}_{m}^{(i+1)})^{\dagger}\hat{A}_{n}^{(i)}(\bar{z}_{1},\bar{z}_{2}|z_{1},z_{2})\right|_{z_{1},\bar{z}_{1}\rightarrow
0,z_{2},\bar{z}_{2}\rightarrow 0}=0\,,$ (42)

for all $n,m$. The next order approximation is then given by

$\displaystyle\hat{A}_{n}^{(i+1)}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{1+|N_{i+1,n}|^{2}}}\left(\hat{A}_{n}^{(i)}+N_{i+1,n}\hat{B}_{n-\lceil\frac{i}{2}\rceil}^{(i+1)}\rangle\right)\,,$
$\displaystyle N_{i+1,n}$ $\displaystyle=$
$\displaystyle\left(\left.\hat{B}^{i+1}_{n-\lceil\frac{i}{2}\rceil}\partial_{t}^{n}F(\bar{z}_{1},\bar{z}_{2},t)\right|_{\bar{z}_{1}\rightarrow
0,\bar{z}_{2}\rightarrow 0,t\rightarrow
0}\right)\left(\left.\hat{A}^{i}_{n}\partial_{t}^{n}F(\bar{z}_{1},\bar{z}_{2},t)\right|_{\bar{z}_{1}\rightarrow
0,\bar{z}_{2}\rightarrow 0,t\rightarrow 0}\right)^{-1}\,.$

When terminated at order $j$ the Lanczos coefficients may thenbe extracted as

$\displaystyle b^{(j)}_{n}$ $\displaystyle=$
$\displaystyle(\hat{A}^{(j)}_{n})^{\dagger}(\hat{A}^{(j)}_{n-1})\left.(\bar{z}_{1},\bar{z}_{2}|H|z_{1},z_{2})\right|_{z_{1},\bar{z}_{1}\rightarrow
0,z_{2},\bar{z}_{2}\rightarrow 0}$ $\displaystyle a^{(j)}_{n}$
$\displaystyle=$
$\displaystyle(\hat{A}^{(j)}_{n})^{\dagger}(\hat{A}^{(j)}_{n})\left.(\bar{z}_{1},\bar{z}_{2}|H|z_{1},z_{2})\right|_{z_{1},\bar{z}_{1}\rightarrow
0,z_{2},\bar{z}_{2}\rightarrow 0}$ (43)

##  4 Jacobi Lanczos coefficients and complexity

We are now in a position to compute the approximate Lanczos coefficients for
the Hamiltonian (29). In the next few sections we will compute these at 0th,
1st and 2nd order for some choices of parameters and compare these to the
exact Lanczos coefficients which we determine numerically. Following this we
will be able to compute the corresponding approximate Krylov complexities. A
nice feature is that, since the terms appearing in the sum are known for all
values of $n$, we have good control over the accuracy of the results.

###  4.1 Lanczos Coefficients

We now proceed to apply the algorithm outlined above up to second order,
setting $h=\frac{1}{4}$ in line with (26). For ease of notation, we set

$w_{n}=\sqrt{|\alpha|^{2}+2Mn|\beta|^{2}}.$ (44)

We find the families of differential operators

$\displaystyle\hat{B}_{n}^{(1)}$ $\displaystyle=$
$\displaystyle\frac{\alpha^{*}}{\alpha\sqrt{M}\sqrt{16h^{2}\alpha^{2}+8nhM\beta^{2}}\sqrt{(n-1)!(2h+1)_{(n-1)}}}\left(2M\beta^{*}\partial_{z_{2}}^{n}-4\alpha^{*}h\partial_{z_{2}}^{n-1}\partial_{z_{1}}\right)$
$\displaystyle A_{n}^{(1)}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{1+|N_{2,n}|^{2}}}\left(\hat{A}_{n}^{(0)}+N_{2,n}\hat{B}^{(1)}_{n-1}\right)$
$\displaystyle B_{n}^{(2)}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{1+|N_{2,n}|^{2}}}\left(N_{2,n}\hat{A}_{n}^{(0)}-\hat{B}^{(1)}_{n-1}\right)$
$\displaystyle A_{n}^{(2)}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{1+|N_{3,n}|^{2}}}\left(\hat{A}_{n}^{(1)}+N_{3,n}\hat{B}^{(2)}_{n-1}\right)$
(45)

with coefficients

$N_{2,k}\\\
=\frac{\sqrt{Mk(k-1)}(4(2k-1)M\beta^{3}-6\alpha^{2}\beta^{*}+3\alpha\beta\gamma))}{6\alpha\sqrt{2k-1}w_{k}w_{k-1}}$

and

$\displaystyle N_{3,k}$ $\displaystyle=$ $\displaystyle\frac{\sqrt{M\
k(k-1)(k-2)}}{180\sqrt{2}\alpha^{3}\beta\sqrt{(1+|N_{2,k-1}|^{2})(1+|N_{2,k-1}|^{2})}\sqrt{(2k-3)(2k-1)}w_{k}w_{k-1}^{2}w_{k-2}}\times$
$\displaystyle\left[(5|\alpha|^{2}+w_{k-1}^{2})(4M\beta^{3}(2k-1)-6\alpha^{2}\beta^{*}+3\alpha\beta\gamma)(4M\beta^{3}(2k-3)-6\alpha^{2}\beta^{*}+3\alpha\beta\gamma)\right.$
$\displaystyle\left.-36\alpha(2\alpha\beta^{*}-\beta\gamma)(2M\beta^{3}(k-1)+3\alpha^{2}\beta^{*}+\alpha\beta\gamma)w_{k-1}^{2}\right]$

These expressions may be readily obtained using a simple Mathematica code
implementing the algorithm of the previous subsection. The Lanczos
coefficients may subsequently be obtained as

$\displaystyle b_{n}^{(2)}$ $\displaystyle=$
$\displaystyle\frac{\sqrt{n(n-\frac{1}{2})}w_{n}\alpha}{\sqrt{(1+|N_{2,n-1}|^{2})(1+|N_{2,n}|^{2})(1+|N_{3,n-1}|^{2})(1+|N_{3,n}|^{2})}w_{n-1}}$
(46) $\displaystyle+$
$\displaystyle\frac{\sqrt{(2n-3)(n-2)}w_{n-1}w_{n-2}N_{2,n-1}N_{2,n}^{*}\alpha+\sqrt{M(n-1)}(4(n-1)M\beta^{3}-2\alpha^{2}\beta^{*}+\alpha\beta\gamma)N_{2,n}^{*}}{\sqrt{2}\sqrt{(1+|N_{2,n-1}|^{2})(1+|N_{2,n}|^{2})(1+|N_{3,n-1}|^{2})(1+|N_{3,n}|^{2})}w_{n-1}^{2}}$
$\displaystyle+$
$\displaystyle\frac{N_{3,n-1}N_{3,n}^{*}\left(\sqrt{(2n-5)(n-3)}w_{n-3}w_{n-2}\alpha+\sqrt{(2n-3)(n-1)}w_{n-1}w_{n-2}N_{2,n-1}N_{2,n-2}^{*}\alpha\right)}{\sqrt{2}\sqrt{(1+|N_{2,n-2}|^{2})(1+|N_{2,n-1}|^{2})(1+|N_{3,n-1}|^{2})(1+|N_{3,n}|^{2})}w_{n-2}^{2}}$
$\displaystyle-$
$\displaystyle\frac{N_{3,n-1}N_{3,n}^{*}\sqrt{M(n-1)}\left(4(n-2)M\beta^{3}-2\alpha^{2}\beta^{*}+\alpha\beta\gamma\right)N_{2,n-2}^{*}}{\sqrt{2}\sqrt{(1+|N_{2,n-2}|^{2})(1+|N_{2,n-1}|^{2})(1+|N_{3,n-1}|^{2})(1+|N_{3,n}|^{2})}w_{n-2}}$
$\displaystyle+$
$\displaystyle\frac{N_{3,n}^{*}(-2\sqrt{(n-2)(2n-3)(n-1)}w_{n-1}w_{n-2}M^{\frac{3}{2}}|\beta|^{2}(\beta-\beta^{*}N_{2,n-1}^{2})+\frac{1}{2}\gamma
w_{n-1}^{2}w_{n-2}^{2}N_{2n-1})}{(1+|N_{2,n-1}|^{2})\sqrt{(1+|N_{3,n-1}|^{2})(1+|N_{3,n}|^{2})}w_{n-1}^{2}w_{n-2}^{2}}$
$\displaystyle+$
$\displaystyle\frac{N_{3,n}^{*}N_{2,n-1}(2M((2n-2)w_{n-2}^{2}-|\alpha|^{2})(\alpha^{*}\beta^{2}+\alpha(\beta^{*})^{2})+\gamma|\alpha|^{2}w_{n-\frac{3}{2}}^{2})}{(1+|N_{2,n-1}|^{2})\sqrt{(1+|N_{3,n-1}|^{2})(1+|N_{3,n}|^{2})}w_{n-1}^{2}w_{n-2}^{2}}$

The expression is organised so that the first line represents our 0th order
approximation, the second line the 1st order approximation and the rest the
2nd order approximation. The 0th order may be extracted from the above by
setting $N_{2,n},N_{3,n}\rightarrow 0$ and the 1st order by setting
$N_{3,n}\rightarrow 0$. These expressions are bulky, reflecting the involved
combinatorics determining the Krylov basis. Nevertheless, these are explicit
expressions that we may test numerically. We may also compare these with an
explicit analytical computation of the first few orders of the algorithm, as
described in appendix C.

(a)

(b)

(c)

(d)

Figure 1: The relative error for the Lanczos coeffcients
$\left|\frac{b_{n}-b_{n}^{i}}{b_{n}}\right|$ for $i=0$ (blue), $i=1$ (orange)
and $i=2$ (green) with
$|\alpha|=1,\gamma=0,\frac{\alpha}{\alpha^{*}}=\frac{\beta}{\beta^{*}}=e^{2i\frac{\pi}{3}}$.

(a)

(b)

(c)

(d)

Figure 2: The relative error for the Lanczos coeffcients
$\left|\frac{b_{n}-b_{n}^{i}}{b_{n}}\right|$ for $i=0$ (blue), $i=1$ (orange)
and $i=2$ (green) with
$|\alpha|=1,\beta=\frac{1}{4},\frac{\alpha}{\alpha^{*}}=e^{i\pi},\frac{\beta}{\beta^{*}}=1$.

In Figs. (1) and (2) we have plotted the relative error of terminating the
algorithm at 0th, 1st and 2nd order. The Lanczos coefficients are computed
numerically and compared with (46). We note that, even for values
$|\beta|\approx|\alpha|$ the leading order approximation provides an error of
less than $10$ percent with greater accuracy at smaller values of $|\beta|$.
At subsequent orders the approximation increases in accuracy, especially for
smaller values of $n$. As advertised, the algorithm perfectly reproduces
$b_{j+1}$ at order $j$.

###  4.2 K-complexity

Given a Krylov basis $\\{|K_{n}\rangle\\}$, we define the Krylov complexity
operator as,

$\hat{K}=\sum_{n}n|K_{n}\rangle\langle K_{n}|.$ (47)

The associated K-complexity of a particular target state is then given by the
expectation value of this operator in this state. In our case, we will use as
the target state the time-evolved reference state (37) which (up to an overall
normalisation) is a Jacobi coherent state with values

$\displaystyle z_{1}(t)$ $\displaystyle=$
$\displaystyle\frac{(2\alpha\beta^{*}-\beta\gamma)(\cosh(\sqrt{\zeta}t)-1)+2i\beta\sqrt{\zeta}\sinh(\sqrt{\zeta}t)}{2\zeta\cosh(\sqrt{\zeta}t)-i\gamma\sqrt{\zeta}\sinh(\sqrt{\zeta}t)}\,,$
$\displaystyle z_{2}(t)$ $\displaystyle=$
$\displaystyle-\frac{2\alpha}{\gamma+2i\sqrt{\zeta}\coth(\sqrt{\zeta}t)}\,,$
(48) $\displaystyle\zeta$ $\displaystyle\equiv$
$\displaystyle|\alpha|^{2}-\frac{\gamma^{2}}{4}\,.$

In terms of the differential operators we write

$K^{(i)}=\sum_{n}n\frac{1}{(\bar{z}_{1}(t),\bar{z}_{2}(t)|z_{1}(t),z_{2}(t))}\left.(A_{n}^{(i)})^{\dagger}(A_{n}^{(i)})\left(\bar{z}_{1},\bar{z}_{2}|z_{1}(t),z_{2}(t)\right)\left(\bar{z}_{1}(t),\bar{z}_{2}(t)|z_{1},z_{2}\right)\right|_{z_{1},z_{2}\rightarrow
0}$ (49)

which is the approximate spread complexity at "ith" order. A useful formula to
compute the above is

$\displaystyle\left.\partial_{\bar{z}_{2}}^{n}\partial_{\bar{z}_{1}}^{2m}(\bar{z}_{1},\bar{z}_{2}|z_{1},z_{2})\right|_{\bar{z}_{1}\rightarrow
0,\bar{z}_{2}\rightarrow 0}$ $\displaystyle=$
$\displaystyle\sum_{k=0}^{n+m}2^{m-k}\left(\begin{array}[]{c}n+m\\\
k\end{array}\right)\frac{\Gamma(n+m+\frac{1}{2})}{\Gamma(k+\frac{1}{2})}z_{1}^{2k}z_{2}^{n+m-k}$
(52)
$\displaystyle\left.\partial_{\bar{z}_{2}}^{n}\partial_{\bar{z}_{1}}^{2m+1}(\bar{z}_{1},\bar{z}_{2}|z_{1},z_{2})\right|_{\bar{z}_{1}\rightarrow
0,\bar{z}_{2}\rightarrow 0}$ $\displaystyle=$
$\displaystyle\sum_{k=0}^{n+m}2^{m-k}\left(\begin{array}[]{c}n+m\\\
k\end{array}\right)\frac{\Gamma(n+m+\frac{3}{2})}{\Gamma(k+\frac{3}{2})}z_{1}^{2k+1}z_{2}^{n+m-k}$
(55)

which may be written formally as a hypergeometric function.

Before proceeding further, it is worthwhile to reflect on the physical
expectations for Krylov complexity. The Hamiltonians we consider are similar
to harmonic oscillator and inverted harmonic oscillator studied previously in
the complexity literature Ali:2019zcj ; Bhattacharyya:2020art ; Haque:2021kdm
(see also the $su(1,1)$ example in Caputa:2021sib ; Haque:2021hyw ). Indeed,
when using the creation / annihilation operator representation, these systems
are related by a shift of the creation / annihilation operators and Bogoliubov
transformation. The frequency of the oscillator is given by
$\omega=\sqrt{\frac{\gamma^{2}}{4}-|\alpha|^{2}}$.

When $\omega$ is real we expect complexity to be a periodic function in time.
Note that this feature is already guaranteed by the parametrisation of the
time-evolved state as a Jacobi coherent state (48). For real $\omega$ the
hyperbolic trigonometric functions become ordinary trigonometric functions. In
contrast, when $\omega$ is imaginary we expect complexity similar to that of
the inverted oscillator which should grow exponentially. This may again be
expected by the form of the parametrisation (48) which saturates at late times
and, furthermore, the large $n$ behavior of the Lanczos coefficients (46) is
$b_{n}^{2}=\alpha n+o(1)$. Put together we expect exponential growth of Krylov
complexity at late time.

(a)

(b)

(c)

Figure 3: The approximate spread complexity, $K^{(i)}$ of the time-evolved
reference state at leading, subleading and subsubleading order in the scheme.
We have taken $\alpha=1$, $\beta=\frac{1}{2}$, $\gamma=1$ and
$b=\frac{\pi}{2}$. The red line is a plot of the expectation value (56) for
comparison.

(a)

(b)

(c)

Figure 4: The approximate spread complexity, $K^{(i)}$ of the time-evolved
reference state at leading, subleading and subsubleading order in the scheme.
We have taken $\alpha=1$, $\gamma=0$ and $a=b=\frac{\pi}{2}$. The red line is
a plot of the expectation value (56) for comparison.

(a)

(b)

(c)

Figure 5: The approximate spread complexity, $K^{(i)}$ of the time-evolved
reference state at leading, subleading and subsubleading order in the scheme.
We have taken $|\alpha|=|\beta|=1$, $\gamma=3$ and $b=\frac{\pi}{2}$. A little
more than one period of the periodic approximate spread complexity is shown.

In Figs. (3), (4) we have plotted the spread complexity for a range of
parameters corresponding to an imaginary frequency while those in Fig. (5)
correspond to real frequency. For comparison we have plotted the expectation
value

$\tilde{C}=\sum_{k}k|\langle h,k,M|z_{1}(t),z_{2}(t)\rangle|^{2}+(k+2)|\langle
h,k+\frac{1}{2},M|z_{1}(t),z_{2}(t)\rangle|^{2}$ (56)

which is always greater than (49) for $i=0,1,2$ for all parameter choices that
we have numerically checked. The physical relevance of this quantity needs
some unpacking. First, note that it is slightly larger than half the value for
the cost function

$C_{B}=\sum_{n}n\langle B_{n}|z_{1}(t),z_{2}(t)\rangle$ (57)

with $|B_{n}\rangle=|h,\frac{n}{2},M\rangle$. This is the basis one would use
with only access to the Heisenberg generators. However, since the time
evolution involves the $su(1,1)$ generators (which ladder up faster than the
Heisenberg generators) it is easy to understand why the complexity is less
than the cost function, $C_{B}$. Indeed, the expectation value (56) may be
interpreted as a cost function with the basis generated by the $su(1,1)$
generators from a "complexity $0$ state" $|h,0,M\rangle$ and a "complexity $2$
state" $|h,\frac{1}{2},M\rangle$. This value of $2$ is allocated since both
$H|h,0,M\rangle$ and $H^{2}|h,0,M\rangle$ are needed to produce both relevant
linear combinations of $|h,\frac{1}{2},M\rangle$ and $|h,1,M\rangle$ present
in the problem555If the state $|h,\frac{1}{2},M\rangle$ is assigned a
”complexity” of $1$ the cost function may be smaller than $K^{(i)}$ at early
times.. The complexity we obtain is less than $\tilde{C}$. One way to think
about this is that the Heisenberg generators provide additional ways to
synthesise the target state. For a choice of parameters corresponding to
imaginary frequencies we observe that the complexity is an increasing function
of time with exponential growth at later times. If, instead, the parameters
are chosen to correspond to a real frequency, the complexity is a periodic
function in time. Both of these results are in line with our physical
expectations.

Another noteworthy feature we would like to highlight is the appearance of an
inflection point in $K^{(i)}$ at intermediate times and for some choices of
parameters, see Figs. (3), (4). Though we do not yet fully understand its
appearance, intuitively, it may be understood as the result of the fact that a
particular target state is obtainable in various ways using the ladder
operators $L_{+}$ and $l_{+}$. As we follow the time-evolved state, its
decomposition into the various Krylov basis components will thus involve many
powers of $z_{1}(t)$ and $z_{2}(t)$, which grow differently as functions of
time. As such, it is possible for the gradient of complexity to change,
depending on the precise coefficients appearing in front of various powers of
$z_{1}(t)$ and $z_{2}(t)$. A similar feature may be noted for a family of
geodesics in the context of Nielsen complextiy Nielsen_2006 ; Nielsen1 ;
Nielsen3 , where the tangent vector to evolution on the space of unitaries may
be an arbitrary combination of $l_{+}$ and $L_{+}$. The geodesic length may
increase differently along these different directions (for a simple example of
this see appendix B).

What may be surprising is that the appearance of the inflection point appears
to be sensitive to the phases $a,b$, even for the leading order $K^{(0)}$.
This is despite the Lanczos coefficient at this order (46) not being sensitive
to the phase. Here it is important to stress that, though the operator
$K^{(i)}$ is an approximate K-complexity operator, the states
$|z_{1}(t),z_{2}(t)\rangle$ are the exact time-evolved states, and the
dependence on the phases enters through it. Indeed, we have noted that, in
order for the inflection point to appear, a necessary condition is that
$|z_{1}(t)|$ is comparable to $|z_{2}(t)|$ at intermediate times.

###  4.3 Generalisations

We would like to emphasise that the recursive method, in vector form in (35),
(36) and in differential form in (42), (3.2), is a general algorithm to
generate the Krylov basis. Provided that a natural guess for a
“$0^{\mathrm{t}h}$ order” basis exists, one may iterate to generate a sequence
of bases that converge on the Krylov basis. We have applied this in the
context of Jacobi coherent states and the resulting expressions capture the
involved combinatorics associated with the Gram-Schmidt procedure. The Jacobi
group provides three simplifications that would not be available in a more
general setting. We would like to comment on the application of our proposed
algorithm in the absence of these simplifications.

Firstly, for the Jacobi group we were able to obtain an explicit realisation
of the time-evolved state in terms of Jacobi coherent states (37). This was
used to obtain a generating function (38) for the overlaps of
$H^{n}|\psi_{r}\rangle$ with the $|A_{n}^{(i)}\rangle$ vectors in the
approximation scheme. For generic Hamiltonians this would not be possible
though the action of $H^{n}$ may always be represented in terms of
differential operators acting on coherent states. These overlaps can still be
obtained, even without the convenient generating functional.

Secondly, for a Hamiltonian from the Jacobi algebra the dimension of the span
of $\left\\{H^{n}|\phi_{0}\rangle\right\\}$ is exactly half that of
$\left\\{|h,n,M\rangle,|h,n-\frac{1}{2},M\rangle\right\\}$. Because of this
the vectors $|B_{n}^{(i)}\rangle$ could be obtained with relative ease at each
step in the recursive algorithm. In general there could be several families of
vectors orthogonal to $|A_{n}{(i)}\rangle$ at each step. In this case one
would have to find appropriate linear combinations of these vectors, similar
to what was performed for the $N_{i+1,n}$ coefficients.

Finally, our Hamiltonian was chosen as an element of some finite-dimensional
algebra. General Hamiltonians will not be of this type. However, it may be
possible to select a reference state which transforms trivially under many
terms present in the Hamiltonian. In such cases, a natural choice for the
“$0^{\mathrm{th}}$” order basis should be possible.

##  5 Discussion

Spread complexity - a measure of the difficulty of preparing a desired target
quantum state from a simple reference state - has emerged as a useful
diagnostic for various quantum phenomena, from chaos to phase transitions. It
naturally involves the spectrum of the system under consideration and requires
only the reference state, target state and system Hamiltonian to be specified.
Finding the Krylov basis associated with a particular choice of Hamiltonian
and reference state is one of the central tasks in computing the complexity.
For simple symmetry groups the Krylov basis can readily be identified in terms
of the ladder operators. In general systems, the Lanczos algorithm provides an
iterative way to implement the required Gram-Schmidt orthogonalisation
procedure.

In this article we have proposed a modification of the Lanczos algorithm which
we anticipate will be well-suited to a variety of problems. Starting from a
basis assembled from from the ladder operators present in the problem, the
basis may be iteratively refined until ultimately converging on the Krylov
basis. A key advantage of this modification is that an approximate basis for
the Krylov subspace is available at each step in the iteration. For an
infinite-dimensional Krylov subspace this is especially useful to approximate
the infinite sums involved in computing the spread complexity in a stable way.

As an illustration, we have applied the algorithm to Jacobi coherent states,
where some degree of analytical control is available. Analytic expressions may
be obtained for the approximate Lanczos coefficients at each step of the
algorithm and we have demonstrated that the approximation is already accurate
at subsubleading order for a range of parameters. With the aid of analytic
expressions, the approximate Krylov may be computed with a simple numerical
algortihm. The behavior we obtain is in accord with physical expectations.
Notably, we have found that the K-complexity for the semi-direct product group
is less than that of its constitutent subgroups. At first sight, this seems
counter intuitive. One way to think about this is by comparison to the
equivalent result in circuit complexity666In Appendix B, we illustrate some
features of the circuit complexity for the Jacobi group. where, in enlarging
the symmetry group we are effectively adding more gate sets thereby creating
new pathways in the space of unitaries with potentially shorter geodesics
between the initial and target states. In this sense then, quantum systems
with this more generic quadratic Hamiltonian are less complex than systems
evolved by the constituent subsystem Hamiltonians. We anticipate that this
result might have interesting practical applications for quantum computations.
We will report on these findings in detail in a future work.

We anticipate that this algorithm may be applied to Hamiltonians drawn from
more general symmetry groups and may even be applicable to problems where
perturbative treatments are unavoidable. One example that would be of great
interest is the $d$-dimensional conformal group $SO(d,2)$. Results for the
Nielsen complexity of conformal primaries are known Chagnet:2021uvi ;
Koch:2021tvp ; Rabambi:2022jwu , and it would be interesting to compare these
results with Krylov complexity. A key question in this context is whether the
spread-complexity of a time-evolved state under an $SO(d,2)$ Hamiltonian
saturates and, if so, how dependent is this property on the choice of
reference state.

It would be useful to develop a better physical intuition for the information
contained in the spread complexity of a quantum state. For example, one might
like to know whether (and under what circumstances) spread complexity is
sensitive to phase transition or integrable-to-chaotic transitions. These have
been investigated for some models, e.g. Jian:2020qpp ; Kim:2021okd ;
Rabinovici:2021qqt ; Rabinovici:2022beu ; Caputa:2022eye ; Caputa:2022yju
(see also Bhattacharjee:2022vlt ), but the list may be extended substantially.
The tools we have developed in this paper are taking a step towards performing
spread complexity computations for these more general cases.

On a more technical level, the qualitative results for the Lanczos
coefficients and spread complexity of the Jacobi group coherent states appear
in good agreement with the (numerically obtained) exact results. The algorithm
would, however, benefit from tools that would allow one to estimate the error.
At least in the context of symmetry groups we are optimistic that this may be
done. If a choice of reference state can be made for which the complexity is
analytically computable, a more general choice can be made by performing a
unitary transformation of either the state or Hamiltonian (but not both). It
may be possible to estimate the error by imposing a suitable metric on the
space of such unitary transformations.

Finally, in this paper we have focused on the example of a time-independent
Hamiltonian. Extending to Hamiltonians that carry an explicit time-dependence
would be an important generalisation of this modified Lanczos method.

###### Acknowledgements.

We would like to thank Pawel Caputa, Bret Underwood and Arpan Bhattacharyya
for discussions and comments. JM would like to acknowledge support from the
ICTP through the Associates Programme and from the Simons Foundation through
grant number 284558FY19. SH would like to thank Brac University for their
hospitality during the final stage of this work. H.J.R.vZ is supported by the
“Quantum Technologies for Sustainable Devlopment" grant from the National
Institute for Theoretical and Computational Sciences (NITHECS).

##  Appendix A Special choice of reference state

In this appendix we would like to highlight a choice of reference state for
which the Krylov complexity is exactly solvable. This choice of reference
state is dependent on the coefficients appearing in the Hamiltonian and
therefore, obviously a special choice. Nevertheless, it is useful to recognize
that there exists reference states for which the problem is solvable. We first
observe that we may identify various $su(1,1)$ subalgebras of the Jacobi
algebra

$\displaystyle\tilde{L}_{+}$ $\displaystyle=$ $\displaystyle
L_{+}+cl_{+}+\frac{c^{2}}{2}l_{0}\,,$ $\displaystyle\tilde{L}_{-}$
$\displaystyle-$ $\displaystyle
L_{-}+c^{*}l_{-}+\frac{(c^{*})^{2}}{2}l_{0}\,,$ (58)
$\displaystyle\tilde{L}_{0}$ $\displaystyle=$ $\displaystyle
L_{0}+\frac{c^{*}}{2}l_{+}+\frac{c}{2}l_{-}+\frac{cc^{*}}{2}l_{0}\,.$

The Hamiltonian (29) may be written (up to an overall central element) as a
linear combination of these generators

$H=\alpha\tilde{L}_{+}+\alpha^{*}\tilde{L}_{-}+\gamma\tilde{L}_{0}-\left(\alpha\frac{c^{2}}{2}+\alpha^{*}\frac{(c^{*})^{2}}{2}+\gamma\frac{cc^{*}}{2}\right)l_{0}\,,$
(59)

with the choice

$c=\frac{4\alpha^{*}\beta-2\beta^{*}\gamma}{4|\alpha|^{2}-\gamma^{2}}\,.$ (60)

If one can find a reference state such that

$\displaystyle\tilde{L}_{-}|\phi_{r}\rangle$ $\displaystyle=$ $\displaystyle
0$ $\displaystyle\tilde{L}_{0}|\phi_{r}\rangle$ $\displaystyle=$
$\displaystyle\Delta|\phi_{r}\rangle$ (61)

the problem reduces to one of $SU(1,1)$ coherent states, for which exact
analytic results are known Caputa:2021sib . This choice of reference state is
given by

$|\phi_{r}\rangle=e^{-c^{*}l_{+}}|h,0,M\rangle$ (62)

for which it can be shown without too much trouble that $\Delta=0$. The exact
Krylov basis is given by

$|K_{n}\rangle=\mathcal{N}_{n}(\tilde{L}_{+})^{2n}e^{-c^{*}l_{+}}|h,0,M\rangle.$
(63)

Note that this choice of reference state for the Hamiltonian is unitarily
equivalent to setting $\beta\rightarrow 0$ with the reference state
$|h,0,M\rangle$ we have considered in the main text.

##  Appendix B Nielsen complexity for Jacobi coherent states

In this appendix we unpack some results for the Nielsen complexity of Jacobi
coherent states. Briefly, for Nielsen complexity one needs to specify a set of
continuous unitary gates, a reference state, target state and a “cost”
function on the space of unitary gates. Our choice of unitary gates is taken
as the Jacobi symmetry group. The accessible target states are those that
differ from the reference state by a Jacobi unitary transformation and are, in
fact, the Jacobi group coherent states. For the lowest weight reference state
$|h,0,M\rangle$ we have

$|z_{1},z_{2}\rangle=N(z_{1},z_{2})e^{z_{1}l_{+}}e^{z_{2}L_{+}}|h,0,M\rangle.$
(64)

As the cost function we will use the Fubini-Study metric, which may be
obtained from the coherent states above as

$g_{a\bar{b}}=\partial_{z_{a}}\partial_{\bar{z}_{b}}\log(\langle
h,0,M|e^{\bar{z}_{1}l_{-}}e^{\bar{z}_{2}L_{-}}e^{z_{1}l_{+}}e^{z_{2}L_{+}}|h,0,M\rangle).$
(65)

We have already computed the coherent state overlap in the main text in
eq(28). Consequently, the Fubini-Study metric is easily computed as

$\displaystyle ds^{2}$ $\displaystyle=$
$\displaystyle\frac{M}{1-|z_{2}|^{2}}dz_{1}d\bar{z}_{1}+\frac{M|z_{1}+\bar{z}_{1}z_{2}|^{2}+2h(1-|z_{2}|^{2})}{(1-|z_{2}|^{2})^{3}}dz_{2}d\bar{z}_{2}$
(66) $\displaystyle+$
$\displaystyle\frac{M(z_{1}+\bar{z}_{1}z_{2})}{2(1-|z_{2}|^{2})^{2}}dz_{1}d\bar{z}_{2}+\frac{M(\bar{z}_{1}+z_{1}\bar{z}_{2})}{2(1-|z_{2}|^{2})^{2}}d\bar{z}_{1}dz_{2}\,.$

The (Nielsen) complexity may now be computed as the minimal geodesic distance
connecting the reference and target states. On the $z_{1}=0$ submanifold the
geodesic may be solved as

$z_{2}(\sigma)=e^{i\phi}\tanh(c\sigma)$ (67)

recovering the results of (Chapman:2021jbh ) for the geodesic distance when
plugging this into

$L=\int_{\sigma_{min}}^{\sigma_{max}}\sqrt{\left|\sum_{a,b}g_{\bar{a}b}\bar{z}^{\prime}_{a}(\sigma)z^{\prime}_{b}(\sigma)\right|}\
d\sigma$ (68)

For general points on the manifold the geodesic equation may be solved
numerically. For real $z_{1},z_{2}$ we find that the geodesic length (i.e.
complexity) is a monotonically increasing function, see Fig (6). The gradient
for the geodesic length has a different slope in the $z_{1}$, $z_{2}$
directions.

Figure 6: The geodesic length connecting the origin to the point $z_{1},z_{2}$
on the manifold with both $z_{1}$ and $z_{2}$ real. We have set
$h=\frac{1}{4}$ and $M=1$. Note that the complexity increases with increasing
$|z_{i}|$ but that the gradient is different for the two directions.

Finally, in the same way that we noted a dependence on the phases $a,b$ in the
context of spread complexity, the Nielsen complexity also displays a
dependence on the phase of $z_{i}$, see Fig. (7)

Figure 7: The geodesic length connecting the origin to the point $z_{1},z_{2}$
on the manifold with $z_{1}=2$ and $z_{2}=0.8e^{i\phi}$. We have set
$h=\frac{1}{4}$ and $M=1$.

##  Appendix C Lanczos Algorithm for our Hamiltonian

In this appendix we demonstrate some of the involved combinatorics that arise
when the Lanczos algorithm is applied to the Jacobi Hamiltonian and highest
weight state considered in the main text. By definition we have as starting
point

$|K_{0}\rangle=|h,0,M\rangle$ (69)

where we set $h=\frac{1}{4}$. Acting on this state with the Hamiltonian gives

$H|K_{0}\rangle=\frac{\alpha}{\sqrt{2}}|h,1,M\rangle+\sqrt{M}\beta|h,\frac{1}{2},M\rangle+\frac{\gamma}{4}|h,0,M\rangle$
(70)

from which we obtain

$\displaystyle|A_{1})$ $\displaystyle=$
$\displaystyle\frac{\alpha}{\sqrt{2}}|h,1,M\rangle+\sqrt{M}\beta|h,\frac{1}{2},M\rangle$
$\displaystyle b_{1}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}\sqrt{|\alpha|^{2}+2M|\beta|^{2}}\,,$ (71)

and $|K_{1}\rangle=\frac{1}{b_{1}}|A_{1})$. The next step in the Lanczos
algorithm then yields

$\displaystyle H|K_{1}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}\sqrt{|\alpha|^{2}+2M|\beta|^{2}}|h,0,M\rangle+\frac{\sqrt{M}(4\alpha\beta^{*}+3\beta\gamma)}{2\sqrt{2}\sqrt{|\alpha|^{2}+2M|\beta|^{2}}}|h,\frac{1}{2},M\rangle+\frac{8M\beta^{2}+5\alpha\gamma}{4\sqrt{|\alpha|^{2}+2M|\beta|^{2}}}|h,1,M\rangle$
(72) $\displaystyle+$
$\displaystyle\frac{2\sqrt{3}\sqrt{M}\alpha\beta}{\sqrt{|\alpha|^{2}+2M|\beta|^{2}}}|h,\frac{3}{2},M\rangle+\frac{\sqrt{3}\alpha^{2}}{\sqrt{|\alpha|^{2}+2M|\beta|^{2}}}|h,2,M\rangle\,,$

from which we extract

$\displaystyle|A_{2})$ $\displaystyle=$
$\displaystyle\frac{4M\beta^{3}-2\alpha^{2}\beta^{*}+\alpha\beta\gamma}{\sqrt{2}(|\alpha|^{2}+2M|\beta|^{2})^{\frac{3}{2}}}\left(2M\beta^{*}|h,1,M\rangle-\sqrt{2}\sqrt{M}\alpha^{*}|h,\frac{1}{2},M\rangle\right)$
$\displaystyle+$
$\displaystyle\frac{\sqrt{3}\alpha}{\sqrt{|\alpha|^{2}+2M|\beta|^{2}}}\left(\alpha|h,2,M\rangle+2\sqrt{M}\beta|h,\frac{3}{2},M\rangle\right)\,,$
$\displaystyle b_{2}$ $\displaystyle=$
$\displaystyle\left(\frac{3|\alpha|^{2}(|\alpha|^{2}+4M|\beta|^{2})}{|\alpha|^{2}+2M|\beta|^{2}}+\frac{M|4M\beta^{3}-2\alpha^{2}\beta^{*}+\alpha\beta\gamma|^{2}}{2(|\alpha|^{2}+2M|\beta|^{2})^{2}}\right)^{\frac{1}{2}}\,,$
(73)

and $|K_{2}\rangle=\frac{1}{b_{2}}|A_{2})$. Note that this vector, compared to
$|K_{1}\rangle$ contains both a linear combination of
$|h,\frac{3}{2},M\rangle$, $|h,2,M\rangle$ as well as the linear combination
of $|h,\frac{1}{2},M\rangle$, $|h,1,M\rangle$ orthogonal to $|K_{1}\rangle$.
Similarly, the next Krylov basis vector will contain a linear combination of
$|h,\frac{5}{2},M\rangle$, $|h,3,M\rangle$ as well as a vector orthogonal to
$|K_{2}\rangle$ and $|K_{1}\rangle$. The specific linear combinations of
vectors become progressively more unwieldy as the algorithm progresses.

## References

  * (1) M.A. Nielsen and I.L. Chuang, _Quantum Computation and Quantum Information_ , Cambridge University Press (2000).
  * (2) J. Watrous, _Quantum computational complexity_ , _Computational Complexity: Theory, Techniques, and Applications_ (2008)  [0804.3401].
  * (3) S. Chapman and G. Policastro, _Quantum computational complexity from quantum information to black holes and back_ , _Eur. Phys. J. C_ 82 (2022) 128 [2110.14672].
  * (4) B. Chen, B. Czech and Z.-z. Wang, _Quantum information in holographic duality_ , _Rept. Prog. Phys._ 85 (2022) 046001 [2108.09188].
  * (5) 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].
  * (6) O. Bohigas, M.J. Giannoni and C. Schmit, _Characterization of chaotic quantum spectra and universality of level fluctuation laws_ , _Phys. Rev. Lett._ 52 (1984) 1.
  * (7) T. Ali, A. Bhattacharyya, S.S. Haque, E.H. Kim, N. Moynihan and J. Murugan, _Chaos and Complexity in Quantum Mechanics_ , _Phys. Rev. D_ 101 (2020) 026021 [1905.13534].
  * (8) A. Bhattacharyya, W. Chemissany, S. Shajidul Haque and B. Yan, _Towards the web of quantum chaos diagnostics_ , _Eur. Phys. J. C_ 82 (2022) 87 [1909.01894].
  * (9) A. Bhattacharyya, W. Chemissany, S.S. Haque, J. Murugan and B. Yan, _The Multi-faceted Inverted Harmonic Oscillator: Chaos and Complexity_ , _SciPost Phys. Core_ 4 (2021) 002 [2007.01232].
  * (10) A. Bhattacharyya, S.S. Haque and E.H. Kim, _Complexity from the reduced density matrix: a new diagnostic for chaos_ , _JHEP_ 10 (2021) 028 [2011.04705].
  * (11) V. Balasubramanian, M. DeCross, A. Kar, Y.C. Li and O. Parrikar, _Complexity growth in integrable and chaotic models_ , _JHEP_ 07 (2021) 011 [2101.02209].
  * (12) P. Caputa, J.M. Magan and D. Patramanis, _Geometry of Krylov complexity_ , _Phys. Rev. Res._ 4 (2022) 013041 [2109.03824].
  * (13) V. Balasubramanian, P. Caputa, J. Magan and Q. Wu, _Quantum chaos and the complexity of spread of states_ , 2202.06957.
  * (14) P. Caputa and S. Liu, _Quantum complexity and topological phases of matter_ , 2205.05688.
  * (15) P. Caputa, N. Gupta, S.S. Haque, S. Liu, J. Murugan and H.J.R. Van Zyl, _Spread Complexity and Topological Transitions in the Kitaev Chain_ , 2208.06311.
  * (16) T. Ali, A. Bhattacharyya, S. Shajidul Haque, E.H. Kim and N. Moynihan, _Post-Quench Evolution of Complexity and Entanglement in a Topological System_ , _Phys. Lett. B_ 811 (2020) 135919 [1811.05985].
  * (17) W.P. Su, J.R. Schrieffer and A.J. Heeger, _Solitons in polyacetylene_ , _Phys. Rev. Lett._ 42 (1979) 1698.
  * (18) A.Y. Kitaev, _Unpaired majorana fermions in quantum wires_ , _Physics-Uspekhi_ 44 (2001) 131.
  * (19) W. Mück and Y. Yang, _Krylov complexity and orthogonal polynomials_ , 2205.12815.
  * (20) B.N. Parlett and D.S. Scott, _The lanczos algorithm with selective orthogonalization_ , _Mathematics of Computation_ 33 (1979) 217\.
  * (21) H.D. Simon, _The lanczos algorithm with partial reorthogonalization_ , _Mathematics of Computation_ 42 (1984) 115.
  * (22) G. Gambolati and M. Putti, _Comparison of lanczos and optimization methods in the partial solution of sparse symmetric eigenproblems_ , _International Journal for Numerical Methods in Engineering_ 37 (1994) 605 .
  * (23) S. Berceanu, _Coherent states associated to the real Jacobi group_ , _AIP Conf. Proc._ 956 (2007) 233.
  * (24) S. Berceanu, _The jacobi group and the squeezed states—some comments_ , _arXiv: Differential Geometry_ 1191 (2009) 21.
  * (25) D. Kressner and C. Tobler, _Krylov subspace methods for linear systems with tensor product structure_ , _SIAM J. Matrix Analysis Applications_ 31 (2010) 1688.
  * (26) K. Adhikari and S. Choudhury, _Cosmological krylov complexity_ , 2022.
  * (27) C. Lanczos, _An iteration method for the solution of the eigenvalue problem of linear differential and integral operators_ , _Journal of research of the National Bureau of Standards_ 45 (1950) 255.
  * (28) V. Viswanath and G. Müller, _The recursion method: application to many-body dynamics_ , vol. 23, Springer Science & Business Media (2008).
  * (29) A.M. Perelomov, _Coherent states for arbitrary Lie group_ , _Communications in Mathematical Physics_ 26 (1972) 222 .
  * (30) R.J. Glauber, _Coherent and incoherent states of the radiation field_ , _Physical Review_ 131 (1963) 2766.
  * (31) S.S. Haque, C. Jana and B. Underwood, _Saturation of thermal complexity of purification_ , _JHEP_ 01 (2022) 159 [2107.08969].
  * (32) S.S. Haque, C. Jana and B. Underwood, _Operator Complexity for Quantum Scalar Fields and Cosmological Perturbations_ , 2110.08356.
  * (33) M.A. Nielsen, M.R. Dowling, M. Gu and A.C. Doherty, _Quantum computation as geometry_ , _Science_ 311 (2006) 1133 [quant-ph/0603161].
  * (34) M.A. Nielsen, _A geometric approach to quantum circuit lower bounds_ , 2005\.  10.48550/ARXIV.QUANT-PH/0502070.
  * (35) M.R. Dowling and M.A. Nielsen, _The geometry of quantum computation_ , 2007\.  10.48550/ARXIV.QUANT-PH/0701004.
  * (36) N. Chagnet, S. Chapman, J. de Boer and C. Zukowski, _Complexity for Conformal Field Theories in General Dimensions_ , _Phys. Rev. Lett._ 128 (2022) 051601 [2103.06920].
  * (37) R.d.M. Koch, M. Kim and H.J.R. Van Zyl, _Complexity from spinning primaries_ , _JHEP_ 12 (2021) 030 [2108.10669].
  * (38) P. Rabambi and H.J.R. Van Zyl, _Complexity for superconformal primaries from BCH techniques_ , _JHEP_ 11 (2022) 095 [2208.05520].
  * (39) 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].
  * (40) J. Kim, J. Murugan, J. Olle and D. Rosa, _Operator delocalization in quantum networks_ , _Phys. Rev. A_ 105 (2022) L010201 [2109.05301].
  * (41) E. Rabinovici, A. Sánchez-Garrido, R. Shir and J. Sonner, _Krylov localization and suppression of complexity_ , _JHEP_ 03 (2022) 211 [2112.12128].
  * (42) E. Rabinovici, A. Sánchez-Garrido, R. Shir and J. Sonner, _Krylov complexity from integrability to chaos_ , _JHEP_ 07 (2022) 151 [2207.07701].
  * (43) B. Bhattacharjee, X. Cao, P. Nandy and T. Pathak, _Krylov complexity in saddle-dominated scrambling_ , _JHEP_ 05 (2022) 174 [2203.03534].