††institutetext: 1Department of Physics, Indian Institute of Technology
Guwahati, North Guwahati-781039, India††institutetext: 2 Faculty of Physics
and Applied Computer Science,
AGH University of Krakow, al. A. Mickiewicza 30, 30-059 Krakow,
Poland††institutetext: 3Department of Physics, National Tsing-Hua University,
Hsinchu 30013, Taiwan††institutetext: 4Department of Physics, Indian Institute
of Technology Kharagpur, Kharagpur, 721302, India

# A study of three butterflies: entanglement wedge method, OTOC and pole-
skipping

Banashree Baishya1 Adrita Chakraborty2,3 Nibedita Padhi4
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

In this work, we investigate two salient chaotic features, namely Lyapunov
exponent and butterfly velocity, in the context of an asymptotically Lifshitz
black hole background with an arbitrary critical exponent. These features are
computed using three methods: entanglement wedge method, out-of-time-ordered
correlator computation and pole-skipping. We present a comparative study of
the aforementioned features where all of these methods yield exactly similar
results for the butterfly velocity and Lyapunov exponent. This establishes an
equivalence between all three methods for probing chaos in the chosen gravity
background. Furthermore, we evaluate the chaos at the classical level by
computing the eikonal phase and Lyapunov exponent from the bulk gravity. These
quantities emerge as nontrivial functions of the anisotropy index. By
examining the classical eikonal phase, we uncover different scattering
scenarios in the near-horizon and near-boundary regimes. We also discuss
potential limitations regarding the choice of the turning point of the null
geodesic in our approach.

###### Keywords:

Chaos, Lifshitz holography

##  1 Introduction

The presence of chaos is robust in various disordered many-body systems, as
evidenced by the study of black holes. Dynamical phenomena such as the growth
of local perturbations due to sufficiently small changes in the initial
conditions, the resulting deviations in the trajectory, and the scrambling of
quantum information, are some immediate consequences of chaos in disordered
systems. The Anti-de Sitter (AdS)/Conformal Field Theory (CFT) correspondence,
established by Maldacena:1997re ; Aharony:1999ti  provides one of the most
competent platforms to realise conformally invariant strongly interacting
quantum field theories through their dual AdS background containing the field
theory on its boundary. Such holographic duality can also be extended as
general gauge/gravity duality where the field theory breaks conformal
invariance and its dual geometry is eventually a non-AdS one. A class of
nonrelativistic Lifshitz field theories are such theories that follow Lifshitz
scaling symmetry

$t\rightarrow\Omega^{\xi}t,\leavevmode\nobreak\ \leavevmode\nobreak\
\vec{x}\rightarrow\Omega\vec{x},$ (1.1)

characterized by the dynamical critical exponent $\xi\in\mathbb{R}$, $\xi\geq
1$. Such QFTs break Lorentz invariance due to the presence of an anisotropic
scale factor along the temporal direction and thus do not adhere to CFT. The
holographic dual of Lifshitz field theories is the nonrelativistic Lifshitz
geometry which shares similar scale symmetry Bala ; Walter ; Juan08 ; Allan08
; Kachru . Exploration of such theories finds its significance in the context
of the correlation between high energy physics and condensed matter physics
due to the fact that these theories, especially for $\xi=2$, are directly
related to strongly correlated electron systems. A series of research for
rigorous understanding of such field theories have been going on over the last
decade by implementing various perspectives Ker2017 ; KORD ; parkChanyong ;
Cheyne:2017bis ; Basak:2023otu . One of the crucial features of such theories
is the appearance of chaos caused by the dynamic quantum critical point.
Previously the holographic Lifshitz bulk spacetime has been found to possess
chaotic behaviour of any probe circular string moving in it Bai:2014wpa . In
addition to this, the nonintegrability of the class of hyperscale violating
Lifshitz theories are studied in Giataganas:2014hma . Apart from classical
chaos, the chaotic growth of marginally irrelevant operators in 2+1D Lifshitz
field theory with $\xi=2$, also famously known as the quantum Lifshitz model,
is previously explored in Plamadeala_2018  by using a perturbative approach.
Holographic studies on scrambling properties have been done Sircar:2016old
for generic gravitational theories including the asymptotically Lifshitz black
hole which reads as the plausible gravity dual of the finite temperature
counterpart of Lifshitz field theoriesBala .

In this article, we aim to study the chaotic properties of such a theory from
the bulk perspective. There are many ways to study quantum chaos such as
Random matrices Cotler:2016fpe , Complexity Susskind:2014rva ; Brown:2015lvg ;
Magan:2018nmu , Ergodicity Pausch_2021 ; ergo  etc. But, in this work, we will
implement three well-known methods in bulk gravity. In a recent
studyDong:2022ucb , the butterfly effects of higher derivative gravity
theories are analysed by two equivalent holographic understandings, to name,
using the method of entanglement wedge reconstruction and incorporating
localized gravitational shock waves at the horizon of the corresponding dual
black hole. Our agenda in this paper is to check the equivalence of the
following three methods of computing quantum chaos for the asymptotically
Lifshitz black hole :

  * •

Entanglement wedge method Mezei:2019zyt ; Fischler:2018kwt ; Dong:2022ucb

  * •

Out of time-ordered correlation function (OTOC) Shenker:2013pqa ;
Maldacena:2015waa ; Shenker:2013yza ; Roberts:2014isa ; Shenker:2014cwa ;
Plamadeala_2018 ; Roberts:2014ifa ; Khetrapal:2022dzy ; Lin_2018 ; Gu:2016hoy
; Das:2022jrr

  * •

Pole-skipping Blake:2017ris ; Blake:2018leo ; Grozdanov:2018kkt ; Ahn:2019rnq
; Blake:2021hjj ; Sil:2020jhr ; Li:2019bgc ; Amano:2022mlu ; Abbasi:2019rhy ;
Baishya:2023mgz ; Baishya:2023ojl ; Yuan:2020fvv

  1. (A)

In AdS/CFT correspondence, entanglement wedge reconstruction entang ;
Wall:2012uf ; Czech:2012bh  depicts that all the information of the
entanglement wedge in the bulk lies in the boundary region. When we perturb
the boundary state by a local operator and let the system evolve, then the
information gets scrambled in the whole space at a late time. This information
propagates outward with a constant velocity. The whole scenario has a bulk
description. Perturbing the boundary means probing a particle close to the
asymptotic boundary, which is falling toward the black hole in the bulk. This
particle then follows a trajectory that reaches inside the extremal surface,
also known as the Ryu-Takayanagi(RT) surface rt . As the trajectory of the
particle changes, the RT surface changes its shape. At late time, this RT
surface reaches up to the near-horizon region of the black hole with a
constant velocity – butterfly velocity ($v_{B}$). This method is very well-
known for extracting the butterfly velocity Mezei:2016wfz . In higher-
derivative gravity theories also, this method is well understood Dong:2022ucb
. In our work, we employ this method to calculate the butterfly velocity in
the Einstein-Proca-type Lifshitz gravity theory. The action of such theories
contains a massive vector field that breaks Lorentz symmetryTaylor:2015glc .
It is very motivating to study the chaotic properties of such a theory.

  2. (B)

Quantum chaos can further be studied by calculating the double commutator of
any two generic operators $V$ and $W$ as,

$\displaystyle C(t)$ $\displaystyle=\langle
V^{\dagger}(0)W^{\dagger}(t)W(t)V(0)+V^{\dagger}(0)W^{\dagger}(t)W(t)V(0)$
$\displaystyle\quad\quad\quad-\underbrace{W^{\dagger}(t)V^{\dagger}(0)W(t)V(0)-V^{\dagger}(0)W^{\dagger}(t)V(0)W(t)\rangle_{\beta}}_{\text{out-
of time-ordered}}\,.$ (1.2)

The out-of-time-ordered correlators (OTOC) show a large value until the
"scrambling time" Sekino:2008he ; Hayden:2007cs  and decrease rapidly
thereafter. The decrease of OTOC leads to an increase of $C(t)$ and signifies
the butterfly effect. OTOC has been rigorously studied in connection to black
holes and quantum chaos Shenker:2013pqa ; Maldacena:2015waa ; Shenker:2013yza
; Roberts:2014isa ; Shenker:2014cwa . To mention the field theory perspective,
OTOC has been computed in 2D CFT Roberts:2014ifa ; Khetrapal:2022dzy , Ising
chain Lin_2018 , rational CFT Gu:2016hoy , driven CFT Das:2022jrr  and
numerically in the quantum Lifshitz model Plamadeala_2018 . Study of the
commutator in holographic systems Shenker:2013pqa ; Shenker:2014cwa
generalizes the expression ((B)) to,

$C(t,x)\approx\exp{\lambda_{L}\left(t-t_{*}-\frac{|x|}{v_{B}}\right)}\,,$
(1.3)

where the operators are separated by a spatial distance $x$. The exponent
$\lambda_{L}$ is called the Lyapunov exponent, $v_{B}$ is the butterfly
velocity and $t_{*}$ is the scrambling time. The scrambling time is defined as
the time at which the commutator square $C(t)$ in ((B)) becomes of order
$\mathcal{O}(1)$. In two-derivative gravity, this Lyapunov exponent and
butterfly velocity is bounded as Maldacena:2015waa ; Mezei:2016wfz ,

$\lambda_{L}\leq\frac{2\pi}{\beta}\,,\hskip 56.9055ptv_{B}\leq
v_{B}^{\text{sch}}=\sqrt{\frac{d}{2(d-1)}}\,,$ (1.4)

where $\beta$ is the inverse temperature and $v_{B}^{\text{sch}}$ is the
butterfly velocity for a $d+1$ \- dimensional AdS-Schwarzschild black brane. A
natural question one can ask is whether a theory with broken Lorentz
invariance can lead to a violation of this Lyapunov exponent and butterfly
velocity bound (1.4). In this paper, we have addressed this question.

  3. (C)

The chaotic nature we discussed in the previous two methods can be captured in
the properties of a two-point energy density correlation function at the pole-
skipping Blake:2017ris ; Blake:2018leo ; Grozdanov:2018kkt ; Ahn:2019rnq
point and those points are holographically related to the near-horizon
properties of the metric perturbation in bulk. Pole-skipping (P-S) occurs when
lines of poles and zeros of retarded Green’s function intersect, i.e., a
would-be pole gets skipped! At these specific values of the frequency and
momentum (P-S point), Green’s function is not uniquely defined Blake:2019otz ;
Natsuume:2019xcy . These specific values of the frequency ($\omega_{*}$) and
momentum ($k_{*}$) in the energy-density two-point function are connected to
the measure of chaos Blake:2017ris ; Blake:2018leo  through the following
relations:

$\omega_{*}=i\lambda_{L},\hskip 56.9055ptk_{*}=\frac{i\lambda_{L}}{v_{B}},$
(1.5)

where $\lambda_{L}$ and $v_{B}$ are the Lyapunov exponent and the butterfly
velocity respectively. In this paper, we have checked whether pole-skipping is
a good diagnostic in a theory with broken Lorentz invariance or not. A
plethora of studies have so far explored the chaotic nature of the holographic
systems utilizing this P-S method in the black hole background
Grozdanov:2017ajz ; Blake:2018leo ; Ahn:2019rnq ; Wu:2019esr , plasma physics
Sil:2020jhr ; Li:2019bgc ; Amano:2022mlu , conformal field theories
Ramirez:2020qer ; Haehl:2019eae ; Das:2019tga , holographic system with chiral
anomaly Abbasi:2019rhy ; Abbasi:2023myj , with stringy corrections
Grozdanov:2018kkt , higher derivative corrections Natsuume:2019vcv ;
Wu:2019esr ; Baishya:2023mgz , brane set-up Baishya:2023ojl  etc. This
phenomenon has been extensively studied in various contexts and explored in
different directions in Ceplak:2021efc ; Kim:2020url ; Grozdanov:2020koi ;
Jeong:2021zhz ; Wang:2022mcq ; Wang:2022xoc ; Yuan:2023tft ; Natsuume:2023lzy
; Baishya:2023xbj .

We calculate the butterfly effect in asymptotically Lifshitz black holes using
three distinct methods: the entanglement wedge, out-of-time-ordered correlator
(OTOC) computation, and pole-skipping. Our analysis compares the results
obtained from these methods, demonstrating their equivalence. In
Sircar:2016old , the authors determined the scrambling time in Kruskal
coordinates for the Lifshitz black hole. However, the explicit forms of the
Lyapunov exponent and butterfly velocity for black holes with arbitrary
anisotropy have yet to be fully explored in the context of OTOCs. Our primary
objective is to calculate these two parameters using all the aforementioned
methods, thereby providing a comprehensive understanding of the chaotic
properties of Lifshitz black holes.

In addition to the comparative analysis of the quantum chaos, we also
elaborate on the eikonal phase in the heavy-heavy-light-light gravitational
scattering and the Lyapunov exponent, from the classical perspective. The
eikonal phase is priorly known to be related to the deflection angle of the
null geodesic of a gravity background Kabat:1992tb ; DAppollonio:2010krb ;
Bjerrum-Bohr:2014zsa . It also finds its own significance in the AdS/CFT
holography, see Cornalba:2006xk ; Cornalba:2006xm ; Brower:2007qh  for review.
To compute for eikonal phase shift, the null geodesic equation plays a
dominant role in determining the impact parameter, often taken as the ratio of
conserved momenta to conserved energy, which eventually affects the nature of
the eikonal phase. For real and imaginary eikonal phases, one must obtain
elastic and completely nonelastic gravitational scattering respectively. Here
in our article, we wish to see how the anisotropy of our chosen gravity
background causes different kinds of eikonal phase and thus the results in
different scatterings when we shift the extrema of the null geodesic from near
boundary to near horizon limit. We verify the eikonal phase using the WKB
approximation of the equation of motion of a scalar field. Specifically, we
examine for a probable matching of the eikonal phase yielded from both the
methods. Furthermore, we wish to briefly exhibit the classical nature of the
Lyapunov exponent in the anisotropic background and its consistency with the
restrictions that we must impose on the choices of the extrema for the
validity of our study.

We organize the paper as follows. In section 2, we give a brief idea of the
Lifshitz black hole. In section 3, we perform the explicit calculation of
Lyapunov exponent and butterfly velocity using all three methods as explained
above. Section 4 is devoted to the computation of the eikonal phase and
Lyapunov exponent via a classical approach. We utilize the notion of geodesic
stability to accomplish our purpose here. In section 5, we conclude with an
elaborate discussion of the implications of our results for both classical and
quantum chaos as well as some future scopes.

##  2 Asymptotically Lifshitz black hole

In this section, we briefly revisit the black hole geometries that asymptote
the nonrelativistic planar Lifshitz background. Such black holes were first
developed in Balasubramanian:2009rx  for $d=2,\leavevmode\nobreak\ \xi=2$.
These black hole solutions are found to consistently satisfy the equations of
motion of the Einstein-Proca type gravity theories with the inclusion of a
massive vector field Taylor_2016 ,

$\begin{gathered}\mathcal{S}=\frac{1}{2\kappa^{2}}\int
d^{d+1}x\sqrt{-g}\left[\mathcal{R}-2\Lambda-\frac{1}{4}F_{\mu\nu}F^{\mu\nu}+\frac{1}{2}m^{2}A_{\mu}A^{\mu}\right]\,,\\\
F_{\mu\nu}=\partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu}\,.\end{gathered}$ (2.1)

Here, for any (d+1) dimensional Lifshitz black hole, the cosmological constant
$\Lambda$ must be a negative quantity that depends on the anisotropy index
$\xi$ as $\Lambda=\frac{d-1}{2(d+1)}\mathcal{R}(\xi)$. Again, $F_{\mu\nu}$ is
the electromagnetic field strength. The appearance of a massive vector field
with mass $m$ in the above action causes the breaking of the Lorentz
invariance in the corresponding metric. The generic form of the $(d+1)$
dimensional Lifshitz black hole solution of the above action can be written
as,

$\begin{gathered}ds^{2}=\frac{R^{2}}{z^{2}}\left[-\frac{R^{2(\xi-1)}}{z^{2(\xi-1)}}f(z)dt^{2}+d\vec{x}_{d-1}^{2}+\frac{dz^{2}}{f(z)}\right]\,,\\\
f(z)=1-\left(\frac{z}{z_{h}}\right)^{d-1+\xi}\,.\end{gathered}$ (2.2)

$R$ is the AdS radius, $z$ is the radial coordinate and $z_{h}$ is the black
hole horizon. This geometry asymptotically reaches the planar Lifshitz
spacetime

$ds^{2}=\frac{R^{2}}{z^{2}}\left[-\frac{R^{2(\xi-1)}}{z^{2(\xi-1)}}dt^{2}+d\vec{x}_{d-1}^{2}+dz^{2}\right]\,.$
(2.3)

Similarly as the planar Lifshitz spacetime, the Lifshitz black hole geometry
also follows anisotropic scaling along the time and space directions

$t\rightarrow\Omega^{\xi}t,\leavevmode\nobreak\ \leavevmode\nobreak\
z\rightarrow\Omega z,\leavevmode\nobreak\ \leavevmode\nobreak\
\vec{x}\rightarrow\Omega\vec{x}\,,$ (2.4)

that breaks the Lorentz invariance. Equation (2.2) describes a one-parameter
family of linearly charged Lifshitz black holes that are thermodynamically
stable and become extremal at the vanishing size Danielsson:2009gi . The
thermodynamics of such an asymptotically Lifshitz black hole can be reproduced
from the holographic renormalization of the gravity theory represented by the
action (2.1) parkChanyong . The legitimate dual field theory that lives on the
boundary of Lifshitz black hole geometry is the finite temperature version of
the non-relativistic Lifshitz field theories. When we take 2D thermal CFT, the
flat boundary that accommodates the CFT is compactified along the time
direction and becomes a cylinder. This compactification is done by taking
$t\rightarrow t+\beta$, where $\beta$ gives the circumference of the
compactifying circle. Thus $\beta$ acquires the dimension of length. Now, let
us consider a similar compactification of the boundary in the case of thermal
LFT for which the dimension of time is [length]ξ due to Lifshitz scaling. Thus
if we take $t\rightarrow t+\tilde{\beta}$ with the analogy of thermal CFT for
the compactification, then $\tilde{\beta}$ should have the dimension of time.
Thus we can assume the inverse temperature of thermal LFT as
$\tilde{\beta}=\beta^{\xi}$, where $\beta$ is the inverse temperature of the
thermal CFT. Hence, in the dual bulk side, the temperature of the Lifshitz
black hole will also be $T=\frac{1}{\beta^{\xi}}$. The temperature of the
black hole is given by

$T=\frac{1}{z_{h}^{\xi}}\frac{d-1+\xi}{4\pi}=\frac{1}{\beta^{\xi}}\equiv\frac{1}{\tilde{\beta}}\,.$
(2.5)

Note that, $\tilde{\beta}$ has dimension of length. For maintaining simplicity
in our further calculations, we will denote $\tilde{\beta}$ as $\beta$. In the
next section, we start with the $(d+1)$-dimensional Lifshitz black hole
geometry given in (2.2) to analyse the chaotic features via a quantum
approach.

##  3 Analyses of quantum chaos

In this section, we will show the explicit calculations of the three methods
that we are interested in to study the chaotic parameters. Our focus is to
compare these results and see the similar dependence of the chaotic parameters
on the anisotropy.

###  3.1 Entanglement wedge method

As discussed in the introduction, we aim in this subsection to find the
butterfly velocity with the entanglement wedge method. For this, we need to
calculate the size of the smallest boundary region whose entanglement wedge
encloses the infalling particle as shown in Figure 1. In this pictorial
description, we can see that $z=0$ is the boundary and $z=1$ is the horizon of
the black hole. A particle (shown by the red dot) is on the RT surface (shown
by the light teal portions) bounded by the boundary region. In general, the RT
surface location is determined by extremizing the holographic entanglement
entropy functional. We can write the entropy functional as,

$S_{\text{EE}}=2\pi\int\differential^{d-1}y\sqrt{\gamma}\,,$ (3.1)

where $\gamma$ is the determinant of the induced metric and $y$ is the set of
coordinates on an appropriate codimension-2 surface.

Figure 1: A particle (red dot) is enclosed in the RT surface (light teal)

Now, we will use background (2.2) for this wedge method. At a constant $t$
hypersurface ($\differential t=0$) and parameterizing $z$ as $z(r)$ with
assumption $r=|x^{i}|$, we can write the induced metric as,

$\displaystyle\gamma_{\alpha\beta}\,\differential x^{\alpha}\differential
x^{\beta}$ $\displaystyle=\frac{R^{2}}{z^{2}}\left[\left(1+\frac{z^{\prime
2}}{f(z)}\right)\differential
r^{2}+r^{2}\differential\Omega_{d-2}^{2}\right].$ (3.2)

Note that in the (2.2), the near-horizon limit is achieved by $z=1$ limit and
the boundary limit is achieved by $z=0$ limit. As we will work with the near-
horizon region ($z=1$), we define the RT surface in the $z=1$ limit by,

$z(r)=1-\epsilon\,y(r)^{2},$ (3.3)

where we parameterize the near-horizon limit($z=1$) with a new function
$y(r)$. The coefficient $\epsilon\,(>0)$ is a very small number and the
function $y(r)$ is the RT profile, which we will solve from the RT equation by
performing a Taylor series expansion around $\epsilon=0$. Now, expanding the
induced metric up to order $\epsilon$, we get,

$\displaystyle\gamma_{\alpha\beta}\,\differential x^{\alpha}\differential
x^{\beta}$ $\displaystyle=R^{2}\left[\left(\frac{4\,\epsilon^{2}y^{2}y^{\prime
2}}{(1-\epsilon\,y^{2})^{2}f_{1}(1-z)}+\frac{1}{(1-\epsilon\,y^{2})^{2}}\right)\differential
r^{2}+\frac{r^{2}}{(1-\epsilon\,y^{2})^{2}}\differential\Omega_{d-2}^{2}\right]\,,$
$\displaystyle=R^{2}\left[\left(1+2\epsilon(y^{2}+\frac{2}{f_{1}}y^{\prime
2})\right)\differential r^{2}+r^{2}(1+2\epsilon
y^{2})\differential\Omega_{d-2}^{2}+\mathcal{O}(\epsilon^{2})\right],$ (3.4)

where we Taylor expand the factor $f(z)$ near the horizon up to the first
order as $f(z)\approx f_{1}(1-z)$. The determinant of the induced metric is
given by,

$\displaystyle\sqrt{\gamma}=R^{2}r^{d-2}\left[1+\epsilon\left((d-1)y(r)^{2}+\frac{2}{f_{1}}y^{\prime}(r)^{2}\right)\right]+\mathcal{O}(\epsilon^{2})\,.$
(3.5)

Now, we can write the RT equation by varying the above equation with respect
to $y(r)$. Keeping terms only up to order $\epsilon$, we get,

$(d-1)y(r)-\frac{2}{f_{1}}\left(y^{\prime\prime}(r)+(d-2)\frac{y^{\prime}(r)}{r}\right)=0\,.$
(3.6)

Solving this second-order differential equation, we get,

$y(r)=r^{n}\left[\text{J}_{n}(\mu r)+\text{Y}_{n}(\mu r)\right]\,,$ (3.7)

where $\text{J}_{n}$ and $\text{Y}_{n}$ are the Bessel functions of first kind
and second kind respectively. In the above equation,

$n=\frac{1}{2}(d-3),\hskip 56.9055pt\mu=\sqrt{\frac{(d-1)(d-1+\xi)}{2}}\,.$
(3.8)

Here, the scaling parameter $\mu$ is behaving as the momentum. From the eq.
(3.3), we can see that near the horizon $z=1$, the function $y(r)$ takes the
form of $0$. So, we will see the form of these Bessel functions near $r=0$. At
this vicinity, both these Bessel functions will behave as $\sim\mu^{n}$.
Following Dong:2022ucb , where the authors have calculated the butterfly
velocity with an exponential ansatz; here in this paper, we have prescribed a
more general result for calculating butterfly velocity without taking any
ansatz.

The surface enclosing the particle is approaching the horizon at a constant
speed, termed as butterfly velocity ($v_{B}$). We assume that at each point in
time, the tip of the RT surface touches the particle which is visualised by
demanding $y(r=0,t)\sim e^{-\frac{2\pi}{\beta}t}$. Therefore, with the value
of $\mu$, we can calculate the butterfly velocity as,

$\displaystyle v_{B}$
$\displaystyle=\frac{2\pi}{\beta\mu}=\frac{2\sqrt{2}\pi}{\beta\sqrt{(d-1)(d-1+\xi)}}\,.$
(3.9)

So, the butterfly velocity $v_{B}$ becomes (at $z_{h}=1$),

$v_{B}=\sqrt{\frac{d-1+\xi}{2(d-1)}}\,.$ (3.10)

This is the butterfly velocity for a $d+1$\- dimensional Lifshitz black hole.
In $\xi\rightarrow 1$ limit, we are getting the planar black hole result
$v_{B}=\sqrt{\frac{d}{2(d-1)}}$. So, the anisotropy affects the butterfly
velocity nontrivially. Butterfly velocity characterizes the propagation of
chaos in a local system. It has been observed previously that the chaotic
feature of the Lifshitz invariant system increases for large $\xi$ values, in
other words, for large anisotropy. From the expression of $v_{B}$ that we get,
it is obvious that $v_{B}$ increases monotonically with $\xi$ for any fixed
dimension of the bulk. Thus, the chaotic feature of the asymptotic Lifshitz
black hole monotonically increases with larger anisotropy which is consistent
from the previous results.

In the next subsection, we will demonstrate the appearance of similar
butterfly velocity in the context of out-of-time-ordered correlators.

###  3.2 Out-of-Time-Ordered correlators

In this subsection, we will study the out-of-time-ordered correlator (OTOC)
for thermal Lifshitz field theory dual to a linearly charged asymptotically
Lifshitz black hole. It is convenient to study the OTOC in the Kruskal-
Szekeres form of the bulk metric, which smoothly covers the globally extended
space-time. In Kruskal geometry, the OTOC can be interpreted as the two-
particle gravitational scattering between particles moving along the two
horizons Shenker:2014cwa ; Maldacena:2015waa ; Shenker:2013yza ;
Roberts:2014isa . We shall compute the OTOC in the specific dual-field theory
in the form given as

$\langle\hat{W}(t_{2},x_{2})\hat{V}(t_{1},x_{1})\hat{W}(t_{2},x_{2})\hat{V}(t_{1},x_{1})\rangle_{\beta}\,,$
(3.11)

where we consider $t_{2}-t_{1}>>\beta$. The particles that participate in the
gravitational scattering process are the operators $\hat{W}$ and $\hat{V}$,
which are inserted in the dual thermofield double (TFD) state. Though the
quantity (3.11) is one-sided, our whole computation of the OTOC has been done
in a two-sided geometry. So far, the OTOC via the gravitational scattering
process has been explored with various black hole backgrounds Shenker:2014cwa
; Maldacena:2015waa ; Shenker:2013yza ; Roberts:2014isa ; Ahn:2019rnq ;
Blake:2021hjj ; Fischler:2018kwt ; Jahnke:2019gxr ; Mezei:2019dfv ;
Reynolds:2016pmi ; Poojary:2018esz ; Banerjee:2019vff ; Saha:2024bpt . In our
work, we have calculated this quantity (3.38) for a linearly charged
asymptotically Lifshitz black hole with arbitrary anisotropy. In the Lifshitz
background, the authors of Sircar:2016old  have calculated the scrambling
time, though the background they worked with is slightly different from ours.
Here, we have provided an explicit calculation of the Lyapunov exponent and
butterfly velocity in the two-sided geometry, which yields interesting
results.

####  3.2.1 Kruskal extension

To study the out-of-time-ordered correlators, let us first construct the
Kruskal extension of the (2+1) dimensional Lifshitz black hole. We will start
with the metric (2.2), but replacing $z=1/r$ and considering the AdS radius to
be 1 for simplicity,

$\differential s^{2}=-r^{2\xi}f(r)\differential
t^{2}+\frac{1}{r^{2}}\frac{\differential
r^{2}}{f(r)}+r^{2}\differential{x}^{2}\,,$ (3.12)

and the blackening factor $f(r)$ is given by,

$f(r)=1-\left(\frac{r_{h}}{r}\right)^{1+\xi}\,.$ (3.13)

$r$ is the radial coordinate and $r_{h}$ is the black hole horizon in the
emblackening factor $f(r)$. In Alvarez:2014pra ; Zhu_2020 , nonlinearly
charged Lifshitz black holes with one and two horizons are obtained with the
Eddington-Finkelstein form of the metric. In this case, nonlinearity appears
in the blackening factors, and it is characterized by using some extra
parameters along with the dynamical exponent $\xi$ of the Lifshitz black hole.
Single horizons and double horizons for such black holes are obtained by using
specific parametric conditions. However, for the metric in (3.12), we will
construct an explicit form of the linearly charged thermal Lifshitz black hole
with double horizons by using the Kruskal-Szekeres coordinates. To work in
these coordinates, we transform,

$u=e^{\frac{\alpha}{2}(r_{*}-t)},\leavevmode\nobreak\ \leavevmode\nobreak\
v=e^{\frac{\alpha}{2}(r_{*}+t)}\,,$ (3.14)

where $\alpha$ is a function of $r_{h}$ that relates with the finite
temperature $T$ of the black hole as $\alpha\equiv 2\pi T$. The tortoise
coordinate $r_{*}$ will play a very crucial role in defining the smooth
horizon in the ($u,v$) coordinate. We define,

$r_{*}(r)=\int\frac{\differential r}{{r}^{1+\xi}f(r)}=\int\frac{\differential
r}{{r}^{1+\xi}\left[1-\left(\frac{r_{h}}{r}\right)^{1+\xi}\right]}\,.$ (3.15)

Integrating equation (3.15), we get,

$r_{*}(r)\approx\frac{r^{-\xi}}{\xi(1+\xi)}\,{{}_{2}}F_{1}\left(1,-\xi,1-\xi,\frac{r}{r_{h}}\right)\,.$
(3.16)

We consider $f(r)=(r-r_{h})f^{\prime}(r_{h})$ while solving the above
equation. Now, taking the near-horizon limit ($r\rightarrow r_{h}$) of the
whole solution, we get,

$r_{*}\approx\frac{1}{r_{h}^{\xi}(1+\xi)}\log\left(\frac{r-r_{h}}{r_{h}}\right)\,.$
(3.17)

As there is a relation between ($u,v$) coordinate and $r_{*}$ (from (3.14)),
we can easily write,

$uv=e^{\alpha r_{*}}\approx\left(\frac{r-r_{h}}{r_{h}}\right)\,.$ (3.18)

This equation allows us to claim that the $uv=0$ is the same limit as the
$r=r_{h}$. This result will be quite helpful for the analysis afterwards.
Writing the metric (3.12) in terms of the null coordinates $u$ and $v$, we get

$\differential
s^{2}=\frac{4r(u,v)^{2\xi}}{\alpha^{2}uv}\,f(r(u,v))\differential
u\,\differential v+r(u,v)^{2}\,\differential x^{2}=2A(u,v)\differential
u\,\differential v+B(u,v)\differential x^{2}\,.$ (3.19)

In this Kruskal geometry, $uv=0$ is defined as the horizon, singularity is at
$uv=1$ and both left/right boundaries are located at $uv=-1$. This unperturbed
background will obey Einstein’s equation which is of the form,

$E_{\mu\nu}=\kappa T_{\mu\nu}\,.$ (3.20)

Here, $E_{\mu\nu}$ is the Einstein tensor, $\kappa=8\pi G_{N}$ is a constant
related to the Newton’s constant $G_{N}$ and $T_{\mu\nu}$ is the stress-energy
tensor which is of the form,

$T_{0}^{\text{matter}}=2T_{uv}\differential u\differential
v+T_{uu}\differential u^{2}+T_{vv}\differential v^{2}+T_{xx}\differential
x^{2}\,.$ (3.21)

It is worth mentioning that this stress tensor accounts for the cosmological
constant and is consistent with the Ricci tensor of the unperturbed
background. For the unperturbed Einstein equation, solving $T_{uu},T_{vv}$ and
$T_{uv}$, we get,

$\displaystyle T_{uu}$
$\displaystyle=\frac{1}{4A(u,v)B(u,v)^{2}}\left[A(u,v)\partial_{u}B(u,v)^{2}+2B(u,v)\left(\partial_{u}A(u,v)\partial_{u}B(u,v)\right.\right.$
$\displaystyle\left.\left.\hskip
170.71652pt-A(u,v)\partial_{u}^{2}B(u,v)\right)\right]\,,$ (3.22a)
$\displaystyle T_{uv}$
$\displaystyle=\frac{1}{4B(u,v)^{2}}\left[2B(u,v)\partial_{u}\partial_{v}B(u,v)-\partial_{v}B(u,v)\partial_{u}B(u,v)\right]\,,$
(3.22b) $\displaystyle T_{vv}$
$\displaystyle=\frac{1}{4A(u,v)B(u,v)^{2}}\left[A(u,v)\partial_{v}B(u,v)^{2}+2B(u,v)\left(\partial_{v}A(u,v)\partial_{v}B(u,v)\right.\right.$
$\displaystyle\left.\left.\hskip
170.71652pt-A(u,v)\partial_{v}^{2}B(u,v)\right)\right]\,.$ (3.22c)

####  3.2.2 Gravitational backreaction due to shock wave

In the previous subsection, we constructed a two-sided Lifshitz black hole
geometry in a (2+1) - dimensional background. This section aims to study the
back-reacted metric due to a null pulse of energy localized along the $v=0$
horizon. In the bulk-boundary picture, this means we have inserted an operator
in the boundary thermal state at some past time $t$ and let it evolve.
Inserting an operator in the boundary leads to releasing a perturbation close
to the boundary, which falls into the bulk. From the point of view of the
$t=0$ slice, the energy of this perturbation is exponentially increasing with
time $t$ Shenker:2013pqa . Thus, we can approximate this perturbation as a
null pulse of energy $E$ localized at the $v=0$ horizon along the
$u$-direction. We must consider the back-reaction of this null pulse in the
two-sided geometry.

$u=0$$v=0$$u$$v$ Figure 2: Left: Penrose diagram of two-sided geometry. Right:
A shock wave (red zigzag line) localised along $v=0$ causes a discontinuity of
the $u=0$ horizon (the blue arrow)

We consider the following form of the stress-tensor localised at $v=0$,

$T_{vv}^{\text{shock}}=E\exp(\frac{2\pi t}{\beta})\delta(v)\delta(x)\,.$
(3.23)

The shock wave localised at $v=0$ horizon is indeed splitting the causal past
(white) and future (green) of $v$. The metric in the green region will change
while in the white region, it will remain the same. In the Penrose diagram, we
can clearly see the shift along the $v$-direction. We can write the form of
the shift in the following way:

$\displaystyle\bar{u}=u+\Theta(v)\eta(u,v,x)\rightarrow\differential{u}=\differential\bar{u}-\eta(u,v,x)\delta(v)\differential
v\,,$ (3.24) $\displaystyle\bar{v}=v\,,\hskip 56.9055pt\bar{x}=x\,.$ (3.25)

Our goal is to calculate the form of the function $\eta(u,v,x)$ – the
shockwave profile. Note that the Heaviside step function $\Theta(v)$ is
present to ensure that only the causal future of the pulse is affected by its
presence. In terms of these new coordinates ($\bar{u},\bar{v},\bar{x}$), we
can write the metric as,

$\differential
s^{2}=2A(\bar{u},\bar{v})\differential\bar{u}\differential\bar{v}+B(\bar{u},\bar{v})\differential\bar{x}^{2}-2A(\bar{u},\bar{v})\eta(u,v,x)\delta(\bar{v})\differential\bar{v}^{2}\,.$
(3.26)

With this back-reacted form of the metric, we need to calculate a form of
$\eta(u,v,x)$ which will obey eq.(3.20). In the above metric,

$\displaystyle A(\bar{u},\bar{v})$
$\displaystyle=\frac{2r(\bar{u},\bar{v})^{2\xi}}{\alpha^{2}\bar{u}\bar{v}}f(r(\bar{u},\bar{v})),\hskip
56.9055ptB(\bar{u},\bar{v})=r(\bar{u},\bar{v})^{2}\,,$ (3.27)

where we have expressed the radial coordinate $r$ in terms of ($u,v$). Now,
for the perturbed background, we wish to calculate the stress tensor. The
stress tensor for the matter part reads as

$\begin{split}T_{\mu\nu}^{\text{matter}}=2\left[T_{\bar{u}\bar{v}}-2T_{\bar{u}\bar{u}}\eta(\bar{u},\bar{v},x)\delta(\bar{v}\right]\differential\bar{u}\differential\bar{v}+\left[T_{\bar{v}\bar{v}}+T_{\bar{u}\bar{u}}\eta(\bar{u},\bar{v},x)^{2}\delta(\bar{v})^{2}\right.\\\
\left.-2T_{\bar{u}\bar{v}}\eta(\bar{u},\bar{v},x)\delta(\bar{v})\right]\differential\bar{v}^{2}+T_{\bar{u}\bar{u}}\differential\bar{u}^{2}+T_{\bar{x}\bar{x}}\differential\bar{x}^{2}\,.\end{split}$
(3.28)

We will drop the bar convention for further analyses so that the result looks
simpler. Now, we will solve the back-reacted Einstein’s equation which is of
the form :

$E_{\mu\nu}=\kappa\left(T_{\mu\nu}^{\text{matter}}+T_{\mu\nu}^{\text{shock}}\right)\,.$
(3.29)

The expressions of $T_{\mu\nu}^{\text{matter}}$ and
$T_{\mu\nu}^{\text{shock}}$ is given by (3.28) and (3.23) respectively.
Solving the $vv$ component of the perturbed Einstein equation, we get,

$\displaystyle\left[\partial_{x}^{2}-\frac{1}{A(u,v)}\partial_{u}\partial_{v}B(u,v)\right]\eta(u,v,x)-\frac{8\pi
G_{N}E\,B(u,v)}{A(u,v)}e^{2\pi t/\beta}\delta(v)\delta(x)=0\,,$
$\displaystyle\Rightarrow\left[\partial_{x}^{2}-\mathcal{M}(u,v)^{2}\right]\eta(u,v,x)=e^{\frac{2\pi}{\beta}(t-t_{*})}\delta(v)\delta(x)\,.$
(3.30)

Where,

$\mathcal{M}(u,v)=\sqrt{\frac{\partial_{u}\partial_{v}B(u,v)}{A(u,v)}}\,,\hskip
56.9055ptt_{*}(u,v)=\frac{\beta}{2\pi}\log\left(\frac{A(u,v)}{8\pi
G_{N}E\,B(u,v)}\right)\,.$ (3.31)

The function $t_{*}(u,v)$ has some physical meaning, which we will discuss in
the upcoming part. After implementing the values of $A(u,v)$ and $B(u,v)$ to
$t_{*}(u,v)$, it will take a constant form. We can see that the shockwave
profile $\eta(u,v,x)$ is obeying an ODE which we can solve easily. In writing
the above equation, we assume, $T_{uu}=0$ near the $v=0$ horizon
Fischler:2018kwt . Our task is to solve the equation (3.30) near the horizon.
Near the horizon ($v=0$ or $r=r_{h}$),

$A(u,v)|_{v=0}=\frac{8}{1+\xi},\hskip
56.9055pt\partial_{u}\partial_{v}B(u,v)|_{v=0}=4r_{h}^{2}\,.$ (3.32)

In taking the derivative of $B(u,v)$ with $u$ and $v$, we have used the
relation (3.18) and performed a chain rule differentiation with $r_{*}$. Now,
implementing the above two expressions into (3.31), we get,

$\mathcal{M}(u,v)|_{v=0}=\sqrt{\frac{1+\xi}{2}}\,,\hskip
56.9055ptt_{*}=\frac{\beta}{2\pi}\log\left(\frac{1}{\pi
G_{N}E(1+\xi)}\right)\,,$ (3.33)

where we have considered $r_{h}=1$. The constant $t_{*}$ is called the
scrambling time, which is dependent on the anisotropy index $\xi$. Near the
horizon, the equation (3.30) becomes,

$\left[\partial_{x}^{2}-\frac{1+\xi}{2}\right]\eta(x)=e^{\frac{2\pi}{\beta}(t-t_{*})}\delta(x)\,.$
(3.34)

To solve this, first, we will solve the homogeneous equation and get,

$\eta(x)=\begin{cases}c_{1}e^{\mathcal{M}x}+c_{2}e^{-\mathcal{M}x}\quad\text{for}\leavevmode\nobreak\
\leavevmode\nobreak\ x>0\,,\\\
c_{3}e^{\mathcal{M}x}+c_{4}e^{-\mathcal{M}x}\quad\text{for}\leavevmode\nobreak\
\leavevmode\nobreak\ x<0\,.\end{cases}$ (3.35)

Now, to solve the non-homogeneous equation, we will check the discontinuity of
the first derivative of the solution and get,

$\eta^{\prime}(\epsilon)-\eta^{\prime}(-\epsilon)=e^{\frac{2\pi}{\beta}(t-t_{*})}\,,$
(3.36)

where, $\epsilon$ is a small parameter $\epsilon\rightarrow 0$. Imposing the
continuity of the solution at $x=0$ and with the discontinuity equation
(3.36), we get a relation between the coefficients as,

$c_{1}-c_{3}=c_{4}-c_{2}=\frac{e^{\frac{2\pi}{\beta}(t-t_{*})}}{\sqrt{2(1+\xi)}}\,.$
(3.37)

We will choose $c_{2}=c_{3}=0$ and the final solution for the shock will be,

$\eta(t,x)=\frac{1}{\sqrt{2(1+\xi)}}e^{\frac{2\pi}{\beta}(t-t_{*}-\frac{\beta\mathcal{M}}{2\pi}x)}\,.$
(3.38)

This shockwave profile will completely determine the four-point out-of-time
ordered function - OTOC. We can identify the commutator $C(t,x)$ in equation
(1.3) with the shockwave profile $\eta(t,x)$ and write,

$\lambda_{L}=2\pi T\,,\hskip 56.9055ptv_{B}=\sqrt{\frac{1+\xi}{2}}\,.$ (3.39)

This Lyapunov exponent agrees with our results of the entanglement wedge
method. For the butterfly velocity, this result is found to be the same as
obtained from the entanglement wedge result for $d=2$ in (3.10). Thus, we can
claim that these two methods are congruent for a reliable study of chaos in
Lifshitz background.

###  3.3 Pole-skipping

In this subsection, we wish to study the retarded energy density correlation
function $G^{\text{R}}_{T_{00}T_{00}}(\omega,k)$ near the pole-skipping point
(1.5). This correlation function can be calculated by perturbing the
gravitational background (2.2) and imposing ingoing boundary conditions at the
black hole horizon. To perform the pole-skipping analysis, it is customary to
work with Eddington-Finkelstein coordinates, defined by

$v=t+r_{*},\hskip 56.9055pt\frac{\differential r_{*}}{\differential
r}=\frac{1}{r^{1+\xi}f(r)}\,,$ (3.40)

where $v$ is the null coordinate and $r_{*}$ is the tortoise coordinate. With
the above coordinate system, the background metric (2.2) is transformed into

$\differential s^{2}=-r^{2\xi}f(r)\differential v^{2}+2r^{\xi-1}\differential
v\differential r+r^{2}\differential x^{2}\,.$ (3.41)

For our current analysis, we choose the gravity background as the (2+1)-D
Lifshitz black hole. Now, to calculate the energy density correlation
function, we must perturb the background metric - specifically the sound modes
(longitudinal modes). We perturb the background metric as,

$g_{\mu\nu}\rightarrow g_{\mu\nu}+\delta g_{\mu\nu}(r)e^{-i\omega v+ikx}\,.$
(3.42)

Here we perform a Fourier transformation to the perturbed mode and choose the
wavenumber $k$ along the $x$-direction. Then, the sound modes will be–

$\delta g_{vv},\delta g_{vx},\delta g_{vr},\delta g_{rr},\delta g_{rx},\delta
g_{xx}\,.$ (3.43)

As the wave points along the $x$\- direction, sound modes are the modes that
are only along the $(v,r,x)$ directions. But, imposing radial gauge condition
$\delta g_{r\mu}=0$ and traceless ness condition $g^{\mu\nu}\delta
g_{\mu\nu}=0$, we get some redundant modes and left with–

$\delta g_{vv},\delta g_{vx}\,.$ (3.44)

The retarded energy density Green’s function is governed by the perturbed
equations that are regular at the horizon in ingoing EF coordinates. So, we
Taylor expand the modes near the horizon as,

$\delta g_{\mu\nu}(r)=\delta g_{\mu\nu}^{(0)}+(r-r_{h})\delta
g_{\mu\nu}^{(1)}+...\,,$ (3.45)

Perturbing the sound modes of the background and imposing the regularity
condition (3.45) near the black hole horizon, we can achieve the expanded form
of Einstein’s equations near the horizon. Our goal is to check whether we get
the results of chaos obtained in the previous subsections. For field $\delta
g_{vv}$, the associated equation of motion is $\delta E_{vv}=\delta T_{vv}$.
But, the perturbation to stress tensor $\delta T_{vv}$ does not vanish
necessarily, while $\delta T^{r}_{v}$ vanishes at the horizon Blake:2018leo ;
Wang:2022mcq . Keeping this in mind, we will set $\delta T^{r}_{v}=0$ to
calculate the pole-skipping point connected to chaos. Expanding this equation
near the black hole horizon, we get,

$k\left[2\omega-ir_{h}^{\xi}(1+\xi)\right]\delta
g_{vx}^{(0)}+\left[k^{2}-ir_{h}^{2-\xi}\omega+(r_{h}^{2}-r_{h}^{1+\xi})(1+\xi)\right]\delta
g_{vv}^{(0)}=0\,.$ (3.46)

This equation is a constraint relation that gives a relation between the near-
horizon coefficients. The above equation is very special at some specific
value of frequency($\omega$) and momentum ($k$). We can verify that, at
$\omega=\omega_{*}=\frac{i}{2}r_{h}^{\xi}(1+\xi)$, coefficient of $\delta
g_{vx}^{(0)}$ is 0. So, at this specific value of $\omega$, we get a relation,

$\left[k^{2}+\frac{1}{2}\left(3r_{h}-2r_{h}^{\xi}\right)\left(1+\xi\right)\right]\delta
g_{vv}^{(0)}=0\,.$ (3.47)

For general $k$, this equation gives rise to $\delta g_{vv}^{(0)}=0$. However,
at $k=k_{*}$, this equation is automatically satisfied. That means at this
specific value $\omega$ and $k$, the $\delta T^{r}_{v}$ is automatically zero
and hence they will not impose any constraint on the near-horizon coefficients
$\delta g_{vv}^{(0)}\,$ and $\delta g_{vx}^{(0)}$. This specific value of
frequency ($\omega_{*}$) and momentum ($k_{*}$) is called the pole-skipping
point. The pole-skipping point we get,

$\omega_{*}=\frac{i}{2}r_{h}^{\xi}(1+\xi)\,,\hskip
56.9055ptk_{*}^{2}=-\frac{1}{2}(3r_{h}-2r_{h}^{\xi})(1+\xi)\,.$ (3.48)

Writing $r_{h}$ in terms of temperature as $r_{h}=\left(\frac{4\pi
T}{1+\xi}\right)^{1/\xi}$, we can express $\omega$ and $k$ in terms of the
temperature of the black hole. Now, we can calculate the Lyapunov exponent and
butterfly velocity from (1.5) as,

$\lambda_{L}=2\pi T\,,\hskip 56.9055ptv_{B}=\sqrt{\frac{1+\xi}{2}}\,.$ (3.49)

Note that, the butterfly velocity ($v_{B}$) is a positive real quantity. So,
we will take the absolute value of momentum while calculating $v_{B}$, i.e.
$v_{B}=\frac{|\omega_{*}|}{|k_{*}|}$. From the above expression, we can see
that both the Lyapunov exponent and butterfly velocity take the same form as
we achieved in our previous subsections. So, we can say that all three methods
to compute the chaotic features are reliable. These methods are giving
consistent results in higher orders too. One can check it explicitly.

##  4 Analysis of classical chaos

Here we elucidate a parallel survey of some of the chaos properties, to be
precise, eikonal phase shift and Lyapunov exponent, at the classical level.
The anatomy of the null geodesics for our specified background and the
corresponding turning points adequately serve our purpose. We present a notion
of emergent dependency of classical chaos in our framework on the values of
these turning points.

###  4.1 Null geodesic

To carry out the eikonal phase shift in the bulk following, here we will first
enunciate the details of the null geodesic for our chosen background. The bulk
phase shift for a particle described by a plane wave is given by manuela ;
Parnachev:2020zbr

$\delta\equiv-p.(\Delta x)\,,$ (4.1)

where, $p$ denotes the components of the conserved momentum vector and $\Delta
x$ gives the deflections in the corresponding bulk directions. For our metric
in equation (2.2) in 2+1 dimension, the conserved energy and momentum are
given by

$p_{t}=\frac{1}{z^{2\xi}}f(z)\frac{\partial t}{\partial
s}\,,\leavevmode\nobreak\ \leavevmode\nobreak\
p_{x}=\frac{1}{z^{2}}\frac{\partial x}{\partial s}\,,$ (4.2)

$s$ being the intrinsic affine parameter. These conserved quantities
correspond to the two killing vectors $\partial_{t}$ and $\partial_{x}$ of our
chosen metric during a heavy-heavy-light-light 2-particle scattering, it is
convenient to consider that the highly boosted particle consists of large
values of conserved momenta. We introduce a parameter
$\gamma=\frac{p_{x}}{p_{t}}$ which often serves as the impact parameter of the
corresponding scattering process. The null geodesic equation is given by

$g_{\mu\nu}\dot{x}^{\mu}\dot{x}^{\nu}=0\,,$ (4.3)

which yields

$\dot{z}^{2}=(p_{t})^{2}z^{2\xi+2}-(p_{x})^{2}z^{4}f(z)\,.$ (4.4)

At the turning point $z=z_{0}$, $z$ attains extreme value for which we can set
$\dot{z}|_{z=z_{0}}=0$. Hence, the values of turning points $z_{0}$ should
satisfy the relation

$z_{0}^{2\xi-2}+\frac{\gamma^{2}}{z_{h}^{\xi+1}}z_{0}^{\xi+1}-\gamma^{2}=0\,.$
(4.5)

The above equation is a polynomial of either integer order or fractional
order, depending on the values of the critical exponent $\xi$. Thus the roots
of this equation which are valid for our approach are all real and lie within
the singularity and the boundary of the bulk geometry. However, these values
of the possible turning points depend on the ratio $\gamma$. At this stage,
some further crucial remarks on the turning point follow immediately:

  * •

Let us take the assumption that the energy $p_{t}$ is much greater than the
momentum $p_{x}$, i.e., $\gamma\gg 1$. With this assumption, taking $z_{0}=1$
equation (4.5) gives,

$\gamma=\sqrt{\frac{z_{h}^{\xi+1}}{z_{h}^{\xi+1}-1}}\,.$ (4.6)

It is obvious from the above expression that for any arbitrary value of
$z_{h}$ in between $z_{h}\rightarrow 0$ and $z_{h}\rightarrow\infty$, the
condition $\gamma\gg 1$ is not satisfied.

  * •

If $\gamma=1$, i.e., $p_{x}=p_{t}$, then the equation (4.5) is not satisfied
with $z_{0}=1$.

  * •

Again, for $z_{h}=1$, $z_{0}=1$ is not a solution of the null geodesic
equation for any arbitrary value of $\gamma$.

  * •

The condition (4.6) clearly tells us that, for $z_{0}=1$ as a solution,
$p_{x}$ cannot be less than $p_{t}$. Hence, with $p_{x}<p_{t}$, we must have
$z_{0}\neq 1$.

  * •

For any arbitrary value of $p_{x}$ and $p_{t}$; $z_{0}=0$ does not satisfy the
geodesic equation.

From the above conditions, it is evident that $z_{0}=1$, i.e., the horizon for
our case and $z_{0}=0$, i.e., the boundary are not good choices of $z_{0}$ to
carry out our desired study.

###  4.2 Eikonal phase shift

With the above analysis, here we will exhibit an explicit computation of
eikonal phase shift at the turning points $z_{0}$ as a nontrivial function of
$z_{0}$ and $\xi$. With $p_{t}$ and $p_{x}$ as the conserved energy and
momentum for our system, the expression of the bulk phase shift takes the form
as

$\delta\equiv p_{t}(\Delta t)-p_{x}(\Delta x)\,.$ (4.7)

The deflections along the time and space directions can be found respectively
from the expressions of $p_{t}$ and $p_{x}$ in (4.2) and the null geodesic
equation (4.4) as

$\Delta t\,=\,2\int_{0}^{z_{0}}\frac{z^{2\xi-1}\differential
z}{f(z)\sqrt{z^{2\xi}-\gamma^{2}z^{2}f(z)}}\,=\,2a^{\frac{3}{2}}\int_{0}^{z_{0}}\frac{z^{2\xi-1}\differential
z}{(a-z^{\xi+1})\sqrt{(az^{2\xi}-a\gamma^{2}z^{2}+\gamma^{2}z^{\xi+3})}}\,,$
(4.8)

and

$\Delta x\,=\,2\int_{0}^{z_{0}}\frac{\gamma z\,\differential
z}{\sqrt{z^{2\xi}-\gamma^{2}z^{2}f(z)}}\,=\,2a^{\frac{1}{2}}\gamma\int_{0}^{z_{0}}\frac{z\,\differential
z}{\sqrt{(az^{2\xi}-a\gamma^{2}z^{2}+\gamma^{2}z^{\xi+3})}}\,,$ (4.9)

where, $a=z_{h}^{\xi+1}$. Because of the intricate structure of the above
integrands, we assume two limiting cases to derive it. Firstly we take the
near-boundary zone where $z\leq z_{h}$. Substituting (4.8) and (4.9) in the
expression of the eikonal shift and simplifying, we get for the near boundary
zone,

$\delta\sim\int_{0}^{z_{0}}\sqrt{z^{2\xi-2}-\gamma^{2}(1-\epsilon)}\differential
z\approx p_{t}\int_{0}^{z_{0}}\sqrt{z^{2\xi-2}-\gamma^{2}}\differential z\,,$
(4.10)

where, $\left(\frac{z}{z_{h}}\right)^{\xi+1}=\epsilon\ll 1$ for $\xi\geq 1$.
Now defining a new integration variable $y=\frac{z}{z_{0}}$ the bulk phase
shift is given by

$\delta\sim\int_{0}^{1}z_{0}\sqrt{z_{0}^{2\xi-2}y^{2\xi-2}-\gamma^{2}}\differential
y\approx i\gamma
z_{0}\,{}_{2}F_{1}\left(-\frac{1}{2},\frac{1}{2\xi-2};\frac{2\xi-1}{2\xi-2};\frac{z_{0}^{2\xi-2}}{\gamma^{2}}\right)\,.$
(4.11)

On the other hand, we take the near horizon regime where $z\approx z_{h}$. In
this regime, we can take two liming cases for the impact parameter $\gamma$.
When the energy $p_{t}$ is much greater than the momentum $p_{x}$, i.e.,
$\gamma<1$, the eikonal phase shift in near horizon case gives

$\delta=\frac{z_{0}^{\xi}}{\xi(\xi+1)}\,_{2}F_{1}\left(1,\xi;\xi+1;\frac{z_{0}}{z_{h}}\right)=\frac{z_{0}^{\xi}}{\xi(\xi+1)}\,_{2}F_{1}\left(1,\xi;\xi+1;z_{0}\right)\,.$
(4.12)

Again for the small energy limit, i.e., $p_{x}\geq p_{t}$, the impact
parameter $\gamma\geq 1$. This, in the near horizon zone, produces,

$\delta=-\frac{4\gamma\left(\sqrt{z_{0}-z_{h}}-\sqrt{-z_{h}}\right)\sqrt{z_{h}}p_{t}}{\sqrt{\xi+1}}=-\frac{4\gamma\left(\sqrt{z_{0}-1}-i\right)p_{t}}{\sqrt{\xi+1}}\,.$
(4.13)

We have assumed $z_{h}=1$ in the above expression. It is well-known in the
backdrop of gravitational scattering that the real and imaginary nature of the
eikonal phase shift respectively represent elastic and inelastic scattering.
Thus our study of the classical eikonal phase clearly reveals that the nature
of the gravitational scattering between the heavy-heavy-light-light particle
collision in the Lifshitz black hole depends on the different regimes taken
for the radial coordinate $z$ along with the choices of the impact parameter.
For the near boundary case and the large impact parameter limit in the near
horizon zone, the imaginary eikonal phase shift implies completely inelastic
collisions whereas the small impact parameter limit in the near horizon zone
gives real values of $\delta$ that stand for elastic scattering. Nevertheless,
all the $\delta$’s derived in these limits inevitably happen to be nontrivial
functions of the anisotropy index $\xi$. This allows one to speculate about
the dependencies of absorption cross-sections of the above scatterings on the
finite anisotropy present in the black hole background.

Again, we use a WKB approximation to verify the consistency of the eikonal
phase. Inserting the expressions of time and space deflections in (4.7), we
can write,

$\delta=2p_{t}\int_{0}^{z_{0}}\frac{\sqrt{z^{2\xi}-\gamma^{2}z^{2}f(z)}}{zf(z)}\,\differential
z\,.$ (4.14)

One may check the emergence of the above eikonal phase via WKB approximation
of the solution of the equation of motion of a scalar field. To do this,
consider a heavy scalar field $\phi$ propagating in the chosen Lifshitz
background. In the background (2.2), the Klein-Gordon equation of the scalar
can be written as

$z^{2}f(z)\partial_{z}^{2}\phi+z^{2}\partial_{x}^{2}\phi+z^{2}f^{\prime}(z)\partial_{x}^{2}\phi-\frac{z^{2\xi}}{f(z)}\partial_{t}^{2}\phi-z^{\xi}f(z)\partial_{z}\phi+m^{2}z^{2+\xi}\phi=0\,.$
(4.15)

Performing a Fourier transformation as

$\phi=\tilde{\phi}(z)\,e^{i(-p_{t}t+p_{x}x)}\,,$

and taking an ansatz,

$\tilde{\phi}(z)=e^{ip_{t}\chi(z)}\,,$ (4.16)

we can solve for $\chi(z)$. With the large energy ($p_{t}>>1$) limit, the EOM
boils down to,

$\left(\partial_{z}\chi(z)\right)^{2}=\frac{z^{2\xi}-\gamma^{2}z^{2}f(z)}{z^{2}f^{2}(z)}\,.$
(4.17)

Now, taking the positive term of $\partial_{z}\chi(z)$, we perform an
integration and get a solution of $\chi(z)$ as

$\chi(z)=\int_{0}^{z_{0}}\frac{\sqrt{z^{2\xi}-\gamma^{2}z^{2}f(z)}}{zf(z)}\,\differential
z\,.$ (4.18)

where the terms of order $\mathcal{O}(p_{t}^{-1})$ can be ignored at large
energy limit. This is indeed the leading term of the solution. We can see that
the WKB approximation of the solution gives us the exact eikonal phase we
recovered in (4.14) at leading order.

###  4.3 Lyapunov exponent

The Lyapunov exponent ($\lambda_{L}$) is one of the key measures in
understanding chaotic systems in classical phase space. It quantifies the rate
at which the nearby trajectories diverge and provides information on the
system’s sensitivity to initial conditions. A positive Lyapunov exponent
signifies a chaotic system. Mathematically,

$\delta X(t)=\delta X(0)e^{\lambda_{L}t}\,,$

where $\delta X$ represents the separation between two nearby trajectories. In
their study Cardoso:2008bp , the authors established a connection between the
Lyapunov exponent and the effective potential in the radial motion for both
massive and massless particles. Building on this framework, we aim to
calculate the Lyapunov exponent for a 2+1-dimensional planar Lifshitz black
hole in this subsection. Recently, Kumara:2024obd  conducted a similar
analysis in the context of equatorial hyperscaling-violating black holes,
treating Lifshitz black holes as a specific case within their study. For a
circular geodesic, the Lyapunov exponent in terms of the second derivative of
the effective potential $V_{eff}$ for radial motion is given by Cardoso:2008bp
,Kumara:2024obd

$\lambda_{L}=\sqrt{-\frac{V_{eff}^{{}^{\prime\prime}}}{2\dot{t}^{2}}}\,.$
(4.19)

Now, for a massless particle, the geodesic equation is given in equation (4.4)
from which we can calculate the effective potential

$V_{eff}(z)=\dot{z}^{2}=(p_{t})^{2}z^{2\xi+2}-(p_{x})^{2}z^{4}f(z)\,.$ (4.20)

Where $t$ is the coordinate time. At the turning point $z=z_{0}$, the
condition on circular orbit, i.e., $V_{eff}(z_{0})=0$ and
$V^{\prime}_{eff}(z_{0})=0$ yield

$\frac{p_{t}^{2}}{p_{x}^{2}}=z_{0}^{2-2\xi}f_{0}\,,$ (4.21)

and

$z_{0}f_{0}^{\prime}=2(\xi-1)f_{0}\,,$ (4.22)

where $f_{0}=f(z_{0})$. Upon substituting equations (4.21) and (4.22) in
equation (4.20) we get

$V_{eff}^{\prime\prime}(z_{0})=-z_{0}^{2}p_{x}^{2}\left(z_{0}^{2}f_{0}^{\prime\prime}-2f_{0}\left(2\xi^{2}-5\xi+3\right)\right)\,.$
(4.23)

Now, utilizing equations (4.2) and (4.23) within equation (4.19), we derive
the Lyapunov exponent for the null circular geodesic

$\lambda_{L}=\sqrt{\frac{f_{0}\left(z_{0}^{2}f_{0}^{\prime\prime}-2f_{0}\left(2\xi^{2}-5\xi+3\right)\right)}{2z_{0}^{2\xi}}}\,.$
(4.24)

Since an unstable circular geodesic meets the condition
$V^{\prime\prime}_{eff}<0$, it follows from equations (4.23) and (4.24) that
$\lambda_{L}$ is real whenever the orbit is unstable. It is again
straightforward to note that the Lyapunov exponent goes to zero as the turning
point $z_{0}$ reaches the horizon and becomes infinity when it reaches at
boundary $z_{0}=0$, both of which are inconsistent with the expected scenario.
Therefore, the above approach of studying the Lyapunov exponent of the theory
under consideration does not remain acceptable when the extrema of the null
geodesics are chosen in the vicinity of both the boundary and the horizon of
the black hole. This also supports our comments in section 4.1 on the possible
restrictions on the values $z_{0}=0$ and $z_{0}=1$. At this stage, it is
faithful to say that one may carry out a polynomial distribution by using a
suitable numerical approach to see what values of $z_{0}$ the equation (4.5)
can accommodate as its valid solution. We keep this computation for future
extensions of our work.

##  5 Conclusion

We undergo a comparative analysis of three different holography-based
approaches to quantum chaos to understand any probable equivalence between
them. The asymptotically Lifshitz black hole with arbitrary anisotropy index
is chosen to be a fertile gravity background for our work due to its chaotic
features which are already explored in various theoretical contexts. The
methods involved in this work for studying quantum chaos are entanglement
wedge reconstruction, out-of-time-ordered correlators and pole-skipping. All
of these approaches are eventually based on the well-known gauge/gravity
duality. In the method of entanglement wedge reconstruction, we construct an
extremal hypersurface as a constant time slice, famously called Ryu-Takayanagi
surface, that extracts out the butterfly velocity at late time. The butterfly
velocity occurs to have nontrivial dependence on the anisotropy index of the
chosen background. For the unit anisotropy index, it expectantly reduces to
that for the planar AdS black hole. For the derivation of OTOCs, we explicitly
construct the double-sided Kruskal geometry corresponding to the
asymptotically Lifshitz black hole. We introduce a gravitational shockwave
backreacting in the Kruskal version of our background. Subsequently, the
functional form of the shockwave assumes an exponential form which shows that
the system is chaotic and this shockwave profile gives rise to the four-point
OTOCs. Such a scenario is somewhat consistent with the previously obtained
scenario for a perturbed two-sided AdS black hole. The butterfly velocity and
Lyapunov exponent obtained from this shockwave profile are exactly similar to
those obtained from the entanglement wedge method. Finally, in the method of
pole-skipping, we again observe the exact same dependency of butterfly
velocity as well as the Lyapunov exponent on the anisotropy. As a consequence
of the exact matching of the salient chaotic features, we claim that all our
chosen methods are equivalent and definitive for the understanding of chaos
for the class of asymptotically Lifshitz black holes at the quantum level.

We demonstrate a further investigation of the eikonal phase shift due to
heavy-heavy-light-light gravitational scattering and computed the Lyapunov
exponent, at the classical level by using the bulk metric both in near-
boundary and near-horizon regimes and a detailed study of the associated null
geodesics. This reveals that both of these quantities appear to be nontrivial
functions of the anisotropy index $\xi$. We further achieve some obvious
restrictions on the possible extrema of the null geodesics where we obtain
inconsistent natures of classical eikonal phase shift and Lyapunov exponent.
An interesting outcome of the presence of arbitrary anisotropy in our case is
that the eikonal phase becomes real or imaginary depending on the near-horizon
and near boundary limits of the black hole and the impact parameter. This may
allow one to explore the absorption cross-section for both elastic and
inelastic scattering in terms of the arbitrary anisotropy index of the
asymptotically Lifshitz black hole. We have further verified the eikonal phase
by using a WKB approximation of the solution of the Klein-Gordon equation. We
understood that in the leading order WKB approximation of the scalar solution,
the phase we recovered is the same as that derived from the analysis of the
null geodesic equation. It is significant to note that, the eikonal phase
approximation from the derivation of the out-of-time-order correlators in
(3.38) for near-horizon limit assumes an exponential form similar to $e^{(at-
bx)}$ with a $\xi$-dependent prefactor, $a$ and $b$ being some $\xi$ and
$\beta$ dependent functions. Whereas, for the classical case, it becomes a
hypergeometric function of $\xi$ for the near-boundary regime. It would be
further interesting to explore the possible connection between these two
results which may eventually lead to classical/quantum correspondence in the
context of chaos in Lifshitz black hole. Moreover, it would be rather
fascinating to extend such a relation between the Lyapunov exponents obtained
from quantum and classical analyses. Another intriguing opening is to
explicitly evaluate the OTOCs from the systematic construction of the dual
thermofield double states of the finite temperature Lifshitz field theory and
subsequently study different chaotic features on the field theory side. We
hope to come back with some of these ideas in the near future.

###### Acknowledgements.

We would like to thank Debaprasad Maity for his valuable suggestions regarding
this work. A special thanks to Diandian Wang and Joydeep Chakravarty for some
helpful suggestions in our latest draft. BB is grateful to Arnab Kundu for
hosting a visit to SINP, Kolkata, and to Pankaj Chaturvedi for facilitating a
visit to NIT, Silchar, where part of this work was conducted. Additionally, BB
would like to thank Akhil Sivakumar for various insightful suggestions on this
work. AC appreciates Bibhas Ranjan Majhi for hosting a visit to IIT Guwahati,
which marked the initiation of this work. AC also acknowledges the financial
support from the NSTC, Taiwan (R.O.C.), under grant number
110-2112-M007-015-MY3, which funded a significant portion of this research.
Finally, AC extends sincere thanks to Chong-Sun Chu, Dimitrios Giataganas,
Himansu Parihar, and Jaydeep Kumar Basak for their insightful discussions on
various aspects of Lifshitz theories, which were extremely beneficial to this
work.

## References

  * (1) J.M. Maldacena, _The Large N limit of superconformal field theories and supergravity_ , _Adv. Theor. Math. Phys._ 2 (1998) 231 [hep-th/9711200].
  * (2) O. Aharony, S.S. Gubser, J.M. Maldacena, H. Ooguri and Y. Oz, _Large N field theories, string theory and gravity_ , _Phys. Rept._ 323 (2000) 183 [hep-th/9905111].
  * (3) K. Balasubramanian and J. McGreevy, _Gravity duals for nonrelativistic conformal field theories_ , _Phys. Rev. Lett._ 101 (2008) 061601.
  * (4) W.D. Goldberger, _Ads/cft duality for non-relativistic field theory_ , _Journal of High Energy Physics_ 2009 (2009) 069.
  * (5) J. Maldacena, D. Martelli and Y. Tachikawa, _Comments on string theory backgrounds with non-relativistic conformal symmetry_ , _Journal of High Energy Physics_ 2008 (2008) 072.
  * (6) A. Adams, K. Balasubramanian and J. McGreevy, _Hot spacetimes for cold atoms_ , _Journal of High Energy Physics_ 2008 (2008) 059.
  * (7) S. Kachru, X. Liu and M. Mulligan, _Gravity duals of lifshitz-like fixed points_ , _Phys. Rev. D_ 78 (2008) 106005.
  * (8) V. Keränen, W. Sybesma, P. Szepietowski and L. Thorlacius, _Correlation functions in theories with lifshitz scaling_ , _Journal of High Energy Physics_ 2017 (2017) .
  * (9) M. Kord Zangeneh, B. Wang, A. Sheykhi and Z. Tang, _Charged scalar quasi-normal modes for linearly charged dilaton-lifshitz solutions_ , _Physics Letters B_ 771 (2017) 257.
  * (10) C. Park, _Massive quasinormal mode in the holographic lifshitz theory_ , _Phys. Rev. D_ 89 (2014) 066003.
  * (11) J. Cheyne and D. Mattingly, _Constructing entanglement wedges for Lifshitz spacetimes with Lifshitz gravity_ , _Phys. Rev. D_ 97 (2018) 066024 [1707.05913].
  * (12) J.K. Basak, A. Chakraborty, C.-S. Chu, D. Giataganas and H. Parihar, _Massless Lifshitz Field Theory for Arbitrary $z$_, 2312.16284.
  * (13) X. Bai, B.-H. Lee, T. Moon and J. Chen, _Chaos in Lifshitz Spacetimes_ , _J. Korean Phys. Soc._ 68 (2016) 639 [1406.5816].
  * (14) D. Giataganas and K. Sfetsos, _Non-integrability in non-relativistic theories_ , _JHEP_ 06 (2014) 018 [1403.2703].
  * (15) E. Plamadeala and E. Fradkin, _Scrambling in the quantum lifshitz model_ , _Journal of Statistical Mechanics: Theory and Experiment_ 2018 (2018) 063102.
  * (16) N. Sircar, J. Sonnenschein and W. Tangarife, _Extending the scope of holographic mutual information and chaotic behavior_ , _JHEP_ 05 (2016) 091 [1602.07307].
  * (17) J.S. Cotler, G. Gur-Ari, M. Hanada, J. Polchinski, P. Saad, S.H. Shenker et al., _Black Holes and Random Matrices_ , _JHEP_ 05 (2017) 118 [1611.04650].
  * (18) L. Susskind, _Computational Complexity and Black Hole Horizons_ , _Fortsch. Phys._ 64 (2016) 24 [1403.5695].
  * (19) 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].
  * (20) J.M. Magán, _Black holes, complexity and quantum chaos_ , _JHEP_ 09 (2018) 043 [1805.05839].
  * (21) L. Pausch, E.G. Carnio, A. Rodríguez and A. Buchleitner, _Chaos and ergodicity across the energy spectrum of interacting bosons_ , _Physical Review Letters_ 126 (2021) .
  * (22) P. Kos, B. Bertini and T.c.v. Prosen, _Chaos and ergodicity in extended quantum systems with noisy driving_ , _Phys. Rev. Lett._ 126 (2021) 190601.
  * (23) X. Dong, D. Wang, W.W. Weng and C.-H. Wu, _A tale of two butterflies: an exact equivalence in higher-derivative gravity_ , _JHEP_ 10 (2022) 009 [2203.06189].
  * (24) M. Mezei and J. Virrueta, _Exploring the Membrane Theory of Entanglement Dynamics_ , _JHEP_ 02 (2020) 013 [1912.11024].
  * (25) W. Fischler, V. Jahnke and J.F. Pedraza, _Chaos and entanglement spreading in a non-commutative gauge theory_ , _JHEP_ 11 (2018) 072 [1808.10050].
  * (26) S.H. Shenker and D. Stanford, _Black holes and the butterfly effect_ , _JHEP_ 03 (2014) 067 [1306.0622].
  * (27) J. Maldacena, S.H. Shenker and D. Stanford, _A bound on chaos_ , _JHEP_ 08 (2016) 106 [1503.01409].
  * (28) S.H. Shenker and D. Stanford, _Multiple Shocks_ , _JHEP_ 12 (2014) 046 [1312.3296].
  * (29) D.A. Roberts, D. Stanford and L. Susskind, _Localized shocks_ , _JHEP_ 03 (2015) 051 [1409.8180].
  * (30) S.H. Shenker and D. Stanford, _Stringy effects in scrambling_ , _JHEP_ 05 (2015) 132 [1412.6087].
  * (31) D.A. Roberts and D. Stanford, _Two-dimensional conformal field theory and the butterfly effect_ , _Phys. Rev. Lett._ 115 (2015) 131603 [1412.5123].
  * (32) S. Khetrapal, _Chaos and operator growth in 2d CFT_ , _JHEP_ 03 (2023) 176 [2210.15860].
  * (33) C.-J. Lin and O.I. Motrunich, _Out-of-time-ordered correlators in a quantum ising chain_ , _Physical Review B_ 97 (2018) .
  * (34) Y. Gu and X.-L. Qi, _Fractional Statistics and the Butterfly Effect_ , _JHEP_ 08 (2016) 129 [1602.06543].
  * (35) S. Das, B. Ezhuthachan, A. Kundu, S. Porey, B. Roy and K. Sengupta, _Out-of-Time-Order correlators in driven conformal field theories_ , _JHEP_ 08 (2022) 221 [2202.12815].
  * (36) M. Blake, H. Lee and H. Liu, _A quantum hydrodynamical description for scrambling and many-body chaos_ , _JHEP_ 10 (2018) 127 [1801.00010].
  * (37) M. Blake, R.A. Davison, S. Grozdanov and H. Liu, _Many-body chaos and energy dynamics in holography_ , _JHEP_ 10 (2018) 035 [1809.01169].
  * (38) S. Grozdanov, _On the connection between hydrodynamics and quantum chaos in holographic theories with stringy corrections_ , _JHEP_ 01 (2019) 048 [1811.09641].
  * (39) Y. Ahn, V. Jahnke, H.-S. Jeong and K.-Y. Kim, _Scrambling in Hyperbolic Black Holes: shock waves and pole-skipping_ , _JHEP_ 10 (2019) 257 [1907.08030].
  * (40) M. Blake and R.A. Davison, _Chaos and pole-skipping in rotating black holes_ , _JHEP_ 01 (2022) 013 [2111.11093].
  * (41) K. Sil, _Pole skipping and chaos in anisotropic plasma: a holographic study_ , _JHEP_ 03 (2021) 232 [2012.07710].
  * (42) W. Li, S. Lin and J. Mei, _Thermal diffusion and quantum chaos in neutral magnetized plasma_ , _Phys. Rev. D_ 100 (2019) 046012 [1905.07684].
  * (43) M.A.G. Amano, M. Blake, C. Cartwright, M. Kaminski and A.P. Thompson, _Chaos and pole-skipping in a simply spinning plasma_ , _JHEP_ 02 (2023) 253 [2211.00016].
  * (44) N. Abbasi and J. Tabatabaei, _Quantum chaos, pole-skipping and hydrodynamics in a holographic system with chiral anomaly_ , _JHEP_ 03 (2020) 050 [1910.13696].
  * (45) B. Baishya and K. Nayek, _Probing pole-skipping through scalar Gauss-Bonnet coupling_ , _Nucl. Phys. B_ 1001 (2024) 116521 [2301.03984].
  * (46) B. Baishya, S. Chakrabarti, D. Maity and K. Nayek, _Pole-skipping and chaos in D3-D7 brane_ , 2312.01829.
  * (47) H. Yuan and X.-H. Ge, _Pole-skipping and hydrodynamic analysis in Lifshitz, AdS 2 and Rindler geometries_, _JHEP_ 06 (2021) 165 [2012.15396].
  * (48) X. Dong, D. Harlow and A.C. Wall, _Reconstruction of bulk operators within the entanglement wedge in gauge-gravity duality_ , _Phys. Rev. Lett._ 117 (2016) 021601.
  * (49) A.C. Wall, _Maximin Surfaces, and the Strong Subadditivity of the Covariant Holographic Entanglement Entropy_ , _Class. Quant. Grav._ 31 (2014) 225007 [1211.3494].
  * (50) B. Czech, J.L. Karczmarek, F. Nogueira and M. Van Raamsdonk, _The Gravity Dual of a Density Matrix_ , _Class. Quant. Grav._ 29 (2012) 155009 [1204.1330].
  * (51) S. Ryu and T. Takayanagi, _Holographic derivation of entanglement entropy from the anti–de sitter space/conformal field theory correspondence_ , _Phys. Rev. Lett._ 96 (2006) 181602.
  * (52) M. Mezei and D. Stanford, _On entanglement spreading in chaotic systems_ , _JHEP_ 05 (2017) 065 [1608.05101].
  * (53) M. Taylor, _Lifshitz holography_ , _Class. Quant. Grav._ 33 (2016) 033001 [1512.03554].
  * (54) Y. Sekino and L. Susskind, _Fast Scramblers_ , _JHEP_ 10 (2008) 065 [0808.2096].
  * (55) P. Hayden and J. Preskill, _Black holes as mirrors: Quantum information in random subsystems_ , _JHEP_ 09 (2007) 120 [0708.4025].
  * (56) M. Blake, R.A. Davison and D. Vegh, _Horizon constraints on holographic Green’s functions_ , _JHEP_ 01 (2020) 077 [1904.12883].
  * (57) M. Natsuume and T. Okamura, _Nonuniqueness of Green’s functions at special points_ , _JHEP_ 12 (2019) 139 [1905.12015].
  * (58) S. Grozdanov, K. Schalm and V. Scopelliti, _Black hole scrambling from hydrodynamics_ , _Phys. Rev. Lett._ 120 (2018) 231601 [1710.00921].
  * (59) X. Wu, _Higher curvature corrections to pole-skipping_ , _JHEP_ 12 (2019) 140 [1909.10223].
  * (60) D.M. Ramirez, _Chaos and pole skipping in CFT 2_, _JHEP_ 12 (2021) 006 [2009.00500].
  * (61) F.M. Haehl, W. Reeves and M. Rozali, _Reparametrization modes, shadow operators, and quantum chaos in higher-dimensional CFTs_ , _JHEP_ 11 (2019) 102 [1909.05847].
  * (62) S. Das, B. Ezhuthachan and A. Kundu, _Real time dynamics from low point correlators in 2d BCFT_ , _JHEP_ 12 (2019) 141 [1907.08763].
  * (63) N. Abbasi and K. Landsteiner, _Pole-skipping as order parameter to probe a quantum critical point_ , 2307.16716.
  * (64) M. Natsuume and T. Okamura, _Pole-skipping with finite-coupling corrections_ , _Phys. Rev. D_ 100 (2019) 126012 [1909.09168].
  * (65) N. Ceplak and D. Vegh, _Pole-skipping and Rarita-Schwinger fields_ , _Phys. Rev. D_ 103 (2021) 106009 [2101.01490].
  * (66) K.-Y. Kim, K.-S. Lee and M. Nishida, _Holographic scalar and vector exchange in OTOCs and pole-skipping phenomena_ , _JHEP_ 04 (2021) 092 [2011.13716].
  * (67) S. Grozdanov, _Bounds on transport from univalence and pole-skipping_ , _Phys. Rev. Lett._ 126 (2021) 051601 [2008.00888].
  * (68) H.-S. Jeong, K.-Y. Kim and Y.-W. Sun, _Bound of diffusion constants from pole-skipping points: spontaneous symmetry breaking and magnetic field_ , _JHEP_ 07 (2021) 105 [2104.13084].
  * (69) D. Wang and Z.-Y. Wang, _Pole Skipping in Holographic Theories with Bosonic Fields_ , _Phys. Rev. Lett._ 129 (2022) 231603 [2208.01047].
  * (70) Y.-T. Wang and W.-B. Pan, _Pole-skipping of holographic correlators: aspects of gauge symmetry and generalizations_ , _JHEP_ 01 (2023) 174 [2209.04296].
  * (71) H. Yuan, X.-H. Ge, K.-Y. Kim, C.-W. Ji and Y. Ahn, _Pole-skipping points in 2D gravity and SYK model_ , _JHEP_ 08 (2023) 157 [2303.04801].
  * (72) M. Natsuume and T. Okamura, _Pole skipping in a non-black-hole geometry_ , _Phys. Rev. D_ 108 (2023) 046012 [2306.03930].
  * (73) B. Baishya, S. Chakrabarti and D. Maity, _Effect of scalar condensation on fermionic Pole-Skipping_ , 2311.05314.
  * (74) D.N. Kabat and M. Ortiz, _Eikonal quantum gravity and Planckian scattering_ , _Nucl. Phys. B_ 388 (1992) 570 [hep-th/9203082].
  * (75) G. D’Appollonio, P. Di Vecchia, R. Russo and G. Veneziano, _High-energy string-brane scattering: Leading eikonal and beyond_ , _JHEP_ 11 (2010) 100 [1008.4773].
  * (76) N.E.J. Bjerrum-Bohr, J.F. Donoghue, B.R. Holstein, L. Planté and P. Vanhove, _Bending of Light in Quantum Gravity_ , _Phys. Rev. Lett._ 114 (2015) 061301 [1410.7590].
  * (77) L. Cornalba, M.S. Costa, J. Penedones and R. Schiappa, _Eikonal Approximation in AdS/CFT: From Shock Waves to Four-Point Functions_ , _JHEP_ 08 (2007) 019 [hep-th/0611122].
  * (78) L. Cornalba, M.S. Costa, J. Penedones and R. Schiappa, _Eikonal Approximation in AdS/CFT: Conformal Partial Waves and Finite N Four-Point Functions_ , _Nucl. Phys. B_ 767 (2007) 327 [hep-th/0611123].
  * (79) R.C. Brower, M.J. Strassler and C.-I. Tan, _On the eikonal approximation in AdS space_ , _JHEP_ 03 (2009) 050 [0707.2408].
  * (80) K. Balasubramanian and J. McGreevy, _An Analytic Lifshitz black hole_ , _Phys. Rev. D_ 80 (2009) 104039 [0909.0263].
  * (81) M. Taylor, _Lifshitz holography_ , _Classical and Quantum Gravity_ 33 (2016) 033001.
  * (82) U.H. Danielsson and L. Thorlacius, _Black holes in asymptotically Lifshitz spacetime_ , _JHEP_ 03 (2009) 070 [0812.5088].
  * (83) V. Jahnke, K.-Y. Kim and J. Yoon, _On the Chaos Bound in Rotating Black Holes_ , _JHEP_ 05 (2019) 037 [1903.09086].
  * (84) M. Mezei and G. Sárosi, _Chaos in the butterfly cone_ , _JHEP_ 01 (2020) 186 [1908.03574].
  * (85) A.P. Reynolds and S.F. Ross, _Butterflies with rotation and charge_ , _Class. Quant. Grav._ 33 (2016) 215008 [1604.04099].
  * (86) R.R. Poojary, _BTZ dynamics and chaos_ , _JHEP_ 03 (2020) 048 [1812.10073].
  * (87) A. Banerjee, A. Kundu and R.R. Poojary, _Rotating black holes in AdS spacetime, extremality, and chaos_ , _Phys. Rev. D_ 102 (2020) 106013 [1912.12996].
  * (88) A. Saha and S. Gangopadhyay, _Quantum chaos in the presence of non-conformality_ , 2401.05814.
  * (89) A. Alvarez, E. Ayón-Beato, H.A. González and M. Hassaïne, _Nonlinearly charged Lifshitz black holes for any exponent $z>1$_, _JHEP_ 06 (2014) 041 [1403.5985].
  * (90) K.-X. Zhu, F.-W. Shu and D.-H. Du, _Holographic complexity for nonlinearly charged lifshitz black holes_ , _Classical and Quantum Gravity_ 37 (2020) 195023.
  * (91) M. Kulaxizi, G.S. Ng and A. Parnachev, _Black holes, heavy states, phase shift and anomalous dimensions_ , _SciPost Phys._ 6 (2019) 065.
  * (92) A. Parnachev and K. Sen, _Notes on AdS-Schwarzschild eikonal phase_ , _JHEP_ 03 (2021) 289 [2011.06920].
  * (93) V. Cardoso, A.S. Miranda, E. Berti, H. Witek and V.T. Zanchin, _Geodesic stability, Lyapunov exponents and quasinormal modes_ , _Phys. Rev. D_ 79 (2009) 064016 [0812.1806].
  * (94) A.N. Kumara, S. Punacha and M.S. Ali, _Lyapunov Exponents and Phase Structure of Lifshitz and Hyperscaling Violating Black Holes_ , 2401.05181.


Methods

We developed our analysis in different steps: first, by considering NLP tools
to transform textual objects into numerical ones, an operation called
embedding. Then by applying some well-known machine learning algorithms to
extract useful information from the dataset.

Natural Language Processing Natural Language Processing, or NLP, is a
computational technique that allows language analysis. During this work, NLP
has been used to analyze Latin texts. In particular, we handled some
examinations that are usually carried out manually in ancient languages. The
most common tools to carry out NLP analysis work on modern languages only. For
this reason, analyzing ancient texts with machine learning is usually a
demanding job. Additionally, another challenge encountered is the
unavailability of ancient works in text formats. Usually, OCR techniques are
needed, and it is not easy to find sources online. However, Latin Library
provides a database of text files. In this way, it has not been necessary to
focus on the translation step.

LatinBERT allowed us to improve our analysis since it was a pre-trained model
to work on ancient Latin. bollmann2018normalization In Figure img:graph, the
steps followed to extract information from the text are summarized.

Preprocessing

The preprocessing phase is conducted to prepare the text for further analysis.
The goal is to clean the raw text and make it consistent, to apply algorithms
for the embedding. Cleaning means, for instance, handling punctuation,
removing uninformative words or correcting spelling errors. In particular, the
preprocessing operation applied in this analysis are the following:

  * •

Removing numbers and uninformative words: numbers do not possess semantic
meaning, so they are removed. The same goes for misspelt words that the
algorithm is not able to correct referring to a standard Latin vocabulary.

  * •

Converting to lower cases to standardize the text.

  * •

Removing punctuation: for analyzing whole documents we also remove the
punctuation, while we retain it when we analyzed the texts at the sentence
level.

###  Word embedding

The computational analysis of written documents requires the transformation of
the text into mathematical objects on which we can perform operations. This is
done by assigning a unique numerical vector to each word or sentence, creating
a projection of the terms in a high dimensional vector space. We refer to this
process as word embedding [turian-etal-2010-word].
Numerous algorithms in the literature can be employed to achieve this goal
[w2v], [d2v], [bert], but little has been done about ancient languages. In the
present work, we make use of LatinBERT, a contextual language model for Latin,
trained on 642.7 million words from various sources from the classical era to
the 21st century (Corpus Thomisticum, Internet Archive, Latin Library,
Patrologia Latina, Perseus, Latin Wikipedia) [bamman2020latin].
LatinBERT, as all BERT models, produces contextual embeddings of the words.
Hence, each term is expressed as a unique vector, and its value depends on the
context in which the word appears. BERT models take all the words appearing in
a sentence and infer the relationship between them. This approach allows the
model to deal with polysemic expressions, associating them with diverse
embeddings when representing different concepts.
BERT models have proven to be more effective than other embedding techniques
because, when they compute the embedding for a word, they take into
consideration all previous and next words in the same sentence. The model was
made available by the authors for future research in the field.

##  Similarities measures

By applying Latin-Bert it is possible to encode documents into vectors and
then exploit their geometrical properties to define similarity between the
native texts. The most intuitive way to extract information about the
closeness of two vectors in space is through the Euclidean distance, which
naturally represents the shortest distance between two points.
An alternative way to quantify such closeness is by computing the cosine
similarity $\phi$ of pairs of vectors, a metric often used in NLP to measure
similarity between words or texts. It is measured by the cosine of the angle
$\theta$ between two vectors $V$ and $W$, and it determines whether the two
vectors point in the same direction [datam]:

$\phi=\cos(\theta)=\frac{\mathbf{v}\cdot\mathbf{w}}{\|\mathbf{v}\|\|\mathbf{w}\|}=\frac{\sum_{i=1}^{n}v_{i}w_{i}}{\sqrt{\sum_{i=1}^{n}v_{i}^{2}}\sqrt{\sum_{i=1}^{n}w_{i}^{2}}}.$
(1)

Here $v_{i}$ and $w_{i}$ are the components of the $n-$dimensional vectors v
and w, in our case standing for the word embedding vectors of the documents.
The value of $\phi$ ranges between -1 and 1. Values of $\phi$ closer to -1
indicate dissimilarity, while values of $\phi$ closer to 1 mean similarity
between the two objects. If $\phi$ is close to zero, then no significant
correlation between the two texts exists.

##  Clustering Analysis

Complementary with computing similarity measures between texts, we can group
them based on their closeness in the multi-dimensional space they lie with
clustering analysis. Clustering algorithms belong to the general class of
Unsupervised ML techniques and apply when the objects divide into natural
groups, referred to as clusters. There exist numerous algorithms that can be
used for this purpose. The one we selected is one of the most popular, the
K-means algorithm. It consists of grouping the vectors by evaluating the
euclidean distance between them, resulting in the assignment of each data
point to the nearest cluster [10.5555/1162264].

##  Visualization

A further step that must be taken into account when working with digital
humanities is making the results intuitive and intelligible. The best way to
deal with such interpretability is to create a graphical representation of the
results.
A way to visualize documents in such a high-dimensional numerical form is
employing dimensionality reduction. Since the visualization of points in
dimensions 4 and above is impossible, it’s necessary to reduce them to be
represented as points in dimensions 2 or 3. In our work, we applied the
Uniform Manifold Approximation and Projection (UMAP) embedding algorithm to
reduce to two the dimension of the embedding vectors. UMAP exploits
geometrical and topological features of the vectors to preserves their global
structure: the embedding is found by searching for a low dimensional
projection of the data that has the closest possible equivalent topological
structure [mcinnes2018umap]. Once we obtain such low-dimensional projection,
we can visualize the vectors in the plane and observe how they cluster and
their relationships.
The best way to represent relationships between objects is by employing
graphs. A graph $G$ is an ordered pair $(V,E)$, where $V$ is a set of nodes
and $E$ is a set of edges. Nodes represent items, while edges represent
relationships between them and can be weighted with geometric measures
capturing the degree of such relationships, such as the cosine similarity or
the euclidean distance previously introduced. The information about the
connectivity patterns and the weights are collected in the connectivity matrix
$W$ of the graph. In the present study, nodes represent authors, while the
links represent the connection between them and are weighted with their degree
of similarity.