††thanks<EMAIL_ADDRESS>
# Probing Causation Dynamics in Quantum Chains near Criticality

Roopayan Ghosh  Department of Physics and Astronomy, University College
London, Gower Street, WC1E 6BT London, United Kingdom.  Bin Yi  Department of
Physics and Astronomy, University College London, Gower Street, WC1E 6BT
London, United Kingdom.  Institute of Fundamental and Frontier Sciences,
University of Electronic Science and Technology of China, Chengdu 610051,
China
Sougato Bose  Department of Physics and Astronomy, University College London,
Gower Street, WC1E 6BT London, United Kingdom.

###### Abstract

Distinguishing causation from correlation is crucial and requires careful
consideration. In this study, we utilize a recent quantum extension of Liang
information to investigate causation in quantum chains across their phase
diagram. Our analysis encompasses two distinct scenarios: (i) the Aubry-André-
Harper model, characterized by a spectrum-wide phase transition, and (ii) the
Anisotropic Transverse Field Ising (ANNNI) model, which demonstrates a ground
state transition. We discern a notable shift in causation behavior across the
critical point with each case exhibiting distinct hallmarks, different from
correlation measures. Especially, in the latter case, we observe maximum
causation in the ordered phase just preceding the critical point.

#### Introduction:

The study of dynamics of quantum many-body systems typically involve the time
evolution of correlation functions and the spreading of entanglement [1].
Quantifying causation dynamics in the quantum realm, unlike classical
theories, has been challenging [2]. Descriptions of causation dynamics have
been mostly limited to quantifying scrambling and relating it to chaos by
modifying classical concepts such the Lyapunov exponents. One such attempt
yielded the Out-of-time-ordered-Correlator (OTOC), a reasonably good
quantifier of scrambling but having limited applications to identify quantum
chaos unlike originally anticipated [3, 4, 5, 6, 7, 8]. Hence very recently
other modes of quantification of quantum chaos has been proposed [9, 10].
Nevertheless, the adage "correlation does not always imply causation" holds
true, and in the context of quantum mechanics, a distinct quantifier for
quantum causation remains less established.

Model Hamiltonians of quantum chains provide the microcosm of quantum effects
prevalent in our universe. Hence, it is of paramount importance to find a
simple detector of causation in these systems, which is easily measurable and
can be appropriately connected to the physical intuition born from classical
systems. This led to efforts which extract quantum causation from quantum
correlations (which are different from the classical correlations), drawing
inspiration from the Liang-Kleeman analysis used in classical systems[11]. In
classical systems, this analysis was originally formulated with respect to
Shannon entropy, and then extended to relative entropy[12], which provides a
better measure of predictability than Shannon entropy in certain cases[13, 14,
15, 16]. For quantum mechanical systems, this formalism was adapted using von
Neumann entropy analogous to Shannon entropy of classical systems.[17] The two
big advantages of the quantum version of Liang information over specific
observable based correlation measures such as OTOC are (i) the usage of
information theoretic tools which makes it more universal, and (ii) it is
intuitively connected to the classical picture, and is practically easy to
implement in expererimental setups as its simplest version requires single
site measurements.

Figure 1: Schematic diagram of how Liang information flow is computed between
two sites by freezing one site (denoted by red cross) and the target site
denoted by the white colour. See text for details.

Furthermore, understanding the behaviour of the quantum causation near
critical points is completely unexplored. This is the void we wish to fill in
this work by using the quantum Liang information formulated with respect to
von Neumann entropy as a causation quantifier (note that the direct usage of
von Neumann entropy just measures quantum correlations, not causation). To
achieve this goal, we conduct simulations of non-trivial dynamics in a quantum
state using model Hamiltonians, typically chosen such that their eigenstates
exclude the initial state—–a scenario known as a quantum quench. While
previous research has extensively explored the time evolution of correlation
functions for quenches across critical points in quantum many-body systems[18,
19], leveraging non-analyticities in these functions to detect phase
transitions[20, 21], to the best of our knowledge, there has been no
significant investigation into causality within such a setup. Although some
studies involving Out-of-time-ordered-Correlator (OTOC) are available in the
literature[22, 23, 24], as mentioned before, their behaviour frequently
depends on the choice of operator and we have to be mindful of conservation
laws on the system. On the other hand, our study using Liang information
suffers from no such drawback and we shall use it to examine the behavior of
causation during a quench across critical points in a chain geometry .

#### Definitions:

For density operator $\rho$, the von Neumann entropy $S(\rho)$ is defined as:

$S(\rho)=-{\rm Tr}[\rho\log\rho]$ (1)

Consider a bipartite state $\rho$ evolve under unitary operator $U(t)$
generated from Hamiltonian function $H$, the reduced density state of system
AB and subsystems A, B are then denoted $\rho_{AB}$, $\rho_{A}$, $\rho_{B}$.
Following the methodology of its classical counterpart[12], quantum Liang
information flow is then given by

$T^{r}_{B\rightarrow
A}=\frac{dS(\rho_{A})}{dt}-\frac{dS(\rho_{A\not{B}})}{dt}$ (2)

the cumulative information flow is then

$\mathbb{T}_{B\rightarrow A}=\int T^{r}_{B\rightarrow A}dt=\Delta
S(\rho_{A})-\Delta S(\rho_{A\not{B}}).$ (3)

Here $\rho_{A\not{B}}$ refers to the density state of A evolving with B
frozen[17].

#### Setup:

The schematic setup for investigating Liang information flow is illustrated in
Fig. 1. We begin by considering a chain with only nearest-neighbor couplings.
While the target site can be either a single site or a group, for simplicity,
we focus on a single site denoted by the violet marker in the schematic figure
throughout this study. In the top three chains, we depict situations with only
nearest-neighbor couplings. Removing a site results in a break in the chain,
and the effective evolution then occurs within a smaller chain. With longer-
range couplings, freezing one site does not break the chain, positioning the
target site consistently within the bulk as shown in the bottom three
chains.111To prevent chain breakage upon freezing site $p$, one must ensure
couplings of range $p+1$ in a chain with open boundary conditions. For
periodic boundary conditions, there is a connection between sites $1$ and $L$,
preventing chain breakage but transforming the system into an open boundary
condition system upon freezing a site. However, we opt for open boundary
conditions in this work to avoid this additional complexity. Although it
should not significantly differ from the results in this work for unitary
evolution, discrepancies may arise in the non-unitary scenario, particularly
where potential skin effects emerge. We plan to explore these differences in a
subsequent study..

#### Localization transition in Aubry-Andre-Harper model

Figure 2: Cumulative Liang info growth for the model in Eq. (4). (a) We show
late time cumulative Liang information to a site $d$ sites away from the
frozen site. $\overline{|\mathbb{T}_{d}|}$ denotes averaged $\mathbb{T}_{d}$
for $t\sim 10^{2}$. Averaging is performed to smoothen the plot. The blue
dashed line denotes the parameter $\lambda=2$ where localization phase
transition occurs. (b) The cross section of (a) marked by the red dashed line
alongside similar data for other system sizes denoted by different colours.
The blue line denotes the critical point. The position of the frozen size is
the closest smaller Fibonacci number $+1$ to $L$, i.e. for $L=2584$ it is
$i=1598$ etc.

Liang information, as a measure of causality, proves to be a highly effective
tool for delineating the cone of influence in quantum systems. To illustrate,
we examine a 1D XX model with a spatially varying on-site magnetic field. The
corresponding Hamiltonian is expressed as follows:

$H=\sum_{j=1}^{L-1}(\sigma^{x}_{j}\sigma^{x}_{j+1}+\sigma^{y}_{j}\sigma^{y}_{j+1})+\sum_{j=1}^{L}\mathcal{B}_{j}\sigma^{z}_{j}.$
(4)

where, $\mathcal{B}_{j}=\lambda\cos(2\pi\beta j)$, and $\beta$ is typically
chosen as the inverse golden ratio, $\beta=\frac{\sqrt{5}-1}{2}$. The system
size $L$ is typically chosen as a Fibonacci number to minimize finite size
effects. This model undergoes a localization-delocalization transition at
$\lambda=2$ across its eigenspectrum [26].

In Fig. 2, we explore the causal influence of a selected site on other sites
across the transition, quantified by cumulative Liang information
$\mathbb{T}_{d}$ flow over time between two sites separated by a distance $d$.
The plot is smoothed through an average over $\mathbb{T}_{d}$ at timescales
beyond the initial transient growth ($t\sim 10^{2}$). While the frozen site
can be arbitrarily chosen for an infinite system, for finite systems as the
one utilized in this work, we select the truncated lattice length as the next
lowest Fibonacci number to avoid spurious finite size effects from the
incommensurate potential. 222To avoid cutting the chain one can consider a
power law decaying hopping model instead with large power $~{}\sim O(10)$,
such that the phase transition point remains near $\lambda=2$. We have
verified (not shown) that this does not affect the results qualitatively.

Figure 3: Variation of absolute value of Liang information for quench to
different $\mathcal{B}$ fields for $\kappa=0$. For (a) and (b) The initial
state is chosen to be ground state of $\mathcal{B}=0.01$ while for (c) and (d)
the initial state is taken to be the ground state of the $\mathcal{B}$ field
at which the Liang information is computed. Fig. (a) shows the value of Liang
information $|\mathbb{T}_{3}|$, $3$ sites away from the frozen site. The red
dashed lines indicate the cross section plotted in Fig. (b), in which the
black line indicates the critical point. In (c) and (d) even though we start
from an eigenstate, the non-trivial evolution occurs due to the removal of the
frozen site, and is used to compute the causation flow. (e) and (f) shows the
spatial profile of Liang information for quenches for $L=250$ at a chosen time
instance. (e) has same setup as (a),(b) and (f) has the same setup as (c),(d).
The red dashed line indicates the critical point $\mathcal{B}_{c}=1$

In Fig. 2(a), we observe that in the delocalized regime, there is equitable
Liang information flow from the frozen site to other sites. This corroborates
the unrestricted transport expected in this regime as a result of spatially
extensive single-particle eigenfunctions. However, as $\lambda$ approaches the
critical point, we note a significant increase in causation effects for nearby
sites compared to those farther away. This behavior stems from the restriction
of information propagation beyond a certain distance, the localization length.
Consequently, local effects strongly influence the evolution of a site. Given
that this phase transition occurs across the eigenspectrum, the energy of the
chosen initial state minimally affects our findings. Hence, our choice of the
Neel state as the initial state for this analysis yields qualitative results
similar to those obtained with other typical initial states (see Appendix A).

In Fig. 2(b), we focus on two specific scenarios: (i) $d=1$, representing a
nearby site to the frozen site, and (ii) $d=15$, representing a site at a
significantly greater distance beyond the localization length of the system.
For the site at $d=1$, we observe the anticipated behavior described in the
preceding paragraph as $\lambda$ increases. However for distant sites the
causation flow becomes interesting near criticality. For small $\lambda$, the
behavior mirrors that of nearby sites i.e. there is a gradual increase in
cumulative Liang information with $\lambda$ until approximately $\lambda\sim
1.5$. This is followed by a plateau before the expected sharp decline in the
localized regime. The intriguing aspect is the occurrence of a peak in
causality for $\lambda<2$ for distant sites, which seems to slowly drift
towards $\lambda\sim 1.5$ for larger values of $d$ as seen in Fig. 2(a), from
the faint violet regions. The larger causation values indicate an already
inequitable flow of information between different sites in the said parameter
regime. This in turn signifies traces of localization in the parts of the
system acting as a herald to the onset of localization across the spectrum for
$\lambda=2$. However, it is not fully clear whether this effect holds in the
thermodynamic limit, even though a naive finite size analysis in Fig. 2(b)
suggests it is not dependent on size. Furthermore, for $\lambda>2$
exponentially localized wavefunctions result in exponentially small
information leakage with $d$ beyond the localization length, whose signature
is given by a rapid but continuous decrease in Liang information with $d$.

#### Ground state phase transitions in Ising models:

The localization phase transition in the previous example demonstrated the
goodness of Liang information to measure causation, which encourages us to
study it in more complicated settings. So we next direct our focus to the
Anisotropic Next Nearest Neighbor Ising (ANNNI) model. This model includes a
next-nearest-neighbor coupling, is non-integrable, and lacks $U(1)$ symmetry
present in the previous example. In certain parameter regimes, the ground
state of this model undergoes a ferromagnetic to paramagnetic Ising phase
transition upon tuning the transverse magnetization strength. In the following
analysis, we will investigate the influence of the middle site of the chain on
other chosen sites across different parameter regimes. Generalizing this study
to encompass the influence over multiple sites is straightforward.

The Hamiltonian of ANNNI chain is given by,

$H_{L}=-\sum_{j}^{L-1}\sigma_{j}^{z}\sigma^{z}_{j+1}+\kappa\sum_{j}^{L-2}\sigma_{j}^{z}\sigma_{j+1}^{z}-\mathcal{B}\sum_{j}^{L}\sigma_{j}^{x},$
(5)

where $\kappa$ represents the strength of the next-nearest neighbor term and
$\mathcal{B}$ is the magnitude of magnetic field applied along the transverse
axis and $L$ denotes the number of sites in the system.

At fixed value of $0<\kappa<0.5$, the ground state of the ANNNI model
undergoes quantum phase transition from the ferromagnetic phase to the
paramagnetic phase when transverse field $\mathcal{B}>0$ exceeds a critical
value $\mathcal{B}_{c}$. The phase transition point has been been computed
before[28, 29, 30, 31, 32, 33] and the critical parameters $\kappa_{c}$ and
$\mathcal{B}_{c}$ are known to satisfy:

$1-2\kappa_{c}=\mathcal{B}_{c}-\mathcal{B}_{c}^{2}\frac{\kappa_{c}}{2-2\kappa_{c}}$
(6)

While the quantum phase transition is traditionally considered an equilibrium
phenomenon occurring in the ground state, previous studies [34, 20, 21, 35]
have demonstrated that signatures of this transition can manifest in non-
equilibrium quantum quenches. In the present work, motivated by the
understanding of heightened quantum sensitivity near phase transitions[36], we
actively seek these signatures within the causation measure Liang information.

Since the model in Eq. (5) is non-integrable for generic parameter values we
will resort to numerical simulations (exact diagonalization, TDVP and DMRG)
for our results. However at the special point $\kappa=0$, it reduces to the
integrable transverse Ising model. We first visit the growth of Liang
information in this integrable chain before discussing the results for
$\kappa>0$.

####  $\kappa=0$ case:

For $\kappa=0$, we can use a canonical transformation
($\sigma^{x}\rightarrow\sigma^{z},\sigma^{z}\rightarrow-\sigma^{x}$) and then
a Jordan-Wigner transformation to map the $H_{L}$ to a system of spinless
fermions,

$H_{L}=-\frac{1}{2}\sum_{i=1}^{L-1}[c_{i}^{\dagger}c_{i+1}+c_{i}^{\dagger}c_{i+1}^{\dagger}+{\rm
h.c.}]-\mathcal{B}\sum_{i=1}^{L}(c_{i}^{\dagger}c_{i}-\frac{1}{2})$ (7)

where $c_{i}[c_{i}^{\dagger}]$ is the spinless JW fermion annihilation
[creation] operator. Since we are interested in the influence on a single
site, we require the one site density matrix which in this special case is
given by[37],

$\rho_{j}=\frac{\mathbb{I}+\langle\sigma^{z}_{j}\rangle\sigma^{z}}{2},$ (8)

where $\langle\sigma^{z}_{j}\rangle=2\langle c_{j}^{\dagger}c_{j}\rangle-1$.
Furthermore since the Hamiltonian is quadratic in fermionic operators, the
time evolution of this quantity can be found semi-analytically. Details of the
computation are provided in Appendix. B

The results for such a system starting from different initial states are shown
in Fig. 3. In Fig. 3(a) we show the growth of Lianginfo three sites away from
the frozen site (chosen to be the middle site of the full lattice) when the
initial state is chosen to be the ground state at $\mathcal{B}=0.01$ and we
quench to different values of $\mathcal{B}$. When $\mathcal{B}$ is small, we
are far away from criticality. The initial state in this regime has large
overlap with the ground state, and the causal influence of one site over
another during evolution is very limited as evolution itself is limited. This
gradually changes as we increase $\mathcal{B}$ towards $\mathcal{B}_{c}$ as we
expect maximum sensitivity near the critical point. On the other extreme for
very large $\mathcal{B}$ the system has only local evolution, so we expect
Liang information to again be small during time evolution. This is the
qualitative behaviour we notice in Fig. 3(b). While the peak of causation is
distinct it appears to be quite shifted from the critical $\mathcal{B}_{c}=1$.
Naively we might identify the reason for this anomaly to be our choice of
initial state, which has finite overlaps several high energy states that do
not posses the desired features at criticality. To prevent this one can choose
the initial state to always be the ground state of the corresponding
$\mathcal{B}$, the results of which is shown in Figs. 3 (c) and (d). However
the peak still occurs for $\mathcal{B}_{p}<\mathcal{B}_{c}$, even though now
$\mathcal{B}_{p}$ shifts to larger values.

The strongest indicator is shown in Fig. 3(e) and (f) where we study how the
frozen site influences sites at different distances. Clearly for quenches at
$\mathcal{B}_{p}$, the chosen site influences other sites at all lengthscales.
This is a hallmark of criticality for second order phase transitions, the
divergence of correlation length, which in turn means correlations and thus
entanglement develop at all length scales. While typically such a behaviour is
seen in equilibrium properties of the ground state, this phenomenon still
shows up in Liang information, a non-equilibrium quantity, remarkably well
irrespective of the initial state.

Figure 4: Variation of absolute value of Liang information for quench to
different $\mathcal{B}$ field with constant $\kappa\neq 0$ for system size of
$L=50,250$ computed via TDVP444parameters- cutoff= $10^{-7}$, $dt=0.01$,
$\chi$ or bond dimension was capped at $5000$.. (a) Liang information
$|\mathbb{T}_{3}|$, $3$ sites away from the frozen site when the initial state
is chosen to be ferromagnetic state
$\ket{\downarrow\downarrow\ldots\downarrow}$, and $\kappa=0.4$. (b)
$|\mathbb{T}_{3}|$ for $\kappa=0.2$ starting from the same initial state. (c)
$|\mathbb{T}_{3}|$ for $\kappa=0.2$ starting from the ground state, similar to
Fig. 3(c). (d) Cross section of (c) at $t=30$ compared with $|\mathbb{T}_{3}|$
for a larger system size, $L=250$. The red/black dashed lines represent
$\mathcal{B}_{c}$ computed from Eq. (6).

####  $\kappa>0$ case:

For $\kappa>0$ one can perform quenches in $\mathcal{B}$ keeping $\kappa$
constant or vice-versa. Both yield qualitatively similar results (see also
Appendix D). In Fig. 4 we demonstrate that features shown by Liang information
in the non-interacting case continue to persist for the non-integrable
interacting model. As in the previous section we study two scenarios: the
ferromagnetic initial state
$\ket{\downarrow\downarrow\downarrow\ldots\downarrow}$ which is the ground
state at $\kappa=\mathcal{B}=0$, and the local quench from the ground state of
corresponding $\mathcal{B}$. We also add a very small longitudinal field to
prevent any degeneracy for small values of $\mathcal{B}$.

In Fig. 4(a) we notice that inspite of starting from the ferromagnetic initial
state for much smaller $L=50$ compared to Fig. 3, Liang information shows a
peak just before criticality. (in fact for $\kappa=0$, this small system size
shows noisy data). The reason for the closeness of the peak to the critical
point compared to the $\kappa=0$ case is that, larger $\kappa$ pushes the
critical point towards $\mathcal{B}=0$. Then even at these system sizes this
initial state maintains sufficient overlap with the ground state of the
Hamiltonian before the critical point to exhibit the expected behaviour.

But as we decrease $\kappa$ to $0.2$, there remains no correlation between the
maxima and the critical point, which is denoted by the red dashed lines. This
is expected as the critical point moves to larger $\mathcal{B}$, the initial
state no longer remains a low energy state and thus cannot respond to phase
transition in ground states. Furthermore since this model is non-integrable,
lack of conservation laws allows exploration of more of Hilbert space compared
to the $\kappa=0$ model leading to greater deviations than that case.

However, we do recover the peak just before criticality phenomena once we
start from the ground state, shown in Figs. 4(c),(d), just as we did for
$\kappa=0$, for even smaller system sizes.(See Appendix D for other
parameters). It is clear that for all $\kappa$ which supports phase
transition, the maxima typically occurs at
$\mathcal{B}_{p}\neq\mathcal{B}_{c}(\mathcal{B}<\mathcal{B}_{c}$ for all the
cases considered here), distinguishing causation from correlation. While
correlation length of the ground state can indeed diverge at the critical
point, the causation dynamics can maximize at slightly different values, as it
is governed by the competing behaviour of the various Hamiltonian terms on the
initial state. Of course the correlation lengthscales still plays a crucial
part in causation dynamics, more so in the integrable model and for the ground
state initial state where such non-analyticities are stronger, shown by the
closeness in$\mathcal{B}_{p}$ and $\mathcal{B}_{c}$ in these cases. Further
analytical insights may be obtained by finding a mean field mapping of
$\kappa\neq 0$ model to the solvable transverse Ising model[34], but it is
left as a future endeavour.

#### Discussion:

In this work we have studied causation behaviour of different sites in quantum
chains by studying the cumulative Liang information flow between them.
Specifically, we study the parameter regime where we find a peak in cumulative
Liang information indicating a causality peak. The causality peaks do not
necessarily correspond to a maxima in correlation of the corresponding sites,
but rather roughly detects the maxima in difference of expectation values of
local operators with and without the site whose influence we want to compute.
We have showed that this measure successfully showcases the democracy of
influence in the fully delocalized regime and near-site causation in localized
regime in one dimensional Aubry Andre Harper model. Furthermore it even
responds to the diverging correlation lengths in ground state second order
phase transitions of the next-nearest neighbour Ising model. Our most crucial
find is that the peak of causation does not occur exactly at the critical
point, but rather slightly towards the ordered side of the transition, for
quenches from the ordered side. Thus causation flow encapsulated by quantum
Liang information between two sites is a new non-equilibrium herald of an
equilibrium phase transition.

While classical Liang information analysis is a well-established approach to
quantify causation in classical networks, its applicability in quantum systems
is an area that is just beginning to be explored. This work represents an
initial exploration, focusing on two specific quantum systems. One can
immediately see that Liang information analysis can provide a natural test bed
for the analysis of many-body localization transitions, for instance, it can
be employed to identify resonances between sets of sites that induce the
avalanche breakdown of localization in different geometries. There are also
open questions regarding its behaviour in the vicinity of topological phase
transitions, and in higher dimensional quantum systems including complex
networks.

###### Acknowledgements.

RG thanks Krishnendu Sengupta, Arnab Das and Madhumita Sarkar for discussions.
RG and SB acknowledge the EPSRC grants Nonergodic quantum manipulation
EP/R029075/1 and Many-Body Phases In Continuous-Time Quantum Computation
EP/Y004590/1 for support .

## References

  * Polkovnikov  _et al._ [2011] A. Polkovnikov, K. Sengupta, A. Silva, and M. Vengalattore, Colloquium: Nonequilibrium dynamics of closed interacting quantum systems, Rev. Mod. Phys. 83, 863 (2011).
  * Brukner [2014] Č. Brukner, Quantum causality, Nature Physics 10, 259 (2014).
  * Maldacena  _et al._ [2016] J. Maldacena, S. H. Shenker, and D. Stanford, A bound on chaos, Journal of High Energy Physics 2016, 106 (2016).
  * Kukuljan  _et al._ [2017] I. Kukuljan, S. Grozdanov, and P. Tomaž, Weak quantum chaos, Phys. Rev. B 96, 060301 (2017).
  * Harrow  _et al._ [2021] A. W. Harrow, L. Kong, Z.-W. Liu, S. Mehraban, and P. W. Shor, Separation of out-of-time-ordered correlation and entanglement, PRX Quantum 2, 020339 (2021).
  * Riddell and Sørensen [2019] J. Riddell and E. S. Sørensen, Out-of-time ordered correlators and entanglement growth in the random-field xx spin chain, Phys. Rev. B 99, 054205 (2019).
  * Notenson  _et al._ [2023] T. Notenson, I. García-Mata, A. J. Roncaglia, and D. A. Wisniacki, Classical approach to equilibrium of out-of-time ordered correlators in mixed systems, Phys. Rev. E 107, 064207 (2023).
  * Xu  _et al._ [2020] T. Xu, T. Scaffidi, and X. Cao, Does scrambling equal chaos?, Phys. Rev. Lett. 124, 140602 (2020).
  * Pandey  _et al._ [2020] M. Pandey, P. W. Claeys, D. K. Campbell, A. Polkovnikov, and D. Sels, Adiabatic eigenstate deformations as a sensitive probe for quantum chaos, Phys. Rev. X 10, 041017 (2020).
  * Lim  _et al._ [2024] C. Lim, K. Matirko, A. Polkovnikov, and M. O. Flynn, Defining classical and quantum chaos through adiabatic transformations (2024), arXiv:2401.01927 [cond-mat.stat-mech] .
  * San Liang [2016] X. San Liang, Information flow and causality as rigorous notions ab initio, Physical Review E 94, 052201 (2016).
  * Liang [2018] X. S. Liang, Causation and information flow with respect to relative entropy, Chaos: An Interdisciplinary Journal of Nonlinear Science 28, 075311 (2018).
  * Schneider and Griffies [1999] T. Schneider and S. M. Griffies, A conceptual framework for predictability studies, Journal of climate 12, 3133 (1999).
  * Kleeman [2002] R. Kleeman, Measuring dynamical prediction utility using relative entropy, Journal of the atmospheric sciences 59, 2057 (2002).
  * Kleeman and Majda [2005] R. Kleeman and A. Majda, Predictability in a model of geophysical turbulence, Journal of the atmospheric sciences 62, 2864 (2005).
  * Kleeman [2008] R. Kleeman, Limits, variability, and general behavior of statistical predictability of the midlatitude atmosphere, Journal of the atmospheric sciences 65, 263 (2008).
  * Yi and Bose [2022] B. Yi and S. Bose, Quantum liang information flow as causation quantifier, Physical Review Letters 129, 020501 (2022).
  * Sengupta  _et al._ [2004] K. Sengupta, S. Powell, and S. Sachdev, Quench dynamics across quantum critical points, Phys. Rev. A 69, 053616 (2004).
  * Basu and Das [2012] P. Basu and S. R. Das, Quantum quench across a holographic critical point, Journal of High Energy Physics 2012, 103 (2012).
  * Haldar  _et al._ [2021] A. Haldar, K. Mallayya, M. Heyl, F. Pollmann, M. Rigol, and A. Das, Signatures of quantum phase transitions after quenches in quantum chaotic one-dimensional systems, Phys. Rev. X 11, 031062 (2021).
  * Dağ  _et al._ [2023a] C. B. Dağ, P. Uhrich, Y. Wang, I. P. McCulloch, and J. C. Halimeh, Detecting quantum phase transitions in the quasistationary regime of ising chains, Phys. Rev. B 107, 094432 (2023a).
  * Shen  _et al._ [2017] H. Shen, P. Zhang, R. Fan, and H. Zhai, Out-of-time-order correlation at a quantum phase transition, Phys. Rev. B 96, 054503 (2017).
  * Sun  _et al._ [2020] Z.-H. Sun, J.-Q. Cai, Q.-C. Tang, Y. Hu, and H. Fan, Out-of-Time-Order Correlators and Quantum Phase Transitions in the Rabi and Dicke Models, Annalen Phys. 532, 1900270 (2020).
  * Heyl  _et al._ [2018] M. Heyl, F. Pollmann, and B. Dóra, Detecting equilibrium and dynamical quantum phase transitions in ising chains via out-of-time-ordered correlators, Phys. Rev. Lett. 121, 016801 (2018).
  * Note [1] To prevent chain breakage upon freezing site $p$, one must ensure couplings of range $p+1$ in a chain with open boundary conditions. For periodic boundary conditions, there is a connection between sites $1$ and $L$, preventing chain breakage but transforming the system into an open boundary condition system upon freezing a site. However, we opt for open boundary conditions in this work to avoid this additional complexity. Although it should not significantly differ from the results in this work for unitary evolution, discrepancies may arise in the non-unitary scenario, particularly where potential skin effects emerge. We plan to explore these differences in a subsequent study.
  * Aubry and André [1980] S. Aubry and G. André, Analyticity breaking and anderson localization in incommensurate lattices, (1980).
  * Note [2] To avoid cutting the chain one can consider a power law decaying hopping model instead with large power $~{}\sim O(10)$, such that the phase transition point remains near $\lambda=2$. We have verified (not shown) that this does not affect the results qualitatively.
  * Selke [1988] W. Selke, The annni model — theoretical analysis and experimental application, Physics Reports 170, 213 (1988).
  * Suzuki  _et al._ [2013] S. Suzuki, J.-i. Inoue, and B. K. Chakrabarti, Annni model in transverse field, in  _Quantum Ising Phases and Transitions in Transverse Ising Models_ (Springer Berlin Heidelberg, Berlin, Heidelberg, 2013) pp. 73–103.
  * Allen  _et al._ [2001] D. Allen, P. Azaria, and P. Lecheminant, A two-leg quantum ising ladder: a bosonization study of the annni model, Journal of Physics A: Mathematical and General 34, L305 (2001).
  * Nagy [2011] A. Nagy, Exploring phase transitions by finite-entanglement scaling of mps in the 1d annni model, New Journal of Physics 13, 023015 (2011).
  * Guimarães  _et al._ [2002] P. R. C. Guimarães, J. a. A. Plascak, F. C. Sá Barreto, and J. a. Florencio, Quantum phase transitions in the one-dimensional transverse ising model with second-neighbor interactions, Phys. Rev. B 66, 064413 (2002).
  * Beccaria  _et al._ [2006] M. Beccaria, M. Campostrini, and A. Feo, Density-matrix renormalization-group study of the disorder line in the quantum axial next-nearest-neighbor ising model, Phys. Rev. B 73, 052402 (2006).
  * Robertson  _et al._ [2023] J. H. Robertson, R. Senese, and F. H. L. Essler, A simple theory for quantum quenches in the ANNNI model, SciPost Phys. 15, 032 (2023).
  * Dağ  _et al._ [2023b] C. B. Dağ, Y. Wang, P. Uhrich, X. Na, and J. C. Halimeh, Critical slowing down in sudden quench dynamics, Phys. Rev. B 107, L121113 (2023b).
  * Zanardi and Paunković [2006] P. Zanardi and N. Paunković, Ground state overlap and quantum phase transitions, Phys. Rev. E 74, 031123 (2006).
  * Osborne and Nielsen [2002] T. J. Osborne and M. A. Nielsen, Entanglement in a simple quantum phase transition, Phys. Rev. A 66, 032110 (2002).
  * Ghosh and Das [2021] R. Ghosh and A. Das, Disorder-induced enhancement of entanglement growth in one dimension: Information leakage at the scale of the localization length, Phys. Rev. B 103, 024202 (2021).
  * Peschel [2003] I. Peschel, Calculation of reduced density matrices from correlation functions, Journal of Physics A: Mathematical and General 36, L205 (2003).
  * Peschel and Eisler [2009] I. Peschel and V. Eisler, Reduced density matrices and entanglement entropy in free lattice models, Journal of Physics A: Mathematical and Theoretical 42, 504003 (2009).
  * Lieb  _et al._ [1961] E. Lieb, T. Schultz, and D. Mattis, Two soluble models of an antiferromagnetic chain, Annals of Physics 16, 407 (1961).
  * Kormos [2017] M. Kormos, Inhomogeneous quenches in the transverse field Ising chain: scaling and front dynamics, SciPost Phys. 3, 020 (2017).
  * Note [3] One can also rewrite them in terms of Majorana fermions and perform the same computation.

##  Appendix A Initial state independence of Aubry-Andre model results

In Fig, 5 we choose two random computational initial states for system size
$L=987$ at approximately half filling $\sum_{i}\sigma^{z}_{i}=1$ and show the
behaviour of smoothed time averaged Liang information
$\overline{\mathbb{T}_{d}}$ with distance $d$ from the frozen site. Evidently,
the qualitative nature of all the cases are the same, but there are a few
quantitative differences in causation flow based on the initial state, which
is picked up by Liang information. Specifically near the localization
transition $\lambda\sim 1.5-2$, the loss of democracy in causality for
different initial states differs slightly from each other. Indeed, all initial
states cannot be expected to relax exactly equivalently, for example,
atyipical states such as the domain wall states have very different
information leakage during evolution compared to Neel states.[38] This shows
up as small quantitative differences in flow of Liang information for states
at intermediate distance from the frozen site as seen in the figure.

Figure 5: Demonstration of cumulative Liang information flow for Aubry Andre
model, similar to Fig. 2

##  Appendix B Deatils of computation for open Ising chain

The nearest neighbour transverse Ising model Hamiltonian,

$H=-\sum_{j}^{L-1}\sigma_{j}^{z}\sigma^{z}_{j+1}-\mathcal{B}\sum_{j}^{L}\sigma_{j}^{x},$
(9)

after the JW transformation is a Hamiltonian quadratic in fermionic operators.

$H=-\frac{1}{2}\sum_{i=1}^{L-1}[c_{i}^{\dagger}c_{i+1}+c_{i}^{\dagger}c_{i+1}^{\dagger}+{\rm
h.c.}]-\mathcal{B}\sum_{i=1}^{L}(c_{i}^{\dagger}c_{i}-\frac{1}{2})$ (10)

Hence the entanglement entropy of a subsystem can be computed just from the
two point correlation functions[39, 40] $\langle c_{i}^{\dagger}c_{j}\rangle$
and $\langle c_{i}^{\dagger}c_{j}^{\dagger}\rangle$. In our case this further
simplifies to just computing $\langle c_{j}^{\dagger}c_{j}\rangle$ since we
are just interested in computing influence on one site i.e. the subsystem
consists of just one site. In fact it can be shown that[37] the one site
density matrix in such a system is given by,

$\rho_{j}=\frac{\mathbb{I}+\langle\sigma^{z}_{j}\rangle\sigma^{z}}{2},$ (11)

where $\langle\sigma^{z}_{j}\rangle=2\langle c_{j}^{\dagger}c_{j}\rangle-1$.
Since we are using open boundary conditions, momentum $k$ is no longer a good
quantum number, hence this Hamiltonian cannot be block diagonalized by the
Fourier transformation. However using the quadratic nature of the Hamiltonian
we can construct the eigenmodes of the system analytically as follows[41, 42].

We define two new fermionic operators555One can also rewrite them in terms of
Majorana fermions and perform the same computation for convenience,

$A_{i}=(c_{i}^{\dagger}+c_{i}),\hskip
14.45377ptB_{i}=(c_{i}^{\dagger}-c_{i}),$ (12)

the Eq. (5) becomes,

$H=\frac{1}{2}\sum_{i}^{L-1}A_{i+1}B+\frac{h}{2}\sum_{i}A_{i}B_{i}$ (13)

and we can compute the two point correlations functions such as
$\langle\sigma^{z}_{j}\rangle=-\langle A_{j}B_{j}\rangle$ appropriately.

Then $H$ in Eq. 10 can be immediately diagonalized to

$H_{D}=\sum_{j}\epsilon_{k}\eta^{\dagger}_{k}\eta_{k}+{\rm const.}$ (14)

where $\epsilon_{k}=\sqrt{1-2\mathcal{B}\cos k+\mathcal{B}^{2}}$, by the real
linear transformation,

$\eta_{k}=\frac{1}{2}\sum_{j}[\phi_{kj}A_{j}-\psi_{kj}B_{j}]$ (15)

with

$\begin{split}\phi_{kj}&=-\mathcal{N}_{k}(-1)^{n}\sin[k(L+1-j)]\\\
\psi_{kj}&=-\mathcal{N}_{k}\sin[kj].\end{split}$ (16)

Here the $k\in[0,\pi]$ modes are obtained from the $n^{{\rm th}}$ (This is the
same $n$ as in Eq. (16)) solution of the equation,

$\frac{\sin[(L+1)k]}{\sin[Lk]}=-\frac{1}{\mathcal{B}}$ (17)

and $\mathcal{N}_{k}$ is the normalization constant given by,

$\mathcal{N}_{k}=\frac{1}{\sqrt{\frac{L}{2}+\frac{\mathcal{B}(\mathcal{B}+\cos
k)}{2\epsilon_{k}^{2}}}}.$

Note that for $\mathcal{B}\rightarrow\infty$, $k\sim\frac{\pi
n}{L+1},n={1,..,L}$. Note that for $\mathcal{B}>1$ Eq. (17) always supports
real valued solutions, but it can have an imaginary solution for
$\mathcal{B}<1$. Clearly the imaginary value of $k$ indicates the wavefunction
is localized, and it happens at the edges. For further details on construction
of this localized solution refer 42. This aspect is not relevant for our
computations. Now we shall compute the time evolution of $\langle
c_{i}^{\dagger}c_{i}\rangle$ to understand the time evolution of the reduced
density matrix. This can be done in the following steps.

  1. 1.

$\langle
A_{j}(t)B_{j}(t)\rangle$=$\langle\cancel{(c_{j}^{\dagger})^{2}}+1-2c_{j}^{\dagger}(t)c_{j}(t)-\cancel{(c_{j})^{2}}\rangle=-\langle\sigma^{z}_{j}(t)\rangle$
where we have used that $c_{j}[c_{j}^{\dagger}]$ are fermionic annihilation
[creation] operators.

  2. 2.

$\langle A_{j}(t)B_{j}(t)\rangle$ can be written in terms of the Bogoliubov
quasiparticles,

$\langle
A_{j}(t)B_{j}(t)\rangle=\sum_{k,l=1}^{L}\phi_{kj}\psi_{lj}\langle(\eta_{k}^{\dagger}(t)+\eta_{k}(t))(\eta_{l}(t)^{\dagger}-\eta_{l}(t))\rangle$

  3. 3.

Now $\eta(t)$ and $\eta^{\dagger}$ being the quasiparticles of the Hamiltonian
has a very simple time evolution,

$\eta_{k}^{\dagger}(t)=e^{i\epsilon_{k}t}\eta_{k}^{\dagger}(0)$

which can be substituted in the previous expression.

  4. 4.

The next step is to rotate the expression back to thye $A,B$ basis to obtain

$\begin{split}\langle
A_{j}(t)B_{j}(t)\rangle&=-\frac{1}{2}\sum_{m,n,k.l}\phi_{kj}\psi_{lj}\phi_{kn}\psi_{lm}\langle
A_{n}(0)B_{m}(0)\rangle(\cos[t(\epsilon_{k}+\epsilon_{l})]+\cos[t(\epsilon_{k}-\epsilon_{l})])\\\
&-\frac{1}{2}\sum_{m,n,k.l}\phi_{kj}\psi_{lj}\psi_{kn}\phi_{lm}\langle
B_{n}(0)A_{m}(0)\rangle(\cos[t(\epsilon_{k}+\epsilon_{l})]-\cos[t(\epsilon_{k}-\epsilon_{l})])\end{split}$
(18)

Here while indices $j,m,n$ can take any integer values in $[1,L]$, $k,l$ label
the eigenmodes obtained by solving Eq. (17). In our case (see Fig.1), since we
remove one site for the system, this amounts to a change in $L$ and $j$ to
$L^{\prime}$ and $j^{\prime}$ of the system. This means upon removal of site
$\langle A_{j}(0)B_{m}(0)\rangle$ will require the available indices
restricted to $L^{\prime}$ with $j$ changing to $j^{\prime}$ reflecting the
new position of the same site on the truncated chain.

##  Appendix C Lightcones in the ANNNI model

Figure 6: (a) The variation of the magnitude of Liang information at a
distance $d$ from the frozen site. $d=0$ denotes the site next to the frozen
site for a quench from the ferromagnetic ground state
$\ket{\downarrow\downarrow\ldots\downarrow}$ to $\mathcal{B}=0.1$ and
$\kappa=0.4$. (b) Plot of $\tau$, the time at which $\mathbb{T}_{d}\geqslant
0.001$ for quenches to different values of $\mathcal{B}$ for $\kappa=0.4$. the
system size chosen is $L=50$.

Let us examine the causality light cones obtained for the next-nearest
neighbor ANNNI model. For this purpose, we choose representative values of
$\mathcal{B}$ and $\kappa$ in Fig. 6 for the quench and study the evolution of
Liang information across different sites in the chain. The initial state is
chosen to be the ferromagnetic state,
$\ket{\downarrow\downarrow\ldots\downarrow}$ As anticipated, on the site
closest to the frozen site ($d=0$), we observe an immediate growth of Liang
information, with the growth occurring at later times as we move to sites with
increasing distance. This characteristic pattern is a hallmark of local
Hamiltonians. Furthermore, there is larger causation flow for
$\mathcal{B}=0.2$ than for $\mathcal{B}=0.1$, since $\mathcal{B}=0.2$ is
closer to the critical point $\mathcal{B}_{c}=0.21$. This is specially
apparent for the $d=3,6$ cases where $|\mathbb{T}_{d}|$ shows values which
differ by several orders of magnitude for the two cases.

In Fig. 6(b) we define $\tau$ such that $|\mathbb{T}_{d}|<0.001$ for $t<\tau$
and plot its variation with $\mathcal{B}$ and distance $d$. We observe an
interesting behaviour that $\tau$ for same $d$ becomes smaller as we increase
$\mathcal{B}$. We have also verified that the even smaller choices of cutoff
does not qualitatively change our results. This result can seem surprising
since we have claimed causation peak is observed for quenches to $\mathcal{B}$
which is just smaller than $\mathcal{B}_{c}$, but we see a decrease in
$\tau(d)$ with $\mathcal{B}$, which seems to suggest the influence of sites
maximize beyond critical quenches. But upon closer look we realize that this
is just a transient effect which is also seen in Fig. 4. Indeed since the
initial state is an eigenstate for $\mathcal{B}=0$, in the initial stages of
quench to intermediate $\mathcal{B}$, the initial undergoes non-trivial
evolution which reflects in the significant influence of different sites on
the chosen site during evolution. Relaxation typically occurs faster for
larger $\mathcal{B}$, which in turn indicates smaller observed values of
$\tau$. But beyond the transient relaxation period the sites exert very little
influence, as the system has completely relaxed. And furthermore large
$\mathcal{B}$ makes the Hamiltonian effectively on-site, so each spin
effectively evolves independently. On the opposite extreme, for small
$\mathcal{B}$, where the initial state has a very large overlap with one of
the eigenstates of the Hamiltonian for both the ordinary and frozen system.
Hence the evolution is slow, and the slope is larger.

##  Appendix D Results for quenches to other parameters

Figure 7: Variation of absolute value of Liang information for quench to
different $\kappa$ fields with constant $B=0.2$. The initial state is chosen
to be ferromagnetic state, i.e. eigenstate of $\mathcal{B}=\kappa=0$. Fig. (a)
shows the value of Liang information $|\mathbb{T}_{3}|$, $3$ sites away from
the frozen site. the black dashed lines indicate the cross sections plotted in
Fig.(c), the red dashed line indicates the cross section plotted in Fig.(d).
Fig.(b) is same as Fig.(a) except the distance is $6$ sites.

###  D.1 Quench in $\kappa$ from ferromagnetic state

In this section, we discuss the effect of the nature of the quench by studying
the variation of Liang information flow while quenching to a specific
$\mathcal{B}$ but varying $\kappa$. The results are shown in Fig. 7, where we
quench to $\mathcal{B}=0.2$ with varying $\kappa$. From Eq. 6, we expect
$\kappa_{c}\sim 0.4$. Hence, using similar arguments as before, we expect
Liang information to show a sharp rise close to $\kappa_{c}$ before showing a
sharp decrease for $\kappa>\kappa_{c}$, a feature evident in Fig. 7(a), (b),
and (d). As expected from the previous section, the liang information of the
farther site, $|\mathbb{T}_{6}|$ shows slightly lower sensitivity to the
critical point than $|\mathbb{T}_{3}|$. Also, one can clearly see the light-
cone effect, which causes the Liang info to rise at later times for the more
distant site.

Figure 8: Variation Liang formation when the quench is for $\kappa=0.45$. (a)
The initial state is the ferromagnetic state. All the rest of the parameters
are same as Fig. 4(a). (b) The initial state is the ground state of the
corresponding $\mathcal{B}$, the rest of the parameters are same as Fig. 4(c)

In Fig. 7(c), we plot the growth for three values of $\kappa$ and see that
remarkably, for values of $\kappa$ close but less than $\kappa_{c}$, we almost
saturate the Liang information growth, indicating that the entanglement growth
becomes fully different between the unfrozen site and frozen site scenarios.
This stems from the still large overlap of the initial state with the ground
state eigenfunction in the unfrozen case, while a significant change occurs in
the frozen case. Comparing to the quench in $\mathcal{B}$, changing $\kappa$
while keeping $\mathcal{B}$ constant causes the evolution of the site-frozen
system to change more than changing $\mathcal{B}$ just before the critical
point. This is possibly because while $\mathcal{B}$ is just an on-site term,
$\kappa$ connects two sites and thus shows an enhanced response when Liang
information is computed.

Figure 9: Variation of Liang formation for the quench starting from the ground
state for top-$\kappa=0.1$ , middle=$\kappa=0.3$ bottom-$\kappa=0.4$. The rest
of the parameters are same as Fig. 4(c)

###  D.2 Further results for quenches in $\mathcal{B}$

In Fig. 8 we plot the scenario for quenches in $\mathcal{B}$ when
$\kappa=0.45$. In (a) we show the case when we start from the ferromagnetic
state, and since the critical point is for even smaller $\mathcal{B}$ than
when $\kappa=0.4$, such an initial state gives good response near the critical
point, comparable to the case when we start from the ground state of
corresponding $\mathcal{B}$ which is shown in (b). Interestingly it seems that
here starting from the ground state gives a weaker signal than starting from
the ferromagnetic state. This is possibly due to the $\mathcal{B}_{c}\sim 0$
making the on site spin flip term of weaker strength. This evolution very slow
in this case and we possibly need to go to much larger timescales to have a
conclusive result which is beyond the scope of ordinary tensor network
algorithms. Finally in Fig. 9 we show more cases of the ‘peak before critical
point’ phenomena computing the Liang information flow between sites separated
by a distance $d=3$ starting from the ground state. These results further
corroborate our conclusions in the main text. Additionally note how the growth
of Liang information near the point of criticality occurs at an earlier
timescale for smaller $\kappa$, which showcases how larger $\mathcal{B}$,
induces spin flip processes at a smaller timescale causing dynamics, if
existent , to speed up, creating a light cone with a greater slope.