# Universality in the tripartite information after global quenches:
spin flip and semilocal charges

Vanja Marić  Université Paris-Saclay, CNRS, LPTMS, 91405, Orsay, France.

###### Abstract

We study stationary states emerging after global quenches in which the time
evolution is under local Hamiltonians that possess semilocal conserved
operators. In particular, we study a model that is dual to quantum XY chain.
We show that a localized perturbation in the initial state can turn an
exponential decay of spatial correlations in the stationary state into an
algebraic decay. We investigate the consequences on the behavior of the
(Rényi-$\alpha$) entanglement entropies, focusing on the tripartite
information of three adjacent subsystems. In the limit of large subsystems, we
show that in the stationary state with the algebraic decay of correlations the
tripartite information exhibits a non-zero value with a universal dependency
on the cross ratio, while it vanishes in the stationary state with the
exponential decay of correlations.

##  I Introduction

This is the third paper in a series of works about the tripartite information
in the stationary states after global quenches, the first two being ref. [1]
and ref. [2], that will be in the following referred to as Paper I and Paper
II respectively. Paper II proves the results announced in Paper I concerning
quantum quenches from ground states of critical Hamiltonians and
bipartitioning protocols. The present paper completes the proofs of the
results announced in Paper I and deals with quenches in which semilocal
conservation laws affect the dynamics.

Quantum quench is a protocol in which the system is prepared in the ground
state of some local Hamiltonian and is suddenly let to evolve with a different
local Hamiltonian. It is perhaps the simplest way to induce non-equilibrium
dynamics and as such it has been thoroughly investigated [3, 4, 5]. In global
quenches, in particular, the two Hamiltonians are macroscopically different.
It has been established that in global quenches quite general isolated quantum
many-body systems locally relax to thermal states, as explained by the
eigenstate thermalization hypothesis [6, 7, 8, 9]. Integrable systems are an
exception [10]. Their stationary properties are instead captured by
generalized Gibbs ensembles (GGE) [11, 12, 13, 14], which carry memory of
additional conserved charges.

While most of the established results on global quenches concern
translationally invariant systems, there has been also a lot of effort to
relax the assumption of translational invariance. Bipartitioning protocols
[15] have been studied substantially, in which the initial state consists of
two macroscopically different parts, for example a domain wall. The theory of
generalized hydrodynamics [16, 17, 18] (GHD) explains that in bipartitioning
protocols the system can at late times be described in terms of locally-quasi-
stationary states, that depend on the ratio of the distance from the
inhomogeneity and the time. At infinite time (in the thermodynamic limit)
around the initial inhomogeneity, in particular, the system relaxes to a
stationary state, usually called non-equilibrium stationary state (NESS).

One of the standard tools to capture the universal properties of systems are
the entanglement entropies. At criticality the entanglement entropy of a
connected block exhibits a simple logarithmic scaling with the subsystem size,
where the prefactor is proportional to the central charge of the underlying
conformal field theory (CFT) [19, 20, 21, 22, 23]. The entanglement entropies
have been thoroughly studied also after global quenches and it has been
established that they saturate to a value that is extensive with the size of
the subsystem [24, 25, 26, 27, 28, 29, 30, 31, 32]. In Papers I and II it has
been pointed out that the entanglement entropy of a connected block after a
global quench can exhibit, both in translationally invariant quench protocols
and in NESS, besides the extensive term, also a subdominant term that grows
logarithmically with the subsystem size (see eq. (6)), similarly to the
universal leading term in CFT. Such subextensive logarithmic terms have been
found also in other works on NESS [33, 34, 35].

Because of the similarity of the subleading logarithmic term to the universal
leading term in CFT, it is desirable to consider a quantity that removes the
extensive and the boundary contributions of the entropies, leaving a
potentially universal quantity. In Papers I and II a special linear
combination of the entropies with such a property, namely the tripartite
information [36] of three adjacent blocks, has been considered for a wide
class of quench protocols in non-interacting spin chains. It has been shown
that quenches from critical points and bipartitioning protocols can exhibit
stationary states with a non-zero tripartite information, which is accompanied
by the subdominant logarithmic term in the entanglement entropy and the
existence of spatial correlation functions that decay algebraically, similarly
to CFT. In 1+1-dimensional CFT the tripartite information of three adjacent
blocks is a model-dependent function with a universal dependency on the cross
ratio [23, 37]. We note that the term “universal” here does not mean the
independence of the parameters of the model, but refers to the property that
the tripartite information of three adjacent blocks depends on their lengths
only through the cross ratio and is independent of the further details of the
configuration. On the lattice this universal property is achieved for large
subsystem lengths. Remarkably, tripartite information in the stationary states
of the aforementioned global quench protocols also exhibits this universality.
A property making it different from CFT is the existence of a nonzero
“residual tripartite information”, introduced in Paper I.

Localized perturbations of initial states are normally not expected to have
any macroscopic effects at large times. For systems of non-interacting
fermions it has been proven under quite general conditions that the
equilibration towards a GGE is resilient to localized perturbations, given
that the initial state has a finite correlation length and that the evolution
Hamiltonian is translationally invariant [38]. However, memory effects can
arise following a local quench in quantum spin chains, related to the non-
locality of the mapping between spins and fermions [39, 40, 41, 42, 43] and to
jammed states [44, 45]. Recently it has been found that the system can keep
the memory of localized perturbations even following a global quench, if it
possesses semilocal conserved operators [46]. In such a system a single spin
flip in the initial product state can induce a subextensive logarithmic term
in the entropy after a global quench [35], which hints at a possible
nontrivial tripartite information.

Semilocal charges are conserved operators whose density does not necessarily
commute with distant localized operators. Sometimes their densities can be
interpreted as semi-inifinte strings. Semilocal charges enable the existence
of string order after global quenches and, in general, they have to be
included in the GGE to capture correctly the properties of the stationary
state [35]. Furthermore, in such systems a single localized perturbation in
the initial state can change the magnetization in the stationary state [46]
and induce time growth of macroscopic entanglement [47]. Similar semilocal
symmetries have also been discussed in the context of bistability of driven-
dissipative fermionic systems [48].

Here we study a system with semilocal charges and show that in such a system a
single spin flip in the initial product state can turn an exponential decay of
spatial correlations in the stationary state into an algebraic decay. More
precisely, we show that, while the stationary state with an exponential decay
of spatial correlations emerges in the quench from a translationally invariant
initial product state, flipping a single spin in the initial state results in
NESS (corresponding to the ray $d/t=0$, where $d$ is the distance from the
initial inhomogeneity and $t$ is time) with an algebraic decay of spatial
correlations. Even though it is possible that the two quench protocols exhibit
different time scales on which the stationary states are reached, we stress
that the inifinite time limit is well defined in both cases. Then we proceed
to studying the consequences on the behavior of the entanglement entropies,
focusing on the tripartite information. We show that NESS arising from the
state with a single spin flip displays the behavior of the tripartite
information found in bipartitioning protocols and quenches from ground states
of critical systems.

The paper is organized as follows. In the remainder of the introduction we
discuss the tripartite information (section I.1) and we introduce the model
and the quench protocols under consideration (section I.2). In section II we
present the results of the paper. Their derivation is presented afterwards, in
section III. Conclusions are drawn in section IV.

###  I.1 Tripartite Information

$A$$B$$C$ Figure 1: We compute the tripartite information of adjacent blocks
$A$, $B$ and $C$ in stationary states after global quenches.

Given a subsystem $A$ described by a reduced density matrix $\rho_{A}$, the
von Neumman entanglement entropy is defined as

$S_{1}(A)=-\mathrm{tr}\left(\rho_{A}\log\rho_{A}\right).$ (1)

It corresponds to the limit $\alpha\to 1$ of the Rényi entanglement entropies

$S_{\alpha}(A)=\frac{1}{1-\alpha}\log\mathrm{tr}\left(\rho_{A}^{\alpha}\right).$
(2)

Typically, Rényi entropies for $\alpha=2,3,\ldots$ are more accessible to
computations than the von Neumann entropy and sometimes the latter can be
obtained from the former using the replica trick [20].

In pure states the entanglement entropies are a measure of quantum
entanglement. In mixed states, such as thermal states and GGEs, this is not
anymore the case, as the entropies carry contributions that are extensive with
the size of the subsystem and that are largely due to classical correlations.
Such extensive contributions are cancelled in the (Rényi-$\alpha$) mutual
information. Given two subsystems, $A$ and $C$, the mutual information [36] is
defined as

$I_{2}^{(\alpha)}(A,C)=S_{\alpha}(A)+S_{\alpha}(C)-S_{\alpha}(AC)\;.$ (3)

Here and in the rest of the paper $AC$ stands for the union $A\cup C$ of two
sets $A,C$. For large disjoint blocks $A,C$ in the configuration of figure 1
mutual information cancels both the extensive and the boundary contributions
of the entropies. Mutual information for $\alpha=1$ is non-negative and is a
measure of total correlations, classical and quantum, between $A$ and $C$
[49]. Moreover, mutual information provides an upper bound for connected
correlation functions [50]. We mention that when $A$ and $C$ constitute the
whole lattice, one being the complement of the other, there are area laws [50,
51, 52, 53, 54, 55] that render mutual information finite in thermal states.

Another quantity that cancels the extensive and the boundary contributions of
the entropies is the (Rényi $\alpha$) tripartite information [36], defined for
three subsystems $A,B,C$ as

$I_{3}^{(\alpha)}(A,B,C)=S_{\alpha}(A)+S_{\alpha}(B)+S_{\alpha}(C)-S_{\alpha}(AB)-S_{\alpha}(AC)-S_{\alpha}(BC)+S_{\alpha}(ABC)\;,$
(4)

which can be expressed in terms of the mutual information as

$I_{3}^{(\alpha)}(A,B,C)=I_{2}^{(\alpha)}(A,B)+I_{2}^{(\alpha)}(A,C)-I_{2}^{(\alpha)}(A,BC)\;.$
(5)

For three adjacent blocks, presented in figure 1, unlike the mutual
information, it has a desirable property of remaining bounded in the limit
$1\ll|B|\ll|A|,|C|$, that will be of central interest. Here and in the
following $|A|$ stands for the size of $A$. Namely, let us suppose that the
entanglement entropy of a (large) connected block scales as

$S_{\alpha}(A)=a_{\alpha}|A|+b_{\alpha}\log|A|+c_{\alpha}\;,$ (6)

for some constants $a_{\alpha},b_{\alpha},c_{\alpha}$. This is the case both
in CFT ($a_{\alpha}=0,b_{\alpha}\neq 0$) and in the stationary states studied
in this work ($a_{\alpha}\neq 0$). Then, for the configuration in figure 1 the
relation between the tripartite information of $A,B,C$ and the mutual
information of $A,C$ is (in the limit of large blocks)

$I_{2}^{(\alpha)}(A,C)=-b_{\alpha}\log(1-x)+I_{3}^{(\alpha)}(A,B,C)\;,$ (7)

where

$x=\frac{|A||C|}{(|A|+|B|)(|B|+|C|)}$ (8)

is the cross ratio. The mutual and the tripartite information can thus be
simply obtained one from another, and the tripartite information describes a
subleading bounded contribution to the mutual information in the limit $x\to
1^{-}$, that corresponds to the limit $1\ll|B|\ll|A|,|C|$. Note that, in
general, the tripartite information can have any sign [56, 57] (even for
$\alpha=1$).

The tripartite information has been studied in many settings. In the context
of topological order in two space dimensions it is usually called simply
“topological entanglement entropy” [58]. It has been studied in quantum field
theory, both in 1+1 dimensions [56, 59, 60, 61, 37, 62, 63, 64, 65, 66, 67,
68, 57, 69, 70, 71, 72], mainly CFT, and in higher dimensions [56, 73, 74,
75]. In holographic theories it has been shown that the tripartite information
is never positive [76]. The tripartite information has also been studied in
continuously monitored chains [77], on Hamming graphs [78] and, with the
partition at different time slices, as a diagnostic tool for quantum
scrambling [79, 80, 81, 82]. Finally, it has also been addressed after global
quenches in closed quantum systems [1, 2, 83, 84].

In this work we focus on the tripartite information of three large adjacent
subsystems embedded in an infinite chain (see figure 1). We refer the reader
to section 2 of paper II for a discussion of the behavior of the tripartite
information in this setting for different systems. Here we note that it is
expected to vanish when the spin correlation functions decay exponentially
with distance, as is the case in thermal states [85] and ground states of
gapped Hamiltonians [86]. An exception are ground states of conformal critical
systems, where the tripartite information is a model dependent function of the
cross ratio in (8) [37, 23]. For example, in the ground state of the XX chain
the Rényi-$2$ tripartite information reads [60, 37, 57]

$I_{3}^{(2)}(A,B,C)=-\log 2+\log\left(1+\sqrt{1-x}+\sqrt{x}\right)\;.$ (9)

Papers I and II identified different global quench protocols that result in a
stationary state in which the spin correlation functions decay algebraically
(in space) and allow for non-zero tripartite information: bipartitioning
protocols and quenches from ground states of critical systems. In the
stationary states of these protocols the entropy of a large single (connected)
block $A$ satisfies scaling (6) with both $a_{\alpha}$ and $b_{\alpha}$
nonzero. The tripartite information of three large adjacent blocks in these
stationary states exhibits a universal dependency on the cross ratio (8),
similarly to critical systems. For example, the quench from the domain wall
state
$\ket{\ldots\uparrow\uparrow\uparrow\downarrow\downarrow\downarrow\ldots}$
under the XY chain Hamiltonian results in the stationary state with the
Rényi-$2$ tripartite information

$I_{3}^{(2)}(A,B,C)=-\log 2+\log\left(1+\sqrt{1-x}\right)\;,$ (10)

which corresponds to (9) with the last term in the logarithm dropped.

Differently from CFT, the tripartite information in the aforementioned
stationary states is not symmetric under the interchange $x\leftrightarrow
1-x$. This property is particularly transparent in the limit $x=1^{-}$, which
corresponds to the limit in which the central subsystem is much smaller than
the others, i.e. the limit $1\ll|B|\ll|A|,|C|$. In this limit the tripartite
information has a nonzero value, termed “residual tripartite information” in
Paper I. The nonzero value is $-\log 2$ and it is common to all studied quench
protocol with non-zero tripartite information. Moreover, this value is
independent of the Rényi index $\alpha$, so it applies also to the von Neumann
tripartite information. It should be contrasted to the zero value found in
equilibrium at any temperature, irrespectively of criticality, or in other
non-equilibrium settings, such as after quenches from ground states of gapped
Hamiltonians. We note that there are also trivial examples with nonzero
residual tripartite information, such as the GHZ state
$\ket{\textrm{GHZ}}=(\ket{\ldots\uparrow\uparrow\uparrow\ldots}+\ket{\ldots\downarrow\downarrow\downarrow\ldots})/\sqrt{2}$,
which has tripartite information $+\log 2$ independently of the size of the
subsystems (as long as the subsystems do not comprise the whole system, in
which case the tripartite information would be zero since the state is pure),
but such examples differ from the systems we are interested in because they
violate clustering.

We note that algebraic decay of some spin correlation functions is not a
sufficient condition for a non-zero tripartite information, as illustrated by
the global quench from the ground state of the Ising chain with a critical
transverse field (see Papers I and II). The latter reaches a stationary state
in which some spin correlation functions decay with distance $r$ as $1/r^{4}$,
yet the tripartite information is zero. This phenomenology seems to be related
to the degree in the power law, since quenches with non-zero tripartite
information posses correlations that decay only as $1/r^{2}$.

###  I.2 Model and Quench Protocols

##### Model

We study the time evolution governed by the dual XY chain [87, 46], given by
the Hamiltonian

$\displaystyle\bm{H}=\sum_{\ell=-\infty}^{\infty}\bm{\sigma}^{x}_{\ell-1}(J_{x}\bm{I}-J_{y}\bm{\sigma}^{z}_{\ell})\bm{\sigma}^{x}_{\ell+1}\;,$
(11)

where we assume $0<|J_{y}|<|J_{x}|$. A particularly interesting property of
the model is the existence of semilocal charges [46, 35], introduced in the
following. We also mention that the special point $J_{x}=J_{y}$ has been
studied for its kinetic constraints [88, 89, 90].

The Hamiltonian is invariant under spin flip $\mathcal{P}_{\sigma}^{z}$, where

$\mathcal{P}_{\sigma}^{z}[\bm{O}]=\lim_{n\to\infty}\left(\prod_{\ell=-n}^{n}\bm{\sigma}_{\ell}^{z}\right)\bm{O}\left(\prod_{\ell^{\prime}=-n}^{n}\bm{\sigma}_{\ell^{\prime}}^{z}\right)$
(12)

for some localized operator $\bm{O}$. Here and in the following we say that a
local operator $\bm{Q}=\sum_{\ell\in\mathbb{Z}}\bm{q}_{\ell}$, where the
density $\bm{q}_{\ell}$ is localized (meaning finite support), is even/odd
under some transformation $\mathcal{P}$ if
$\mathcal{P}[\bm{q}_{\ell}]=\pm\bm{q}_{\ell}$ respectively. We will use the
bold notation exclusively for operators defined on the whole, infinite, chain.

A charge $\bm{Q}$ is an operator commuting with the Hamiltonian
($[\bm{Q},\bm{H}]=0$). While for local charges
$\bm{Q}=\sum_{\ell\in\mathbb{Z}}\bm{q}_{\ell}$ the density $\bm{q}_{\ell}$ is
localized around site $\ell$, semilocal charges of this model are
characterised by density $\bm{q}_{\ell}$ with support on all sites on one side
of $\ell$. Namely, defining the operators $\bm{\Pi}^{z}_{\sigma,+}(\ell)$ by

$\bm{\Pi}^{z}_{\sigma,+}(\ell)\bm{\sigma}_{j}^{x,y}=\begin{cases}\bm{\sigma}_{j}^{x,y}\bm{\Pi}^{z}_{\sigma,+}(\ell)&j<\ell\\\
-\bm{\sigma}_{j}^{x,y}\bm{\Pi}^{z}_{\sigma,+}(\ell)&j\geq\ell\end{cases},\qquad\left[\bm{\Pi}^{z}_{\sigma,+}(\ell),\bm{\sigma}_{j}^{z}\right]=0,\qquad\left(\bm{\Pi}^{z}_{\sigma,+}(\ell)\right)^{2}=\bm{1},$
(13)

which can be thought of as semi-infinite strings
$\bm{\sigma}_{\ell}^{z}\bm{\sigma}_{\ell+1}^{z}\bm{\sigma}_{\ell+2}^{z}\ldots$,
an example of a semilocal charge of the model in (11) is

$\bm{Q}^{(0,-)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\bm{\sigma}_{\ell-1}^{x}\left(\bm{I}-\bm{\sigma}_{\ell}^{z}\right)\bm{\sigma}_{\ell+1}^{y}\bm{\Pi}_{\sigma,+}^{z}(\ell+2)\;.$
(14)

A complete set of one-site shift invariant semilocal charges of the model is
given in section III.1.1. Their construction is directly related to the
Kramers-Wannier duality. Semilocal charges are particularly interesting
because they allow the sytem to keep the memory of localized perturbations in
the initial state. Namely, while localized perturbations can affect the
expectation value of localized charge densities $\bm{q}_{\ell}$ at most for
several sites $\ell$, the effect on semilocal charges can be drastic.

##### Quench protocols

In this work we study and compare two different global quench protocols:

  1. 1)

all-spin-up quench protocol: time evolution of a translationally invariant
product state, $\ket{\Psi(t)}=e^{-i\bm{H}t}\ket{\Uparrow}$;

  2. 2)

flipped-spin quench protocol: time evolution of the state in 1) with a flipped
spin at a single site,
$\ket{\Psi(t)}=e^{-i\bm{H}t}\bm{\sigma}_{0}^{x}\ket{\Uparrow}$.

Here and in the rest of the paper by $\Uparrow$ we denote an inifinite string
of $\uparrow$. We can thus also write
$\bm{\sigma}_{0}^{x}\ket{\Uparrow}=\ket{\Uparrow\downarrow_{0}\Uparrow}$. In
the flipped-spin protocol the system will for large times reach quasi-
stationary states at different rays $\zeta=d/t$, where $d$ is the distance
from the site with the flipped spin. When speaking about the stationary state
we will always refer to the NESS reached at infinite time around the initial
inhomogeneity, corresponding to the ray $\zeta=0$.

##### Kramers-Wannier transformation

We consider a transformation that differs from the standard Kramers-Wannier
duality map, responsible for the self-duality of the transverse-field Ising
model, just in an additional rotation,

$\bm{\tau}_{j}^{x}=\bm{\sigma}_{j-1}^{x}\bm{\sigma}_{j}^{x}\,,\qquad\bm{\tau}_{j}^{z}\bm{\tau}_{j+1}^{z}=\bm{\sigma}_{j}^{z}\,,$
(15)

that we will also refer to as Kramers-Wannier transformation. We will refer to
the representations of the theory in terms of the $\bm{\tau}$ operators as the
dual picture.

The transformation preserves the algebra of Pauli matrices, i.e. the
$\bm{\tau}$ operators satisfy the same algebra as the $\bm{\sigma}$ ones.
Relation (15) specifies the transformation only for localized operators even
under $\mathcal{P}_{\sigma}^{z}$, which are transformed to operators that are
localized also in the dual picture and that are, moreover, even under
$\mathcal{P}^{x}_{\tau}$, where

$\mathcal{P}^{x}_{\tau}[\bm{O}]=\lim_{n\to\infty}\left(\prod_{\ell=-n}^{n}\bm{\tau}_{\ell}^{x}\right)\bm{O}\left(\prod_{\ell^{\prime}=-n}^{n}\bm{\tau}_{\ell^{\prime}}^{x}\right)\;.$
(16)

In this work, the full mapping is needed only for the discussion of semilocal
charges. Details are given in section III.1.

Kramers-Wannier transformation maps the dual XY chain (11) to the XY chain

$\bm{H}=\sum_{\ell=-\infty}^{\infty}\left(J_{x}\bm{\tau}_{\ell}^{x}\bm{\tau}_{\ell+1}^{x}+J_{y}\bm{\tau}_{\ell}^{y}\bm{\tau}_{\ell+1}^{y}\right)\;,$
(17)

a very well known model mappable to free fermions [91], with substantially
developed techniques for quantum quenches [92, 93, 94, 95, 96]. The state of
all spin up maps into itself, while the state with a spin flip maps into a
domain wall state. Namely, we have the identification
$\ket{\Uparrow}^{(\sigma)}=\ket{\Uparrow}^{(\tau)}$ and
$\ket{\Uparrow\downarrow_{0}\Uparrow}^{(\sigma)}=\ket{\Uparrow\Downarrow}^{(\tau)}$
with the domain wall between sites $0$ and $1$. Note that these are not unique
identifications and one could, for example, identify
$\ket{\Uparrow}^{(\sigma)}$ with $\ket{\Downarrow}^{(\tau)}$ instead of
$\ket{\Uparrow}^{(\tau)}$ or with any linear combination of the two. In a
finite system [46] there is a unique choice, related to the boundary
conditions at the edge of the system. Since we let the boundaries to infinity,
for our purposes different choices are equivalent. In any case, a single spin
flip changes the state macroscopically in the dual picture. We note that these
aspects of the mapping are very similar to the hidden symmetry-breaking
picture of the symmetry-protected topological phases [97, 98, 99].

##### Reduced density matrices and tripartite information

We stress it is a highly non-trivial question what are the effects of the
Kramers-Wannier transformation on the tripartite information. For example, the
Jordan-Wigner transformation is also a duality transformation and there are
important differences between the entanglement entropy of disjoint blocks of
spins and fermions [68, 100], and therefore also in the tripartite
information. In fact, this difference is crucial for the phenomenology
discovered in Papers I and II. However, it turns out that the Kramers-Wannier
transformation in the studied stationary states does not affect the tripartite
information of three large adjacent blocks.

For a system in a state $\ket{\Psi}$, the reduced density matrix for subsystem
$X$ is given by

$\rho_{X}=\frac{1}{2^{|X|}}\sum_{\gamma_{\ell}\in\\{0,x,y,z\\},\ell\in
X}\bra{\Psi}\prod_{\ell\in
X}\bm{\sigma}^{\gamma_{\ell}}_{\ell}\ket{\Psi}\bigotimes_{\ell\in
X}\sigma^{\gamma_{\ell}}\;.$ (18)

Here the sum is over an orthogonal basis of operators on $X$, given by all
possible products of Pauli matrices, and we use the standard convention
$\sigma^{0}\equiv\mathbb{I}$. While the sites $\ell$ of some physical
subsystem $X$ are associated to the indices of the operators
$\bm{\sigma}_{\ell}^{\gamma}$, in the dual picture we associate the notion of
subsystem to the indices of the $\bm{\tau}^{\gamma}_{\ell}$ operators.
Accordingly, in the dual picture it is natural to consider the density matrix

$\rho_{X}^{\tau}=\frac{1}{2^{|X|}}\sum_{\gamma_{\ell}\in\\{0,x,y,z\\},\ell\in
X}\bra{\Psi}\prod_{\ell\in
X}\bm{\tau}^{\gamma_{\ell}}_{\ell}\ket{\Psi}\bigotimes_{\ell\in
X}\sigma^{\gamma_{\ell}}\;.$ (19)

For the studied quench protocols the density matrix in (19) can be assessed
using the techniques developed for the model in (17). However, it is a non-
trivial question how are density matrices (18) and (19) related, or how are
the corresponding entanglement entropies related. These questions are some of
the main problems we tackle with in this work. As already mentioned, for a
single block of spins $X$ this problem has been addressed in [35]. The result
is that density matrices (18) and (19) give different entanglement entropies
in general, but the difference is bounded by a constant independent of the
subsystem size $|X|$. However, from this result we cannot conclude about the
tripartite information (4), which has also a contribution from the entropy of
disjoint blocks and which is, moreover, bounded itself. Thus, technically, the
main problem of this work is to study the density matrix (19) when $X$
consists of disjoint blocks $A,C$, as presented in section III. In the end we
are able to conclude that in the studied stationary states the differences in
the entropies corresponding to density matrices (18) and (19) cancel in the
tripartite information of three large adjacent blocks.

##  II Results

Figure 2: Rényi-$2$ mutual information between subsystems $A$ and $C$ (see
figure 1) for $|A|=|C|=16$ and varying distance $|B|$, in the stationary
states following quench protocols 1 and 2. In the stationary state of protocol
1 mutual information decays exponentially with distance ($y$-axis is in the
log scale), while in the stationary state of protocol 2 it decays
algebraically.

###  II.1 Correlation functions

Computing the spin-correlation functions in the stationary states of protocols
1 and 2 is a straightforward application of the Kramers-Wannier duality and
the Jordan-Wigner transformation for the dual model, as commented in section
III.2.

The connected correlation functions of the $z$-components of spins in the
stationary state following quench protocol 1 are exactly zero,

$\textrm{protocol \ref{quench protocol
1}:}\qquad\braket{\bm{\sigma}_{0}^{z}\bm{\sigma}_{r}^{z}}-\braket{\bm{\sigma}_{0}^{z}}\braket{\bm{\sigma}_{r}^{z}}=0,\qquad
r\geq 1,$ (20)

while in the stationary state of protocol 2 they exhibit algebraic decay with
distance,

$\textrm{protocol \ref{quench protocol
2}:}\qquad\braket{\bm{\sigma}_{0}^{z}\bm{\sigma}_{r}^{z}}-\braket{\bm{\sigma}_{0}^{z}}\braket{\bm{\sigma}_{r}^{z}}\simeq\frac{16}{\pi^{4}}\frac{1}{r^{4}},$
(21)

where $\simeq$ means asymptotically equal. The connected spin-correlation
functions of other spin components are exactly zero (for large enough
distance) both in the stationary state following quench protocol 1 and the one
following protocol 2,

$\textrm{protocols \ref{quench protocol 1} and \ref{quench protocol
2}:}\qquad\braket{\bm{\sigma}_{0}^{\gamma}\bm{\sigma}_{r}^{\gamma}}-\braket{\bm{\sigma}_{0}^{\gamma}}\braket{\bm{\sigma}_{r}^{\gamma}}=0,\qquad\gamma=x,y,\
r\geq 5.$ (22)

These results deal only with the correlation functions of operators with
support on one site, which already establish that there are algebraically
decaying correlations in the stationary state of quench protocol 2. For the
stationary state of quench protocol 1 we provide evidence that there are no
connected correlation functions decaying algebraically by computing the
Rényi-$2$ mutual information. The results are presented in figure 2. In the
stationary state of protocol 1 the mutual information $I_{2}^{(2)}(A,C)$
decays exponentially with the distance between $A$ and $C$, while it decays
only algebraically in the stationary state of protocol 2. Assuming that the
von Neumann mutual information $I_{2}^{(1)}(A,C)$ follows an analogous
behavior, we can conclude that all connected correlation functions decay
exponentially in protocol 1 [50].

We note that in the stationary state of quench protocol 2 it is possible to
construct operators whose connected correlation functions decay only as
$1/r^{2}$, as opposed to the $1/r^{4}$ decay of the spin correlation functions
in eq. (21). This is the case for the operator
$\bm{\xi}_{\ell}\equiv\bm{\sigma}_{\ell-1}^{x}\bm{\sigma}_{\ell}^{z}\bm{\sigma}_{\ell+1}^{x}$,
with the correlation functions

$\begin{split}\textrm{protocol \ref{quench protocol
2}:}\qquad\braket{\bm{\xi}_{0}\bm{\xi}_{2r-1}}-\braket{\bm{\xi}_{0}}\braket{\bm{\xi}_{2r-1}}\simeq-\frac{4}{\pi^{2}r^{2}},\qquad\braket{\bm{\xi}_{0}\bm{\xi}_{2r}}-\braket{\bm{\xi}_{0}}\braket{\bm{\xi}_{2r}}=0,\;r\geq
1.\end{split}$ (23)

Figure 3: The analytical prediction for the Rényi-$\alpha$ tripartite
information of three large adjacent blocks (see figure 1) in the stationary
state of the flipped-spin quench protocol, given by (29), as a function of the
cross ratio (8), for different values of $\alpha$. The curves do not differ
much. The Rényi-$\alpha$ tripartite information for any $\alpha$ is equal to
zero for $x=0^{+}$ and equal to $-\log 2$ for $x=1^{-}$. By replica trick the
result applies also to the von Neumann tripartite information. A nonzero value
in the limit $x=1^{-}$ was termed “residual tripartite information” in ref.
[1].

###  II.2 String order

As discussed in ref. [35] in the context of the all-spin-up quench protocol,
string order can survive a global quench. Accordingly, the order in the
stationary state was termed “non-equilibrium symmetry-protected topological
order”. To be more precise, a string of adjacent $\bm{\sigma}^{z}_{\ell}$
operators has a non-zero expectation value in the stationary state even in the
limit of infinite length of the string. This is related to the fact that the
expectation value of a string of $\bm{\sigma}^{z}$ in the dual picture becomes
a two-point correlation function (not the connected one) of localized
operators. We note that in quench protocol 2 this string order parameter
vanishes. Explicitly, we have

$\displaystyle\textrm{protocol \ref{quench protocol
1}:}\qquad\lim_{r\to\infty}\braket{\prod_{\ell=-r}^{r}\bm{\sigma}^{z}_{\ell}}=\frac{1}{4}\left(1+\frac{J_{y}}{J_{x}}\right)^{2},$
(24) $\displaystyle\textrm{protocol \ref{quench protocol
2}:}\qquad\lim_{r\to\infty}\braket{\prod_{\ell=-r}^{r}\bm{\sigma}^{z}_{\ell}}=0.$
(25)

Here we point out that there is a string order also in the stationary state of
protocol 2, but given by a different string order parameter, which, on the
other hand, vanishes in protocol 1. Namely, we have

$\displaystyle\textrm{protocol \ref{quench protocol
1}:}\qquad\lim_{r\to\infty}\braket{\bm{\sigma}_{-r-3}^{x}\bm{\sigma}_{-r-2}^{z}\bm{\sigma}_{-r-1}^{y}\left(\prod_{\ell=-r}^{r}\bm{\sigma}^{z}_{\ell}\right)\bm{\sigma}^{y}_{r+1}\bm{\sigma}_{r+2}^{z}\bm{\sigma}_{r+3}^{x}}=0\;,$
(26) $\displaystyle\textrm{protocol \ref{quench protocol
2}:}\qquad\lim_{r\to\infty}\braket{\bm{\sigma}_{-r-3}^{x}\bm{\sigma}_{-r-2}^{z}\bm{\sigma}_{-r-1}^{y}\left(\prod_{\ell=-r}^{r}\bm{\sigma}^{z}_{\ell}\right)\bm{\sigma}^{y}_{r+1}\bm{\sigma}_{r+2}^{z}\bm{\sigma}_{r+3}^{x}}=-\frac{1}{\pi^{2}}\left(1+\frac{J_{y}}{J_{x}}\right)^{2}\;.$
(27)

The two string operators differ only at the edges, while in the bulk they are
given by the on-site action of the symmetry of the model (as in the symmetry-
protected topological order at zero temperature [98, 101]). The computation of
the string order parameters is commented in section III.2 and appendix A.

###  II.3 Tripartite information

Figure 4: Rényi-$2$ tripartite information of three adjacent blocks $A,B,C$
(see figure 1) of equal length, that we vary. In the stationary state of the
all-spin-up quench protocol the tripartite information vanishes exponentially
with the size of the subsystems, while in the stationary state of the flipped-
spin quench protocol it tends to a non-zero value, for which we have
analytical prediction, given by eq. (10) (solid line).

Based on the representation of the reduced density matrix in the dual picture,
derived in section III.1.3, we argue in section III.1.4 that the tripartite
information in the stationary states of the studied global quenches is not
affected by the Kramers-Wannier transformation, in the limit of large
subsystems. In this way we reduce the problem of computing the tripartite
information in the stationary state of quench protocols 1 and 2 to computing
the tripartite information in the stationary state of the quench protocol in
which the time evolution is with the XY model and the initial state is,
respectively, 1) the state of all spin up $\ket{\Uparrow}$, 2) the domain wall
state $\ket{\Uparrow\Downarrow}$. These problems have already been studied in
details in Papers I and II so we copy the results derived there, and compare
them with the exact numerical results for the second Rényi entropy, obtained
using the methods developed in this work.

In this way we find the following tripartite information in the stationary
states of the studied quench protocols, in the limit of large subsystems:

$\displaystyle\textrm{protocol \ref{quench protocol 1}:}\qquad
I_{3}^{(\alpha)}(A,B,C)=0,$ (28) $\displaystyle\textrm{protocol \ref{quench
protocol 2}:}\qquad
I_{3}^{(\alpha)}(A,B,C)=\frac{1}{\alpha-1}\log\Bigl{[}\sum_{\delta_{j}\in\\{0,\frac{1}{2}\\}\atop
j=1,\dots,\alpha-1}\Bigl{(}\tfrac{\Theta(\vec{\delta}|\hat{\tau}_{x})}{\Theta(\vec{0}|\hat{\tau}_{x})}\Bigr{)}^{2}\Bigr{]}-\log
2\;.$ (29)

Here $\hat{\tau}_{x}$ is the $(\alpha-1)\times(\alpha-1)$ period matrix of the
Riemann surface $\mathcal{R}_{\alpha}$ with elements [37]

$[\hat{\tau}_{x}]_{\ell
n}=\frac{2i}{\alpha}\sum_{k=1}^{\alpha-1}\sin(\tfrac{\pi
k}{\alpha})\cos(\tfrac{2\pi
k(\ell-n)}{\alpha})\tfrac{P_{(k/\alpha)-1}(2x-1)}{P_{(k/\alpha)-1}(1-2x)}\,,$
(30)

where $P_{\mu}(z)$ denotes the Legendre functions and
$\Theta(\vec{z},M)=\sum_{\vec{m}\in\mathbb{Z}^{\alpha-1}}e^{i\pi\vec{m}^{t}M\vec{m}+2\pi
i\vec{m}\cdot\vec{\delta}}$ is the Siegel theta function. The formula (29) is
plotted in figure 3 for several values of $\alpha$. For $\alpha=2$ it reduces
simply to eq. (10). We note that eq. (29) was obtained in Paper I by
establishing a correspondence between some contributions to the entanglement
entropy of disjoint blocks in the expansion of ref. [68] and already known CFT
resuts [102], while the results for $\alpha=2$ and some simpler approximate
formulas for higher $\alpha$ were obtained by a direct computation in Paper
II. Note also that the result (29) is independent of the parameters of the
model (11).

The tripartite information of three large adjacent blocks is zero in the
stationary state of the all-spin-up quench protocol, while in the stationary
state of the flipped-spin quench protocol it is nonzero. This conclusion is
confirmed by numerical results for the second Rényi entropy, presented in
figure 4. There it can be seen that in the stationary state of the all-spin-up
quench protocol the tripartite information goes to zero exponentially fast
with the subsystem size, while in the stationary state of protocol 2 the
tripartite information reaches the asymptotic result (10).

The tripartite information in the stationary state of the flipped-spin
protocol exhibits a universal dependency on the cross ratio (8), as confirmed
by numerical results in figure 5. The tripartite information is a function
only of the cross ratio, notwithstanding that the initial state and the time
evolution Hamiltonian are not related to CFT. Note that the tripartite
information is zero in the limit $x=0^{+}$, as in CFT. Differently from CFT,
we find nonzero residual tripartite information

$\textrm{protocol \ref{quench protocol 2}:}\qquad
I_{3}^{(\alpha)}(A,B,C)\stackrel{{\scriptstyle x\to 1^{-}}}{{=}}-\log 2\;,$
(31)

which corresponds to the limit of small separation between $A$ and $C$ with
respect to their size, i.e. to the limit $1\ll|B|\ll|A|,|C|$. The result is
the same for all $\alpha$ so it holds also for the von Neumann tripartite
information (by replica trick).

Figure 5: The Rényi-$2$ tripartite information in the stationary state of the
flipped-spin quench protocol for different configurations with fixed value of
the cross ratio $x$, defined by (8). From top to bottom the configurations are
$(|A|,|B|,|C|)\propto(5,5,1),(1,1,1),(2,1,3),(4,1,7),(10,1,10)$. Increasing
the size of the subsystems the tripartite information reaches the analytical
prediction (10), given by the solid line.

##  III Methods

###  III.1 Kramers-Wannier transformation

The Kramers-Wanniers transformation is fully specified, as discussed in [35],
by the relations

$\displaystyle\bm{\Pi}^{x}_{\tau,-}(\ell)=$
$\displaystyle\bm{\sigma}_{\ell}^{x}\,,$ (32)
$\displaystyle\bm{\tau}_{\ell}^{y}=$
$\displaystyle\bm{\sigma}_{\ell-1}^{x}\bm{\sigma}_{\ell}^{y}\bm{\Pi}_{\sigma,+}^{z}(\ell+1)\,,$
$\displaystyle\bm{\tau}_{\ell}^{z}=$
$\displaystyle\bm{\Pi}_{\sigma,+}^{z}(\ell)\,,$

where the operator $\bm{\Pi}^{z}_{\sigma,+}(\ell)$ is defined in eq. (13) and
$\bm{\Pi}^{x}_{\tau,-}(\ell)$ is defined by

$\bm{\Pi}^{x}_{\tau,-}(\ell)\bm{\tau}_{j}^{y,z}=\begin{cases}-\bm{\tau}_{j}^{y,z}\bm{\Pi}^{x}_{\tau,-}(\ell)&j\leq\ell\\\
\bm{\tau}_{j}^{y,z}\bm{\Pi}^{x}_{\tau,-}(\ell)&j>\ell\end{cases},\qquad\left[\bm{\Pi}^{x}_{\tau,-}(\ell),\bm{\tau}_{j}^{x}\right]=0,\qquad\left(\bm{\Pi}^{x}_{\tau,-}(\ell)\right)^{2}=\bm{1}$
(33)

The operators $\bm{\Pi}^{z}_{\sigma,+}(\ell)$ and
$\bm{\Pi}^{x}_{\tau,-}(\ell)$ are semilocal(ized) in the $\sigma$ and $\tau$
basis respectively. They can be thought of as semi-infinite strings
$\bm{\sigma}_{\ell}^{z}\bm{\sigma}_{\ell+1}^{z}\bm{\sigma}_{\ell+2}^{z}\ldots$
and $\ldots\bm{\tau}_{\ell-2}^{x}\bm{\tau}_{\ell-1}^{x}\bm{\tau}_{\ell}^{x}$
respectively. In particular, operators $\bm{O}$ that are localized in the
$\tau$-representation and odd under $\mathcal{P}_{\tau}^{x}$
($\mathcal{P}_{\tau}^{x}[\bm{O}]=-\bm{O}$) are in the $\sigma$ representation
even under $\mathcal{P}_{\sigma}^{z}$
($\mathcal{P}_{\sigma}^{z}[\bm{O}]=\bm{O}$) and semilocal.

####  III.1.1 Semilocal charges

In this section we give a complete set of one-site shift invariant semilocal
charges of the dual XY chain (11). A complete set of one-site shift invariant
charges of the quantum XY chain (17) is given by (see e.g. [95, 103])

$\displaystyle\bm{Q}^{(0,+)}=\bm{H},$ (34)
$\displaystyle\bm{Q}^{(1,+)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\left[J_{x}\bm{\tau}_{\ell}^{x}\bm{\tau}^{x}_{\ell+2}+J_{y}\bm{\tau}_{\ell}^{y}\bm{\tau}^{y}_{\ell+2}-(J_{x}+J_{y})\bm{I}\right]\bm{\tau}_{\ell+1}^{z},$
(35)
$\displaystyle\bm{Q}^{(n,+)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\left[\left(J_{x}\bm{\tau}_{\ell}^{x}\bm{\tau}_{\ell+n+1}^{x}+J_{y}\bm{\tau}_{\ell}^{y}\bm{\tau}_{\ell+n+1}^{y}\right)\prod_{j=1}^{n}\bm{\tau}_{\ell+j}^{z}+\left(J_{x}\bm{\tau}_{\ell}^{x}\bm{\tau}_{\ell+n-1}^{x}+J_{y}\bm{\tau}_{\ell}^{y}\bm{\tau}_{\ell+n-1}^{y}\right)\prod_{j=1}^{n-2}\bm{\tau}_{\ell+j}^{z}\right],$
(36) $\displaystyle\qquad\textrm{for}\quad n=2,3,4,\ldots\quad,$
$\displaystyle\bm{Q}^{(0,-)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\left(\bm{\tau}_{\ell}^{x}\bm{\tau}_{\ell+1}^{y}-\bm{\tau}_{\ell}^{y}\bm{\tau}_{\ell+1}^{x}\right),$
(37)
$\displaystyle\bm{Q}^{(n,-)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\left[\left(\bm{\tau}_{\ell}^{x}\bm{\tau}_{\ell+n+1}^{y}-\bm{\tau}_{\ell}^{y}\bm{\tau}_{\ell+n+1}^{x}\right)\prod_{j=1}^{n}\bm{\tau}_{\ell+j}^{z}\right],\qquad
n=1,2,3,\ldots$ (38)

These charges all mutually commute ($[\bm{Q}^{(n,+)},\bm{Q}^{(m,\pm)}]=0$ for
$n,m\in\mathbb{N}_{0}$). We note that quantum XY chain possesses also charges
that are two-site, and not one-site, shift invariant, which are non-Abelian in
general [103]. These charges are not relevant for the quench protocols
considered in this work and will not be discussed. We also note that in
interacting systems in general one has to go beyond local charges to consider
quasilocal charges [104], but this is not the case in the studied non-
interacting system.

Applying the duality transformation to the charges, we obtain the charges of
the dual XY chain (11). While the charges even under $\mathcal{P}_{\tau}^{x}$
are local also in the $\sigma$-representation, a simple example being provided
by the Hamiltonian, the charges odd under $\mathcal{P}_{\tau}^{x}$ are instead
semilocal, as discussed in ref. [35] (which labels charges in a different
way). Specifically, the charges odd under $\mathcal{P}_{\tau}^{x}$ are charges
$\bm{Q}^{(n,+)}$ for odd $n$ and charges $\bm{Q}^{(n,-)}$ for even $n$. These
charges form a complete set of one-site shift invariant semilocal charges of
the model in (11). The $\sigma$-representation of the charge $\bm{Q}^{(0,-)}$
has been already given in eq. (14). The remaining one-site shift invariant
semilocal charges of the model in (11) read

$\displaystyle\bm{Q}^{(1,+)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\left[\left(-J_{x}\bm{\sigma}_{\ell}^{x}\bm{\sigma}_{\ell+1}^{y}+J_{y}\bm{\sigma}_{\ell}^{y}\bm{\sigma}_{\ell+1}^{x}\right)\bm{\sigma}_{\ell-1}^{x}\bm{\sigma}_{\ell+2}^{y}-(J_{x}+J_{y})\bm{\sigma}_{\ell+1}^{z}\bm{\sigma}_{\ell+2}^{z}\right]\bm{\Pi}^{z}_{\sigma,+}(\ell+3)$
(39)
$\displaystyle\bm{Q}^{(n,+)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\Bigg{[}\left(-J_{x}\bm{\sigma}_{\ell}^{x}\bm{\sigma}_{\ell+n}^{y}\prod_{j=1}^{(n-1)/2}\bm{\sigma}_{\ell+2j-1}^{z}+J_{y}\bm{\sigma}_{\ell}^{y}\bm{\sigma}_{\ell+n}^{x}\prod_{j=1}^{(n-1)/2}\bm{\sigma}_{\ell+2j}^{z}\right)\bm{\sigma}_{\ell-1}^{x}\bm{\sigma}_{\ell+n+1}^{y}+$
$\displaystyle\left(-J_{x}\bm{\sigma}_{\ell}^{x}\bm{\sigma}_{\ell+n-2}^{y}\prod_{j=1}^{(n-3)/2}\bm{\sigma}_{\ell+2j-1}^{z}+J_{y}\bm{\sigma}_{\ell}^{y}\bm{\sigma}_{\ell+n-2}^{x}\prod_{j=1}^{(n-3)/2}\bm{\sigma}_{\ell+2j}^{z}\right)\bm{\sigma}_{\ell-1}^{x}\bm{\sigma}_{\ell+n-1}^{y}\bm{\sigma}_{\ell+n}^{z}\bm{\sigma}_{\ell+n+1}^{z}\Bigg{]}\bm{\Pi}^{z}_{\sigma,+}(\ell+n+2)$
(40) $\displaystyle\qquad\mathrm{for}\quad n=3,5,7,\ldots\quad,$
$\displaystyle\bm{Q}^{(n,-)}=\frac{1}{2}\sum_{\ell=-\infty}^{\infty}\left(\bm{\sigma}_{\ell}^{x}\bm{\sigma}_{\ell+n}^{x}\prod_{j=1}^{n/2}\bm{\sigma}_{\ell+2j-1}^{z}+\bm{\sigma}_{\ell}^{y}\bm{\sigma}_{\ell+n}^{y}\prod_{j=1}^{n/2-1}\bm{\sigma}_{\ell+2j}^{z}\right)\bm{\sigma}_{\ell-1}^{x}\bm{\sigma}_{\ell+n+1}^{y}\bm{\Pi}_{\sigma,+}^{z}(\ell+n+2)$
(41) $\displaystyle\qquad\mathrm{for}\quad n=2,4,6,\ldots\quad,$

where we use the convention that products of the form $\prod_{j=1}^{0}$ are
equal to identity.

####  III.1.2 Semilocal Generalized Gibbs ensemble

As pointed out in [35], the stationary state of quench protocol 1 is a
semilocal generalized Gibbs ensemble, meaning that it includes semilocal
charges. Here we point out that the non-equilibrium stationary state emerging
from the flipped-spin quench protocol is also described by a semilocal GGE.

In both cases the GGE can be obtained directly by working in the dual picture
and following the approach of ref. [95], as explained in appendix B. We find
that the GGE

$\bm{\rho}=\frac{1}{\mathcal{Z}}\exp\left[-\sum_{n=0}^{\infty}\left(\lambda^{(n,+)}\bm{Q}^{(n,+)}+\lambda^{(n,-)}\bm{Q}^{(n,-)}\right)\right]\;,$
(42)

where $\mathcal{Z}$ is the normalization, is specified by Lagrange multipliers

protocol 1:
$\displaystyle\lambda^{(n,-)}=0,\qquad\lambda^{(n,+)}=\begin{cases}0,&\textrm{$n$
even}\\\
\frac{8}{\pi}\int_{0}^{\pi/2}\mathrm{artanh}\left[\frac{2(J_{x}+J_{y})\cos
k}{\varepsilon(k)}\right]\frac{1}{\varepsilon(k)}\cos(nk)\ dk,&\textrm{$n$
odd}\end{cases},$ (43) protocol 2:
$\displaystyle\lambda^{(n,-)}=\begin{cases}\mathrm{sgn}(J_{x}J_{y})\frac{8}{\pi}\int_{0}^{\pi/2}\mathrm{artanh}\left[\frac{2(J_{x}+J_{y})\cos
k}{\varepsilon(k)}\right]\sin[(n+1)k]\ dk,&\textrm{$n$ even}\\\ 0,&\textrm{n
odd}\end{cases},\quad\lambda^{(n,+)}=0,\quad$ (44)

for $n\in\mathbb{N}_{0}$.

The Lagrange multipliers are completely different in the two quench protocols,
but in both cases only those associated to semilocal charges are nonzero. In
the dual picture these are the charges odd under $\mathcal{P}_{\tau}^{x}$,
which will be important, in section III.1.4, for arguing that the tripartite
information is not influenced by the Kramers-Wannier transformation. We note
that the nonzero Lagrange multipliers decay rather slowly with $n$, i.e. only
algebraically, which is related to the logarithmic singularity of the function
$\mathrm{artanh}[2(J_{x}+J_{y})/\varepsilon(k)]$, appearing under the
integral, at $k=0$.

(a)$A$

(b)$\tilde{A}$

(c)$\tilde{A}$$P=\tau_{\mathrm{min}\
\\!\\!\tilde{A}}^{z}$$P=\tau_{\mathrm{max}\
\\!\\!\tilde{A}}^{z}$$P=\prod\limits_{\ell\in\tilde{A}}\tau_{\ell}^{x}$

(d)$A$$C$$B$

(e)$\tilde{A}$$\tilde{C}$$\tilde{B}$

(f)$\tilde{A}$$\tilde{C}$$P=\tau_{\mathrm{min}\
\\!\\!\tilde{A}}^{z}$$P=\tau_{\mathrm{max}\
\\!\\!\tilde{A}}^{z}$$P=\prod\limits_{\ell\in\tilde{A}}\tau_{\ell}^{x}$$P=\tau_{\mathrm{min}\
\\!\\!\tilde{C}}^{z}$$P=\tau_{\mathrm{max}\
\\!\\!\tilde{C}}^{z}$$P=\prod\limits_{\ell\in\tilde{C}}\tau_{\ell}^{x}$

Figure 6: (a) A connected block $A$. (b) Block $\tilde{A}$ we work with in the
dual picture, related to the configuration in a). (c) The transformations
$P_{j}$ appearing in eq. (49). (d) Disjoint blocks $A$ and $C$. (e) Blocks
$\tilde{A},\tilde{B},\tilde{C}$ we work with in the dual picture, related to
the configuration in d). (f) The transformations $P_{j}$ appearing in eq.
(54).

####  III.1.3 The reduced density matrix in the dual picture

In this section we find a convenient representation of the reduced density
matrix (18) in view of comparing it with (19). For this task we introduce the
Kramers-Wannier transformation for a finite system. Let us consider a system
of size $L$, described by Pauli matrices
$\sigma^{\gamma}_{\ell}\equiv\mathbb{I}^{\otimes(\ell-1)}\otimes\sigma^{\gamma}\otimes\mathbb{I}^{\otimes(L-\ell)}$
for $\gamma=0,x,y,z$ and $\ell=1,2,\ldots,L$, where $\sigma^{0}=\mathbb{I}$ is
a $2\times 2$ unit matrix. The transformation is given by

$\tau_{\ell}^{x}=\sigma_{\ell-1}^{x}\sigma_{\ell}^{x}\,,\qquad\tau_{\ell^{\prime}}^{z}\tau_{\ell^{\prime}+1}^{z}=\sigma_{\ell^{\prime}}^{z}\,,$
(45)

for $\ell=2,3,\ldots,L$ and $\ell^{\prime}=1,2,\ldots,L-1$. It is in
correspondence with the one in (15) for the infinite system. Eq. (45) does not
specify the transformation of operators $O$ that are odd under
$\Pi^{z}_{\sigma}\equiv\prod_{\ell=1}^{L}\sigma_{\ell}^{z}$
($[O,\Pi^{z}_{\sigma}]=-O$) nor the transformation at the boundaries. The full
mapping can be found in [46]. For our purposes it is just important that the
mapping preserves the algebra of Pauli matrices and that for operators even
under $\Pi^{z}_{\sigma}$ that are not at the boundaries the mapping is
analogous to the one for the infinite system. Note that to describe a
connected block $A$ we need to consider the mapping (45) for $L\geq|A|+1$, to
avoid the peculiarities of the mapping at the boundaries .

For a system in state $\ket{\Psi}$, the reduced density matrix for subsystem
$X$ is given by (18). In order to define a duality transformation in
correspondence with the one in the infinite system we consider an enlarged
space $\tilde{X}=X\cup X^{\prime}$, for some suitable choice of $X^{\prime}$
that will be discussed later. Note that we can extend $\rho_{X}$ by identity
to the enlarged space, defining

$\bar{\rho}_{X}\equiv\rho_{X}\otimes\left(\bigotimes\limits_{\ell\in
X^{\prime}}\frac{\mathbb{I}}{2}\right)\;.$ (46)

To obtain the entanglement entropy from the enlarged density matrix (46) we
just need to remove a constant related to the size of the enlargement
$X^{\prime}$,

$S_{\alpha}(X)=\frac{1}{1-\alpha}\log\mathrm{tr}(\bar{\rho}_{X}^{\alpha})-|X^{\prime}|\log
2\;.$ (47)

Having the duality transformation defined in the enlarged space, it is
convenient to consider the density matrices

$\tilde{\rho}_{Y}=\frac{1}{2^{|\tilde{X}|}}\sum_{\gamma_{\ell}\in\\{0,x,y,z\\},\ell\in
Y}\bra{\Psi}\prod_{\ell\in
Y}\bm{\tau}^{\gamma_{\ell}}_{\ell}\ket{\Psi}\prod_{\ell\in
Y}\tau_{\ell}^{\gamma_{\ell}}\;.$ (48)

for subsystems $Y\subseteq\tilde{X}$. Clearly, the density matrix (48) is in
correspondence with the density matrix $\rho^{\tau}_{Y}$, defined by (19), but
they are different mathematical objects. Instead of comparing the density
matrices (18) and (19) directly we will find it more convenient to compare the
density matrices (46) and (48), for a suitable choice of $Y$, concluding about
the former afterwards.

##### Single block

When the subsystem $X$ is a connected block $A$ (figure 6a) we consider the
enlarged subsystem $\tilde{A}=A\cup\\{\mathrm{max}(A)+1\\}$ (figure 6b). There
is a simple representation for
$\bar{\rho}_{A}=\rho_{A}\otimes\frac{\mathbb{I}}{2}$ in the dual picture,

$\bar{\rho}_{A}=\frac{1}{2^{3}}\sum_{j_{1},j_{2},j_{3}=0}^{1}\left(\prod_{m=1}^{3}{P_{m}^{j_{m}}}\right)\tilde{\rho}_{\tilde{A}}\left(\prod_{m=1}^{3}{P_{m}^{j_{m}}}\right)^{\dagger},$
(49)

with

$P_{1}=\prod_{\ell\in\tilde{A}}\tau_{\ell}^{x},\qquad
P_{2}=\tau_{\mathrm{min}(\tilde{A})}^{z},\qquad
P_{3}=\tau_{\mathrm{max}(\tilde{A})}^{z}.$ (50)

This result has already been obtained in [35] and we revisit it here for
comparison with the case of disjoint blocks. The density matrix
$\bar{\rho}_{A}$ is obtained from $\tilde{\rho}_{\tilde{A}}$ by successive
transformations of the form $\rho\to(\rho+P\rho P)/2$, for three different
Hermitian involutions $P$. Their effect is to remove those terms from
$\tilde{\rho}_{\tilde{A}}$ which act non-trivially on the site
$\mathrm{max}(A)+1$ (in the $\sigma$-representation), leaving only terms with
support on $A$. Note that transformation $P_{1}$ has support on the whole
system $\tilde{A}$ and acts by flipping all $\tau_{\ell}^{z}$, while the
remaining transformations affect only the boundaries. All three
transformations are represented graphically in figure 6c.

The representation (49) follows from the following properties, that are direct
consequences of the Kramers-Wannier transformation (15):

  1. 1)

Any operator $\bm{O}_{A}$ localized (in the $\sigma$ representation) in $A$
and even under $\mathcal{P}_{\sigma}^{z}$ is mapped by the Kramers-Wannier
transformation to an operator that in the $\tau$ representation is localized
in $\tilde{A}$ and is even under $\mathcal{P}_{\tau}^{x}$
($\mathcal{P}_{\tau}^{x}[\bm{O}_{A}]=\bm{O}_{A}$). More technically, if
$\bm{O}_{A}=\prod_{\ell\in A}\bm{\sigma}_{\ell}^{\gamma_{\ell}}$, where
$\gamma_{\ell}\in\\{0,x,y,z\\}$, is even under $\mathcal{P}_{\sigma}^{z}$ then
we have
$\bm{O}_{A}=s\prod_{\ell^{\prime}\in\tilde{A}}\bm{\tau}_{\ell^{\prime}}^{\gamma^{\prime}_{\ell^{\prime}}}$,
for some $\gamma_{\ell^{\prime}}^{\prime}\in\\{0,x,y,z\\}$, with $s\in\\{\pm
1\\}$ allowing for a minus sign.

  2. 2)

Any operator
$\bm{O}_{\tilde{A}}=\prod_{\ell\in\tilde{A}}\bm{\tau}_{\ell}^{\gamma_{\ell}}$,
where $\gamma_{\ell}\in\\{0,x,y,z\\}$, that is even under
$\mathcal{P}_{\tau}^{x}$ and which in the $\sigma$ representation has support
outside $A$, anticommutes with $\bm{\sigma}_{\min(A)-1}^{z}$
($\\{\bm{O}_{\tilde{A}},\bm{\sigma}_{\min(A)-1}^{z}\\}=0$) or
$\bm{\sigma}_{\max(A)+1}^{z}$
($\\{\bm{O}_{\tilde{A}},\bm{\sigma}_{\max(A)+1}^{z}\\}=0$), or both.

A simple example of property 2 are operators
$\bm{\tau}_{\min(\tilde{A})}^{x}=\bm{\sigma}_{\min(A)-1}^{x}\bm{\sigma}_{\min(A)}^{x}$
and
$\bm{\tau}^{x}_{\max(\tilde{A})}=\bm{\sigma}^{x}_{\max(A)}\bm{\sigma}^{x}_{\max(A)+1}$.

##### Disjoint blocks

When the subsystem $X$ consists of disjoint blocks $A$ and $C$ (figure 6d) we
don’t have such a simple relation as (49). The additional complication is that
the duality transformation can introduce strings between the blocks. For
example,

$\bm{\sigma}_{\max(A)}^{x}\bm{\sigma}_{\min(C)}^{x}=\prod_{\ell=\max(A)+1}^{\min(C)}\bm{\tau}_{\ell}^{x}.$
(51)

These strings are similar to the ones appearing in expressing the reduced
density matrix of disjoint blocks in terms of Jordan-Wigner fermions [68] and
we tackle the problem in a similar way.

We have to consider the enlargement of $X=A\cup C$ to the whole subsystem
between $A$ and $C$, and to the site adjacent to $C$ on the right. We thus
work with
$X^{\prime}=\\{\max(A)+1,\max(A)+2,\ldots,\min(C)-1\\}\cup\\{\max(C)+1\\}$. It
is also convenient to define the subsystems $\tilde{A}\equiv
A\cup\\{\max(A)+1\\}$, $\tilde{B}\equiv B-\\{\max(B)\\}$ and $\tilde{C}\equiv
C\cup\\{\max(C)+1\\}$ (figure 1e). Here we assume $|B|\geq 2$. Similarly to
the single block case, we consider the density matrix (48) with
$Y=\tilde{A}\tilde{C}$, which can be written as

$\tilde{\rho}_{\tilde{A}\tilde{C}}=\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{O_{\tilde{A}},O_{\tilde{C}}}\bra{\Psi}\bm{O}_{\tilde{A}}\bm{O}_{\tilde{C}}\ket{\Psi}O_{\tilde{A}}O_{\tilde{C}}\;,$
(52)

where the sum is over all possible products
$O_{\tilde{A}}=\prod_{\ell\in\tilde{A}}\tau_{\ell}^{\gamma_{\ell}}$ with
$\gamma_{\ell}\in\\{0,x,y,z\\}$ and $O_{\tilde{C}}$, defined analogously for
subsystem $\tilde{C}$. The problem of density matrix (52) is that it cannot
describe operators such as (51), that include strings between the blocks in
the dual picture. To describe them we introduce the traceless operator

$\omega_{\tilde{A}\tilde{C}}=\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{O_{\tilde{A}},O_{\tilde{C}}}\bra{\Psi}\bm{O}_{\tilde{A}}\bm{S}^{x}\bm{O}_{\tilde{C}}\ket{\Psi}(O_{\tilde{A}}O_{\tilde{C}})S^{x},$
(53)

which is similar to $\tilde{\rho}_{\tilde{A}\tilde{C}}$, but includes the
string $S^{x}\equiv\prod_{\ell\in\tilde{B}}\tau_{\ell}^{x}$
($\bm{S}^{x}\equiv\prod_{\ell\in\tilde{B}}\bm{\tau}_{\ell}^{x}$).

We find that the reduced density matrix
$\bar{\rho}_{AC}=\rho_{AC}\otimes(\bigotimes_{\ell\in
B\cup\\{\max(C)+1\\}}\mathbb{I}/2)$ is equal to

$\begin{split}\bar{\rho}_{AC}=&\frac{1}{2^{6}}\sum_{j_{1},j_{2},\ldots,j_{6}=0}^{1}\left(\prod_{m=1}^{6}{P_{m}^{j_{m}}}\right)\tilde{\rho}_{\tilde{A}\tilde{C}}\left(\prod_{m=1}^{6}{P_{m}^{j_{m}}}\right)^{\dagger}\\\
+&\frac{1}{2^{6}}\sum_{j_{1},j_{2},\ldots,j_{6}=0}^{1}(-1)^{j_{4}+j_{5}}\left(\prod_{m=1}^{6}{P_{m}^{j_{m}}}\right)\omega_{\tilde{A}\tilde{C}}\left(\prod_{m=1}^{6}{P_{m}^{j_{m}}}\right)^{\dagger}\end{split}$
(54)

where

$\begin{split}&P_{1}=\prod_{\ell\in\tilde{A}}\tau_{\ell}^{x},\qquad
P_{2}=\tau_{\mathrm{min}(\tilde{A})}^{z},\qquad
P_{3}=\tau_{\mathrm{max}(\tilde{A})}^{z},\\\
&P_{4}=\prod_{\ell\in\tilde{C}}\tau_{\ell}^{x},\qquad
P_{5}=\tau_{\mathrm{min}(\tilde{C})}^{z},\qquad
P_{6}=\tau_{\mathrm{max}(\tilde{C})}^{z}.\end{split}$ (55)

The transformations are completely analogous to the single block case and are
graphically represented in figure 6f.

The derivation of representation (54) is more pedantic than in the single
block case, but similar, so we cover just the main ingredients. It is derived
by applying properties 1 and 2 to blocks $A$ and $C$ separately. The
additional complication is that the non-locality of the transformation can
introduce strings between the blocks, as in (51). To tackle it we notice that
any product of Pauli matrices $\bm{O}_{X}$ localized (in the $\sigma$
representation) in $X=A\cup C$ and even under $\mathcal{P}_{\sigma}^{z}$
($\mathcal{P}_{\sigma}^{z}[\bm{O}_{X}]=\bm{O}_{X}$) can be written as a
product $\bm{O}_{X}=\bm{O}_{A}\bm{O}_{C}$ where $\bm{O}_{A}$ and $\bm{O}_{C}$
are localized in $A$ and $C$ respectively, and have the same parity under
$\mathcal{P}_{\sigma}^{z}$, i.e. we have that either $\bm{O}_{A}$ and
$\bm{O}_{C}$ are both even under $\mathcal{P}_{\sigma}^{z}$
($\mathcal{P}_{\sigma}^{z}[\bm{O}_{A,C}]=\bm{O}_{A,C}$) or they are both odd
($\mathcal{P}_{\sigma}^{z}[\bm{O}_{A,C}]=-\bm{O}_{A,C}$). The even operators
$\bm{O}_{A}$ and $\bm{O}_{C}$ are, by property 1, localized in $\tilde{A}$ and
$\tilde{C}$ respectively so an analogous application of property 2 to the
single block case results in the first term in (54). The second term is a
contribution from the terms where both $\bm{O}_{A}$ and $\bm{O}_{C}$ are odd.
To treat such terms we notice that they can be rewritten as
$\bm{O}_{A}=(\bm{O}_{A}\bm{\sigma}_{\ell_{A}}^{x})\bm{\sigma}_{\ell_{A}}^{x}$
and
$\bm{O}_{C}=(\bm{O}_{C}\bm{\sigma}_{\ell_{C}}^{x})\bm{\sigma}_{\ell_{C}}^{x}$,
where $\ell_{A},\ell_{C}$ is an arbitrary site belonging to $A,C$
respectively. The operators in the brackets are even under
$\mathcal{P}_{\sigma}^{z}$ so they can be dealt with by applying properties 1
and 2, while the product
$\bm{\sigma}_{\ell_{A}}^{x}\bm{\sigma}_{\ell_{C}}^{x}=\bm{\tau}_{\ell_{A}+1}^{x}\bm{\tau}_{\ell_{A}+2}^{x}\ldots\bm{\tau}_{\ell_{C}}^{x}$
is a string, common to all odd $\bm{O}_{A}$, $\bm{O}_{C}$. This gives the
second term in (54). A subtlety is that the sign factor $(-1)^{j_{4}+j_{5}}$
appears, because the operators
$\bm{\sigma}^{z}_{\max(A)}=\bm{\tau}^{z}_{\max(\tilde{A})}\bm{\tau}^{z}_{\max(\tilde{A})+1}$
and
$\bm{\sigma}^{z}_{\min(\tilde{C})-1}=\bm{\tau}^{z}_{\min(\tilde{C})-1}\bm{\tau}^{z}_{\min(\tilde{C})}$
arising from property 2 anticommute with the string $\bm{S}^{x}$.

####  III.1.4 Invariance of the tripartite information

Based on representations (49) and (54) for the reduced density matrices we now
argue that the tripartite information in the stationary states of the studied
protocols is not influenced by the Kramers-Wannier transformation.

The first thing to notice is that, in the model (17), that does not have
semilocal charges, the string order does not survive the quench (see [35] for
a general discussion), i.e. in the $\tau$ representation there is no string
order after the quench. Thus, since $\omega_{\tilde{A}\tilde{C}}$, defined in
(53), includes strings between the blocks (in the $\tau$ representation), in
the limit of large subsystems it is expected to be negligible with respect to
$\tilde{\rho}_{\tilde{A}\tilde{C}}$. Therefore, the second line in the
expression (54) for the reduced density matrix of disjoint blocks in the dual
picture can be neglected, leaving us with

$\bar{\rho}_{AC}\sim\frac{1}{2^{6}}\sum_{j_{1},j_{2},\ldots,j_{6}=0}^{1}\left(\prod_{m=1}^{6}{P_{m}^{j_{m}}}\right)\tilde{\rho}_{\tilde{A}\tilde{C}}\left(\prod_{m=1}^{6}{P_{m}^{j_{m}}}\right)^{\dagger},$
(56)

where $\sim$ stands for the limit $|A|,|B|,|C|\gg 1$.

The transformations $P$ acting on the boundaries act in the same way in the
single block case and in the case of disjoint blocks (compare Fig. 6c and Fig.
6f). Moreover, when computing different entropies appearing in the definition
(4) of the tripartite information, these boundary transformations always act
(in the dual picture) on the sites $\ell=\min
A,\max{\tilde{A}},\min{\tilde{C}},\max{\tilde{C}}$. For example, when
computing $S_{\alpha}(B)$ we have transformations that act on the sites
$\ell=\min B=\max\tilde{A}$ and $\ell=\max(B)+1=\min{\tilde{C}}$. Since the
tripartite information is constructed in such a way to cancel the
contributions of the boundaries, these transformations are not expected to
affect it. Therefore, the same tripartite information should be obtained by
starting from the density matrix

$\bar{\rho}_{A}^{\prime}\equiv\frac{1}{2}\sum_{j=0}^{1}P_{3}^{j}\tilde{\rho}_{\tilde{A}}P_{3}^{j}$
(57)

for the reduced density matrix of a single block and

$\bar{\rho}_{AC}^{\prime}\equiv\frac{1}{2^{2}}\sum_{j_{3},j_{6}=0}^{1}\left(\prod_{m=3,6}{P_{m}^{j_{m}}}\right)\tilde{\rho}_{\tilde{A}\tilde{C}}\left(\prod_{m=3,6}{P_{m}^{j_{m}}}\right)$
(58)

for the reduced density matrix of disjoint blocks.

The stationary state of both of the studied quench protocols can be written as

$\bm{\rho}_{\mathrm{GGE}}=\frac{1}{\mathcal{Z}}e^{-\bm{Q}},$ (59)

where $\mathcal{Z}$ is the normalization and $\bm{Q}$ includes the relevant
charges with appropriate Lagrange multipliers, given by eq. (43) and eq. (44).
The operator $\bm{Q}$ can be interpreted as a Hamiltonian and its eigenvalues
$E_{i}$ as energies, that come in pairs $\pm E_{i}$. Since the GGE entropy is
extensive with the system size, the energies $E_{i}$ are also expected to be
extensive. It follows that for any integer $M$ such that $1\leq M\leq\alpha-1$
the ratios

$\frac{\mathrm{tr}\left[e^{-(\alpha-2M)\bm{Q}}\right]}{\mathrm{tr}\left[e^{-\alpha\bm{Q}}\right]}=\frac{\sum_{E_{i}\geq
0}\cosh[(\alpha-2M)E_{i}]}{\sum_{E_{i}\geq 0}\cosh(\alpha E_{i})}$ (60)

are exponentially suppressed with the system size, while for $M=0,\alpha$ they
are equal to unity.

The reduced density matrices $\tilde{\rho}_{\tilde{A}}$ and
$\tilde{\rho}_{\tilde{A}\tilde{C}}$ are obtained by tracing out the GGE
density matrix in the dual picture. Let us focus on the single block case.
When the block $\tilde{A}$ is large enough the reduced density matrix
$\tilde{\rho}_{\tilde{A}}$ is expected to be similar to the GGE in the bulk.
For thermal states the entanglement Hamiltonian is at large temperatures
approximated well by the Hamiltonian of the subsystem [105]. Similarly, we
expect that the properties of the entanglement Hamiltonian in the GGE can be
captured well by the restriction of $\bm{Q}$ to the subsystem. The consequence
is that for $j_{1},\ldots,j_{\alpha-1}\in\\{0,1\\}$, such that $M$ of them are
non-zero, with $1<M<\alpha$, we have

$\frac{\mathrm{tr}\left[\rho\left(P_{m_{1}}^{j_{1}}\rho
P^{j_{1}}_{m_{1}}\right)\left(P^{j_{2}}_{m_{2}}\rho
P^{j_{2}}_{m_{2}}\right)\ldots\left(P^{j_{\alpha-1}}_{m_{\alpha-1}}\rho
P^{j_{\alpha-1}}_{m_{\alpha-1}}\right)\right]}{\mathrm{tr}\rho^{\alpha}}\sim
0\;,$ (61)

where $m_{1},\ldots,m_{\alpha-1}=3$. Here we have used the property that
$\bm{Q}$ includes only charges odd under $\mathcal{P}^{\tau}_{x}$, realization
that the transformation $\bullet\to P_{3}\bullet P_{3}$ is a finite-size
analogue of transformation $\mathcal{P}^{\tau}_{x}$, which changes the sign in
front of the odd charges, and the suppression of the ratios (60). In the case
of disjoint blocks, the transformations $P_{3},P_{6}$, appearing in (58), act
only on a part of the subsystem, but their support is still extensive with the
size of one of the subsystems and we expect the ratios in (61) with
$m_{1},\ldots,m_{\alpha-1}\in\\{3,6\\}$ to be suppressed for similar reasons.

From (61) it follows

$\mathrm{tr}\left[\left(\bar{\rho}_{A}^{\prime}\right)^{\alpha}\right]\sim
2^{1-\alpha}\mathrm{tr}\left[\left(\tilde{\rho}_{\tilde{A}}\right)^{\alpha}\right],\qquad\mathrm{tr}\left[\left(\bar{\rho}_{AC}^{\prime}\right)^{\alpha}\right]\sim
2^{2(1-\alpha)}\mathrm{tr}\left[\left(\tilde{\rho}_{\tilde{A}\tilde{C}}\right)^{\alpha}\right]\;,$
(62)

where in the disjoint blocks case the squared prefactor appears because in eq.
(58) we have two transformations. Now, since the prefactors in (62) cancel in
the tripartite information it follows that the same tripartite information is
obtained if we use $\tilde{\rho}_{\tilde{A}}$ instead of $\bar{\rho}_{A}$
(analogously for other single block contributions to the tripartite
information in (4)) and $\tilde{\rho}_{\tilde{A}\tilde{C}}$ instead of
$\bar{\rho}_{AC}$. Since
$\tilde{\rho}_{\tilde{A}},\tilde{\rho}_{\tilde{A}\tilde{C}}$ are in
correspondence with the density matrices
$\rho_{\tilde{A}}^{\tau},\rho_{\tilde{A}\tilde{C}}^{\tau}$, defined by (19),
we can simply use the latter as well. The conclusion is that for the
stationary states of the studied quench protocols the tripartite information
is not influenced by the Krammers-Wannier transformation: the differences in
the Rényi entropies arising from the differences between the density matrices
(18) and (19) cancel in the tripartite information of three large adjacent
blocks. The same conclusion is expected to hold quite generically for
stationary states described by semilocal Gibbs ensembles. For arbitrary states
the conclusion does not need to hold. For example, for states with non-
negligible string contributions in (54) (in the $\tau$ representation) the
presented argument breaks right at the start.

We note that the argument presented here is largely intuitive, despite being
also quite technical. However, exact numerical results for the second Rényi
entropy, obtained using fermionic techniques, provide an implicit check for
the argument’s validity.

###  III.2 Jordan-Wigner transformation

The XY chain (17) is a quadratic form in Majorana fermions, defined as

$\bm{a}_{2\ell-1}\equiv\bm{a}^{x}_{\ell}=\left(\prod_{j<\ell}\bm{\tau}^{z}_{j}\right)\bm{\tau}^{x}_{\ell},\qquad\bm{a}_{2\ell}\equiv\bm{a}^{y}_{\ell}=\left(\prod_{j<\ell}\bm{\tau}^{z}_{j}\right)\bm{\tau}^{y}_{\ell}\;.$
(63)

The Majorana fermions are self-adjoint and satisfy the algebra

$\\{\bm{a}_{\ell},\bm{a}_{n}\\}=2\delta_{\ell n}\bm{I}\;.$ (64)

The XY chain Hamiltonian can due to translational invariance be written as

$\bm{H}=\frac{1}{4}\sum_{\ell,n=-\infty}^{\infty}\begin{pmatrix}\bm{a}^{x}_{\ell}&\bm{a}^{y}_{\ell}\end{pmatrix}\int_{-\pi}^{\pi}\tfrac{{\rm
d}k}{2\pi}e^{i(\ell-n)k}\mathcal{H}(k)\begin{pmatrix}\bm{a}^{x}_{n}\\\
\bm{a}^{y}_{n}\end{pmatrix},$ (65)

where the $2\times 2$ matrix $\mathcal{H}(k)$, called the Hamiltonian symbol
[106], generates the couplings through its Fourrier coefficients. The symbol
is given by

$\displaystyle\mathcal{H}(k)=2(J_{x}-J_{y})\sin(k)\,\sigma^{x}-2(J_{x}+J_{y})\cos(k)\,\sigma^{y}.$
(66)

The positive eigenvalues of the symbol, given by

$\varepsilon(k)=2\sqrt{J_{x}^{2}+J_{y}^{2}+2J_{x}J_{y}\cos(2k)}$ (67)

are the energies of the quasiparticle excitations of the model.

In a translationally invariant state $\ket{\Psi}$ the correlation matrix

$\Gamma_{\ell,n}=\delta_{\ell,n}I-\bra{\Psi}\bm{a}_{\ell}\bm{a}_{n}\ket{\Psi}$
(68)

can be expressed in terms of its symbol $\Gamma(k)$, which is a $2\times 2$
matrix, as

$\Gamma_{2\ell+i,2n+j}=\int_{-\pi}^{\pi}\tfrac{{\rm
d}k}{2\pi}e^{i(\ell-n)k}\Gamma_{ij}(k)\;,\qquad\ell,n\in\mathbb{Z},\
i,j\in\\{1,2\\}.$ (69)

In the dual picture the all-spin-up quench protocol is a homogenous quench
protocol, while the flipped-spin protocol is a bipartitioning protocol.
Accordingly, the correlation matrix symbol describing their stationary states
has already been derived in the literature (see e.g. [35, 2]). For the
stationary state of the all-spin-up quench protocol it reads

$\textrm{protocol \ref{quench protocol
1}:}\qquad\Gamma(k)=\frac{\mathrm{tr}[\sigma^{y}\mathcal{H}(k)]\mathcal{H}(k)}{2\varepsilon^{2}(k)}\,,$
(70)

which can be obtained simply as a time average. The non-equilibrium stationary
state following the quench from the domain wall state is also translationally
invariant, notwithstanding that the initial state is not, and its correlation
matrix symbol reads

$\textrm{protocol \ref{quench protocol
2}:}\qquad\Gamma(k)=-\mathrm{sgn}[v(k)]\frac{\mathrm{tr}[\sigma^{y}\mathcal{H}(k)]\mathbb{I}}{2\varepsilon(k)}\,,$
(71)

where $v(k)=\varepsilon^{\prime}(k)$ is the velocity of the quasiparticle
excitations. The symbol (71) is derived using the generalized hydrodynamic
equation, which provides the connection between stationary states at different
rays $\zeta=d/t$, where $d$ is the distance from the position of the domain
wall in the initial state.

The spin correlation functions of section II.1 and the string order parameters
of section II.2 are computed straightforwardly from the correlation matrix by
expressing the operators in question in terms of fermions (63), where we use
also some simplifying properties of the correlation matrix in question, as
discussed in appendix A. The Rényi entanglement entropy of a single block is
obtained, as already done in [35], by expressing the density matrix (49) in
terms of fermionic Gaussians, i.e. exponentials of quadratic forms of
fermionic operators. Similarly, the entanglement entropy of disjoint blocks is
obtained by expressing the density matrix (54) in terms of Gaussians. The
procedure is pedantic, but the methods essentially do not go beyond ref. [68]
dealing with the entanglement entropy of disjoint blocks for the model (17).
The procedure of expressing the density matrices in terms of fermionic
Gaussians is reported in appendix C, while the exact formulas for computing
the second Rényi entropy are given in appendix D.

##  IV Conclusions

We have shown that the stationary state of the translationally invariant all-
spin-up quench protocol and the stationary state reached in the flipped-spin
quench protocol (NESS) have different properties. The former exhibits
exponential decay of spatial correlations and zero tripartite information of
three large adjacent subsystems, while in the latter there are correlations
that decay only algebraically and the tripartite information is non-zero. On
the other hand, string order is present in both cases. It would be interesting
to understand further the meaning of different string order parameters.
Importantly, a single spin flip in the initial state is responsible for
differences in the behavior of the spatial correlations and the tripartite
information between the two stationary states. This result complements the
recents findings that a spin flip in the initial state affects the
magnetization [46] and the subleading term in the single block entanglement
entropy [35] in the stationary state. The mechanism responsible for all this
phenomenology is the existence of semilocal charges, whose expectation value
can be strongly affected by localized perturbations. We have also derived
explicitly the expressions for generalized Gibbs ensembles describing the
stationary state of the all-spin-up and the flipped-spin quench protocol. The
set of semilocal charges associated with nonzero Lagrange multipliers is
completely different for the two protocols.

We remark that if the initial state contains several flipped spins instead of
just one the properties of the stationary state will strongly depend on the
parity of the number of flipped spins. Flipping an even number of spins in the
initial state is expected to yield a stationary state with the same properties
as in the all-spin-up quench protocol, while an odd number of flipped spins is
expected to yield the discovered phenomenology of the flipped-spin quench
protocol. The reason is that flipping an even number of spins does not excite
semilocal charges, while an odd number of flips does.

The tripartite information in the stationary state of the flipped-spin quench
protocol is interesting on its own. Although the initial state is a product
state and the time evolution is with a gapped Hamiltonian, the tripartite
information in the stationary state shares universal properties with conformal
field theory. It depends on the lengths of the subsystems only through the
cross ratio (8). It would be nice to have an explanation for this phenomenon
independent of the details of the computation. Let us note that, although the
obtained formula for the tripartite information, given by eq. (29), is
independent of the model parameters, such an independence would be lost in the
non-interacting generalization of the model, the dual generalized XY chain
introduced in ref. [35], for which the methods developed in this work still
apply. Based on the results of Paper II, we conclude, on the other hand, that
universality, i.e. the dependence on the configuration only through the cross
ratio, should still hold. Another aspect of our result is that, as in
protocols studied in Papers I and II, the stationary state exhibits nonzero
residual tripartite information. Namely, in the limit of small size of
subsystem $B$ with respect to the size of subsystems $A,C$ (but still much
larger than the lattice spacing) the tripartite information is nonzero and
equal to $-\log 2$. This value should be contrasted to the zero value found in
equilibrium at any temperature, irrespectively of criticality, or in other
non-equilibrium settings, such as after quenches from ground states of gapped
Hamiltonians. It is an interesting open question whether universality would
survive the generalization to interacting models and what would be the value
of the residual tripartite information in this case. Finally, we mention that
the tripartite information is negative, which in other contexts [84, 76] has
been interpreted as an indication that quantum entanglement dominates over
classical correlations. Maybe computing other entanglement measures, such as
entanglement negativity [107, 108, 109], could shed light on the nature of
this phase.

##  V Acknowledgments

This work was supported by the European Research Council under the Starting
Grant No. 805252 LoCoMacro. I thank Maurizio Fagotti for providing me the idea
for this work, collaboration on related topics and useful discussions. I also
thank Saverio Bocini for useful discussions.

##  Appendix A Correlation functions of Majorana fermions

Majorana correlation functions, determined by the symbols (70) and (71), have
some simplifying properties. We assume $0<|J_{y}|<|J_{x}|$. In the stationary
state of quench protocol 1 the correlations are

$\begin{split}&\braket{\bm{a}_{\ell+r}^{x}\bm{a}_{\ell}^{y}}\\\
&=\begin{cases}0,&r<0\textrm{ or }r\in 2\mathbbm{Z}+1\\\
\frac{i}{2}\big{(}1+\frac{J_{y}}{J_{x}}\big{)},&r=0\\\
\frac{i}{2}\left(1-\frac{J_{y}^{2}}{J_{x}^{2}}\right)\left(-\frac{J_{y}}{J_{x}}\right)^{\frac{r}{2}-1},&r>0\textrm{
and }r\in 2\mathbb{Z}\end{cases}\;,\end{split}$ (72)

as obtained by expanding (70) in the Fourier series, where it is convenient to
exploit the formula for the geometric series. The correlations vanish for all
odd $r$ and for negative even $r$. Moreover the correlations decay
exponentially with $r$, related to the fact that symbol (70) is smooth. Also,
trivially, the correlations
$\braket{\bm{a}_{\ell+r}^{\alpha}\bm{a}_{\ell}^{\beta}}$ for $\alpha=\beta$
vanish for all $r\neq 0$.

In the stationary state of quench protocol 2 the correlations
$\braket{\bm{a}_{\ell+r}^{\alpha}\bm{a}_{\ell}^{\beta}}$ vanish for
$\alpha\neq\beta$, while for $\alpha=\beta$ they vanish for nonzero even $r$.
The remaining ones decay to zero algebraically with $r$, since (71) is
discontinuous. Using partial integration we get

$\braket{\bm{a}_{\ell+2r-1}^{x}\bm{a}_{\ell}^{x}}=\braket{\bm{a}_{\ell+2r-1}^{y}\bm{a}_{\ell}^{y}}\simeq-\mathrm{sgn}(J_{y})\frac{i}{\pi
r}.$ (73)

In computing the string order parameter we will also use the exact result for
$r=1$

$\braket{\bm{a}_{\ell+1}^{x}\bm{a}_{\ell}^{x}}=\braket{\bm{a}_{\ell+1}^{y}\bm{a}_{\ell}^{y}}=-\mathrm{sgn}(J_{y})\frac{i}{\pi}\left(1+\frac{J_{y}}{J_{x}}\right)\;,$
(74)

obtained by basic integral manipulations.

The spin correlation functions in section II.1 and the string order parameters
in section II.2 are computed straightforwardly using the Kramers-Wannier
transformation (15) and the listed properties of Majorana correlators. In
particular, the string operators read

$\displaystyle\prod_{\ell=-r}^{r}\bm{\sigma}^{z}_{\ell}=\bm{\tau}^{z}_{-r}\bm{\tau}^{z}_{r+1}=-\bm{a}_{-r}^{x}\bm{a}_{-r}^{y}\bm{a}_{r+1}^{x}\bm{a}_{r+1}^{y}\;,$
(75)
$\displaystyle\bm{\sigma}_{-r-3}^{x}\bm{\sigma}_{-r-2}^{z}\bm{\sigma}_{-r-1}^{y}\left(\prod_{\ell=-r}^{r}\bm{\sigma}^{z}_{\ell}\right)\bm{\sigma}^{y}_{r+1}\bm{\sigma}_{r+2}^{z}\bm{\sigma}_{r+3}^{x}=\bm{\tau}_{-r-2}^{y}\bm{\tau}_{-r-1}^{x}\bm{\tau}_{r+2}^{x}\bm{\tau}_{r+3}^{y}=\bm{a}_{-r-2}^{x}\bm{a}_{-r-1}^{x}\bm{a}_{r+2}^{y}\bm{a}_{r+3}^{y}\;.$
(76)

##  Appendix B Procedure for finding the Lagrange multipliers

In this section we show how the Lagrange multipliers in the generalized Gibbs
ensembles describing the studied stationary states are obtained, based on
methods developed in [95, 103].

The one-site shift invariant charges of the quantum XY chain (17) given in
section III.1.1 can be written as a quadratic form in Majorana fermions

$\bm{Q}^{(n,\pm)}=\frac{1}{4}\sum_{\ell,n=-\infty}^{\infty}\begin{pmatrix}\bm{a}^{x}_{\ell}&\bm{a}^{y}_{\ell}\end{pmatrix}\int_{-\pi}^{\pi}\tfrac{{\rm
d}k}{2\pi}e^{i(\ell-n)k}\mathcal{Q}^{(n,\pm)}(k)\begin{pmatrix}\bm{a}^{x}_{n}\\\
\bm{a}^{y}_{n}\end{pmatrix},$ (77)

specified by the symbol $\mathcal{Q}^{(n,\pm)}(k)$, which is a $2\times 2$
matrix function. The charges are constructed so that the symbol
$\mathcal{Q}^{(n,\pm)}(k)$ commutes with the Hamiltoian symbol
$\mathcal{H}(k)$, because this implies that charges commute with the
Hamiltonian ($[\bm{Q}^{(n,\pm)},\bm{H}]=0$). Moreover, the symbol is chosen to
have a finite number of Fourier coefficients, to ensure the locality of the
charge. In particular, the charges of section III.1.1 are given by the symbols

$\displaystyle\mathcal{Q}^{(n,+)}(k)=\cos(nk)\mathcal{H}(k)\;,$ (78)
$\displaystyle\mathcal{Q}^{(n,-)}(k)=\sin[(n+1)k]\mathbb{I}\;,$ (79)

for $n\in\mathbb{N}_{0}$.

A translationally invariant fermionic Gaussian can be written as

$\bm{\rho}=\frac{1}{\mathcal{Z}}\exp\left[\frac{1}{4}\sum_{\ell,n=-\infty}^{\infty}\begin{pmatrix}\bm{a}^{x}_{\ell}&\bm{a}^{y}_{\ell}\end{pmatrix}\int_{-\pi}^{\pi}\frac{{\rm
d}k}{2\pi}e^{i(\ell-n)k}\mathcal{W}(k)\begin{pmatrix}\bm{a}^{x}_{n}\\\
\bm{a}^{y}_{n}\end{pmatrix}\right],$ (80)

where $\mathcal{Z}$ is the normalization and $\mathcal{W}(k)$ is related to
the symbol of the correlation matrix (see e.g. (86)) through

$\mathcal{W}(k)=2\ \mathrm{artanh}\left[\Gamma(k)\right]\;.$ (81)

On the other hand, $\mathcal{W}(k)$ is related to the Lagrange multipliers in
the GGE (42) through

$\mathcal{W}(k)=-\sum_{n=0}^{\infty}\left(\lambda^{(n,+)}\mathcal{Q}^{(n,+)}(k)+\lambda^{(n,-)}\mathcal{Q}^{(n,-)}(k)\right)\;.$
(82)

Thus, Lagrange multipliers $\\{\lambda^{(n,+)},\lambda^{(n,-)}\\}_{n}$ are
obtained from the correlation matrix symbol, given by (70) and (71) for the
studied quench protocols, by comparing (81) and (82). In this way the Lagrange
multipliers in (43) and (44) are obtained, where the given integrals are
simplifications of the ones from $0$ to $2\pi$ arising from the Fourier
coefficients.

##  Appendix C Reduced density matrix in the fermionic picture

The usefulness of representation (49) and (54) is that they enable us to go to
the fermionic language. The goal of this section is to use these
representations to express the density matrices in a useful form for computing
the Rényi entropies, which is the subject of appendix D. Assuming $\min A=1$
without loss of generality, the mapping starts by defining the Majorana
fermions restricted to $\tilde{A}\tilde{B}\tilde{C}$,

$a_{2\ell-1}\equiv
a^{x}_{\ell}=\left(\prod_{j=1}^{\ell-1}\tau^{z}_{j}\right)\tau^{x}_{\ell},\qquad
a_{2\ell}\equiv
a^{y}_{\ell}=\left(\prod_{j=1}^{\ell-1}\tau^{z}_{j}\right)\tau^{y}_{\ell}\;,$
(83)

for $\ell=1,2,\ldots,|\tilde{A}\tilde{B}\tilde{C}|$, which are analogous to
those defined on the whole chain and satisfy analogous algebra
$\\{a_{\ell},a_{n}\\}=2\delta_{\ell n}I$.

##### Gaussians of fermionic operators

In the $2^{d}$ dimensional Fock space associated to Majorana fermions
$a_{1},a_{2},a_{3},\ldots,a_{2d}$ (with anticommutation relations
$\\{a_{\ell},a_{n}\\}=2\delta_{\ell,n}$) a generic normalized Gaussian
operator $\rho(\Gamma)$ reads

$\rho(\Gamma)=\frac{1}{\mathcal{Z}}\exp\left(\frac{\vec{a}^{\dagger}W\vec{a}}{4}\right)\,,$
(84)

where $\vec{a}^{\dagger}=(a_{1},a_{2},\ldots,a_{2d})$ and $W$ is a $2d\times
2d$ dimensional antisymmetric matrix with complex entries. The Gaussian is
completely specified by its two-point correlation matrix, given by

$\Gamma_{j,\ell}=\delta_{j,\ell}-\mathrm{tr}\left[\rho(\Gamma)\
a_{j}a_{\ell}\right],\qquad 1\leq j,l\leq 2d,$ (85)

which is related to matrix $W$ in the exponent through the formula [68]

$\Gamma=\tanh\left(\frac{W}{2}\right).$ (86)

The normalization is given by

$\mathcal{Z}=\pm\left[\det\left(e^{\frac{W}{2}}+e^{-\frac{W}{2}}\right)\right]^{\frac{1}{2}}=\pm\left[\det\left(\frac{I-\Gamma}{2}\right)\right]^{-\frac{1}{2}},$
(87)

where the sign in front of the square root is not specified by the formula and
depends on further details of $\Gamma$. When $\rho$ is a density matrix we
have that $\Gamma$ and $W$ are Hermitian so the sign is positive. We will
encounter non-Hermitian $\Gamma$ and $W$ only in the case of disjoint blocks,
analogously to ref. [68].

##### Single block

In the case of a single block, it is an established fact [110, 111, 22, 112]
that the reduced density matrices such as $\tilde{\rho}_{\tilde{A}}$ in (49),
defined by (48), are Gaussians, i.e. exponentials of quadratic forms of
Majorana fermions, determined by the correlation matrix. As discussed in [35]
and in appendix D, each term in the sum in (49) is still a Gaussian, with the
correlation matrix depending on the transformations $P$. The Gaussian
structure can then be exploited for the computation of the Rényi entropies.

##### Disjoint blocks

In the case of disjoint blocks we have derived representation (54) for the
reduced density matrix. Each term in the sum in (54) can be written as a sum
of Gaussians. The Gaussian structure can then be exploited to compute the
Rényi entropies.

Expressing $\tilde{\rho}_{\tilde{A}\tilde{C}}$, appearing in (54), in terms of
Gaussians is a procedure identical to expressing the reduced density matrix in
the XY chain in terms of Gaussians, as done in ref. [68]. We review it briefly
because we use the same techniques for $\omega_{\tilde{A}\tilde{C}}$. We have

$\begin{split}\tilde{\rho}_{\tilde{A}\tilde{C}}&=\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{\mathrm{even}\;F_{\tilde{A}},F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}\right)^{\dagger}\\\
&+\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{\mathrm{odd}\;F_{\tilde{A}},F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}S^{z}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}S^{z}\right)^{\dagger}\;,\\\
\end{split}$ (88)

where the first (second) sum is over all $F_{\tilde{A}}$ and $F_{\tilde{C}}$
that are both a product of an even (odd) number fermions in $\tilde{A}$ and
$\tilde{C}$ respectively, i.e. $F_{\tilde{A}}$ and $F_{\tilde{C}}$ in the
first (second) sum commutes (anticommutes) with
$\Pi_{\tau}^{z}\equiv\prod_{\ell\in A}\tau_{\ell}^{z}$. Explicitly,

$F_{\tilde{A}}=\prod_{\ell\in\tilde{A}}(a_{\ell}^{x})^{j^{x}_{\ell}}(a_{\ell}^{y})^{j^{y}_{\ell}},\qquad
F_{\tilde{C}}=\prod_{\ell\in\tilde{C}}(a_{\ell}^{x})^{j^{x}_{\ell}}(a_{\ell}^{y})^{j^{y}_{\ell}}\;,$
(89)

where indices $j_{k}^{x,y}$ are either $0$ or $1$ and the order of operators
in the product does not matter in (88). In the second term in (88) the string
$S^{z}$, where we denote
$S^{\gamma}=\prod_{\ell\in\tilde{B}}\tau^{\gamma}_{\ell}$ for
$\gamma=0,x,y,z$, appears in order to cancel the string that arises due to the
non-locality of the Jordan-Wigner transformation. The two different sums can
be expressed conveniently using the (anti)commutation relation with the string
over the block $\tilde{A}$

$P_{0}\equiv S_{\tilde{A}}^{z}=\prod_{\ell\in\tilde{A}}\tau_{\ell}^{z}\;.$
(90)

We have

$\tilde{\rho}_{\tilde{A}\tilde{C}}=\frac{\rho_{0}+S_{\tilde{A}}^{z}\rho_{0}S_{\tilde{A}}^{z}}{2}+S^{z}\frac{\rho_{z}-S_{\tilde{A}}^{z}\rho_{z}S_{\tilde{A}}^{z}}{2},$
(91)

where we denote

$\rho_{\gamma}=\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{F_{\tilde{A}},F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}S^{\gamma}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}\right)^{\dagger}$
(92)

for $\gamma=0,x,y,z$. Here the sum is over all possible products of fermions
in $\tilde{A}$ and $\tilde{C}$.

The operator $\omega_{\tilde{A}\tilde{C}}$, appearing in (54), can be
expressed in terms of fermions in a similar way. For even $|\tilde{B}|$ it is
completely analogous to $\tilde{\rho}_{\tilde{A}\tilde{C}}$,

$\begin{split}\omega_{\tilde{A}\tilde{C}}&=S^{x}\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{\mathrm{even}\;F_{\tilde{A}},F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}S^{x}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}\right)^{\dagger}\\\
&+S^{x}\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{\mathrm{odd}\;F_{\tilde{A}},F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}S^{z}S^{x}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}S^{z}\right)^{\dagger}\;,\\\
\end{split}$ (93)

where the first (second) sum is again over all $F_{\tilde{A}}$ and
$F_{\tilde{C}}$ that are both a product of an even (odd) number of fermions
and $S^{z}$ stems from the non-locality of the Jordan-Wigner transformation.
For odd $|\tilde{B}|$ it is slightly different,

$\begin{split}\omega_{\tilde{A}\tilde{C}}&=S^{x}\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{\mathrm{odd}\;F_{\tilde{A}},\;\mathrm{even}\;F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}S^{x}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}\right)^{\dagger}\\\
&+S^{x}\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{\mathrm{even}\;F_{\tilde{A}},\;\mathrm{odd}\;F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}S^{z}S^{x}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}S^{z}\right)^{\dagger}\;.\\\
\end{split}$ (94)

Now the number of fermions in $F_{\tilde{A}}$ and $F_{\tilde{C}}$ has opposite
parities. The string $S^{z}$ has to be included when $F_{\tilde{C}}$ consists
of an odd number of fermions or, consequently, when $F_{\tilde{A}}$ consists
of an even number of fermions.

Overall, we have

$\omega_{\tilde{A}\tilde{C}}=S^{x}\frac{\rho_{x}+(-1)^{|\tilde{B}|}S_{\tilde{A}}^{z}\rho_{x}S_{\tilde{A}}^{z}}{2}+S^{y}\frac{\rho_{y}-(-1)^{|\tilde{B}|}S_{\tilde{A}}^{z}\rho_{y}S_{\tilde{A}}^{z}}{2},$
(95)

where the string $S^{y}$ is a result of multiplying $S^{x}$, that comes from
the non-local character of the Kramers-Wannier transformation, with $S^{z}$,
that comes from the non-locality of the Jordan-Wigner transformation. Finally,
(95) can be expressed in terms of Gaussians, since the operators (92) can be
expressed in terms of Gaussians, as covered in the following.

##### Operators with strings as Gaussians

Let us consider operators defined in (92). We note that they satisfy

$S^{\gamma}\rho_{\gamma}=\frac{1}{2^{|\tilde{A}\tilde{B}\tilde{C}|}}\sum_{F_{\tilde{A}},F_{\tilde{C}}}\bra{\Psi}F_{\tilde{A}}F_{\tilde{C}}\mathcal{S}^{\gamma}\ket{\Psi}\left(F_{\tilde{A}}F_{\tilde{C}}\mathcal{S}^{\gamma}\right)^{\dagger},$
(96)

where we have the freedom of replacing the string between the blocks
$S^{\gamma}$ by
$\mathcal{S}^{\gamma}=F^{\prime}_{\tilde{A}}F^{\prime}_{\tilde{C}}S^{\gamma}$
for any fixed $F_{\tilde{A}}^{\prime}$ and $F_{\tilde{C}}^{\prime}$, that are
a product of fermionic operators in $\tilde{A}$ and $\tilde{C}$ respectively.
The exact choice is discussed in appendix E. This replacement is allowed since
it corresponds simply to changing the dummy summation index. Because of this
freedom we can choose $\mathcal{S}^{\gamma}$ to be a product of an even number
of Majorana fermions (irrespectively of the parity of $|\tilde{B}|$) and
therefore a Gaussian. Moreover, choosing $\mathcal{S}^{\gamma}$ such that it
has a non-zero expectation value $\braket{\mathcal{S}^{\gamma}}$, the
expectation values

$\braket{F_{\tilde{A}}F_{\tilde{C}}\mathcal{S}^{\gamma}}=\braket{\mathcal{S}^{\gamma}}\mathrm{tr}\left[F_{\tilde{A}}F_{\tilde{C}}\rho(\Gamma_{\gamma})\right],$
(97)

are given by the Gaussian $\rho(\Gamma_{\gamma})$ specified by the correlation
matrix

$\left(\Gamma_{\gamma}\right)_{j,\ell}=\delta_{j\ell}-\frac{\braket{a_{f(j)}a_{f(\ell)}\mathcal{S}^{\gamma}}}{\braket{\mathcal{S^{\gamma}}}}\;,\qquad
j,\ell=1,2,\ldots,2|\tilde{A}\tilde{C}|$ (98)

where

$f(j)=\begin{cases}&j,\qquad 1\leq j\leq|2\tilde{A}|\\\
&2[\min{(\tilde{C})}-|\tilde{A}|]+j-2,\qquad 2|\tilde{A}|+1\leq j\leq
2|\tilde{A}\tilde{C}|\end{cases}$ (99)

is just a compact notation to work with the indices related to disjoint
blocks. Here we have used the Wick theorem, the property that the product of
two Gaussians ($\mathcal{S}^{\gamma}$ and $\rho(\Gamma)$) is a Gaussian [68]
and we have divided by the expectation value of the string
$\braket{\mathcal{S}^{\gamma}}$ to ensure the normalization
$\mathrm{tr}\left[\rho(\Gamma_{\gamma})\right]=1$. The correlation matrix (98)
and, accordingly, the matrix $W$ in (84), are Hermitian for $\gamma=0$, but in
general this is not the case. In general they are complex antisymmetric
matrices. We note that Wick theorem holds also in this more general case. It
is desirable to have a proof of Wick theorem as such directly in the Majorana
fermions formalism we work with so we provide it in appendix G.

We have thus shown that $S^{\gamma}\rho_{\gamma}$ can be expressed in terms of
a normalized Gaussian,

$S^{\gamma}\rho_{\gamma}=\left(\mathcal{S}^{\gamma}\right)^{\dagger}\braket{\mathcal{S}^{\gamma}}\rho(\Gamma_{\gamma})\;.$
(100)

The correlation matrix (98) could, in principle, be evaluated from its
definition, but since this procedure would require the evaluation of the
pfaffian for each matrix element, in appendix F we derive an alternative
expression (a similar trick was used in [68] to avoid singularities). Suppose
that the string is given by $\mathcal{S}^{\gamma}=a_{i_{1}}a_{i_{2}}\ldots
a_{i_{m}}$ for some positive even integer $m$. We find

$\Gamma_{\gamma}=\Gamma_{0}-UN^{-1}U^{\mathrm{T}},$ (101)

where

$U_{\ell,n}=\braket{a_{f(\ell)}a_{i_{n}}}\qquad\ell=1,2,\ldots,2|\tilde{A}\tilde{C}|,\;n=1,2,\ldots
m$ (102)

and

$N_{\ell,n}=\braket{a_{i_{\ell}}a_{i_{n}}}\;,\qquad\ell,n=1,2,\ldots,m\;.$
(103)

Note that the squared absolute value of the expectation value of the string is
given by

$|\braket{\mathcal{S}^{\gamma}}|^{2}=|\det N|.$ (104)

##  Appendix D The second Rényi entropy

In the previous section we have derived fermionic representations of reduced
density matrices. In this section we use these expressions to find a way for
computing the Rényi entropies exactly. For simplicity, we focus only on the
second Rényi entropy.

##### Products of Gaussians

First we cover some basic formulas. For the computation of the second Rényi
entropy we only need the knowledge of the trace of the product of two
normalized Gaussians, defined in (84). It has been argued in [68] that for two
normalized Gaussians $\rho(\Gamma_{1})$ and $\rho(\Gamma_{2})$, with
correlation matrices $\Gamma_{1}$ and $\Gamma_{2}$, the trace is given by

$\mathrm{tr}\left[\rho(\Gamma_{1})\rho(\Gamma_{2})\right]\equiv\mathcal{G}(\Gamma_{1},\Gamma_{2})\equiv\prod_{\\{\lambda\\}/2}\lambda\;,$
(105)

where the product is over all eigenvalues $\lambda$ of the matrix
$(I+\Gamma_{1}\Gamma_{2})/2$ with degeneracy, that is always even, reduced by
half. Note that up to an unresolved sign the following formula holds (see also
[113])

$\mathcal{G}(\Gamma_{1},\Gamma_{2})=\pm\sqrt{\det\left(\frac{I+\Gamma_{1}\Gamma_{2}}{2}\right)}\;.$
(106)

##### Single block

The expression for computing exactly the second Rényi entropy of a single
block has already been derived in [35]. We revisit it here since we need the
result for computing the tripartite information, and as a warm-up for the case
of disjoint blocks. In the case of a single block the density matrix
$\tilde{\rho}_{\tilde{A}}$ in (49) is a Gaussian,

$\tilde{\rho}_{\tilde{A}}=\rho(\Gamma)$ (107)

where the correlation matrix is given by
$\Gamma=\delta_{\ell,n}-\braket{a_{\ell}a_{n}}$ for
$\ell,n=1,2,\ldots,2|\tilde{A}|$. The transformations $P_{j}$ in (49) preserve
the Gaussianity of the state but modify the correlation matrix. The latter can
be found simply by studying the effects on the two-point products of Majorana
fermions. Namely, the operator $P_{j}\rho(\Gamma)P_{j}$, for $j=1,2,3$, is the
(normalized) Gaussian $\rho(V_{j}\Gamma V_{j})$, where $V_{j}$ are diagonal
matrices of size $2|\tilde{A}|$ defined by

$\left(V_{1}\right)_{\ell\ell}=\begin{cases}-1,&\ell=1,2\\\
1,&\mathrm{otherwise}\end{cases}\;,\qquad\left(V_{2}\right)_{\ell\ell}=\begin{cases}-1,&\ell=2|\tilde{A}|-1,2|\tilde{A}|\\\
1,&\mathrm{otherwise}\end{cases}\;,\qquad\left(V_{3}\right)_{\ell\ell}=(-1)^{\lfloor\frac{\ell}{2}\rfloor}\;,$
(108)

where $\lfloor.\rfloor$ stands in the standard way for the floor function. The
Rényi entropy now follows by applying the formula (105) for the trace of a
product of two Gaussians. Introducing the compact notation

$\mathcal{V}[j_{1},j_{2},j_{3}](\Gamma)=\left(V_{1}^{j_{1}}V_{2}^{j_{2}}V_{3}^{j_{3}}\right)^{\mathrm{T}}\Gamma\left(V_{1}^{j_{1}}V_{2}^{j_{2}}V_{3}^{j_{3}}\right)\;,$
(109)

we have

$S_{2}(A)=2\log
2-\log\left[\sum_{j_{1},j_{2},j_{3}=0}^{1}\mathcal{G}\Big{(}\Gamma,\mathcal{V}[j_{1},j_{2},j_{3}](\Gamma)\Big{)}\right],$
(110)

where $\mathcal{G}$ is given by (106) with the $+$ sign. Note that the summand
$2\log 2$ comes from the second term in (47) (for $|X^{\prime}|=1$) and from
the factor $1/2^{3}$ in (49).

##### Disjoint blocks

The second Rényi entropy of disjoint blocks can be computed using eq. (54),
relating the $\sigma$ and $\tau$ representation, together with eq. (91), (95),
that allow us to express the density matrix as a sum of Gaussian fermionic
operators. In order to cancel the string operators without additional phase
factors, using
$\left(\mathcal{S}^{\gamma}\right)^{\dagger}\mathcal{S}^{\gamma}=1$, we are
going to start from the trivial relation
$\mathrm{tr}(\bar{\rho}_{AC}^{2})=\mathrm{tr}(\bar{\rho}_{AC}\bar{\rho}_{AC}^{\dagger})$.
Using the property

$\mathrm{tr}\left(\frac{\rho_{1}\pm P\rho_{1}P}{2}\frac{\rho_{2}\pm
P\rho_{2}P}{2}\right)=\mathrm{tr}\left(\rho_{1}\frac{\rho_{2}\pm
P\rho_{2}P}{2}\right),$ (111)

that holds for any operators $\rho_{1},\rho_{2}$ given some Hermitian
involution $P$, we also have the simplification

$\mathrm{tr}\left(\bar{\rho}_{AC}^{2}\right)=\mathrm{tr}\left(\tilde{\rho}_{\tilde{A}\tilde{C}}\bar{\rho}_{AC}^{\dagger}+\omega_{\tilde{A}\tilde{C}}\bar{\rho}_{AC}^{\dagger}\right)$
(112)

It is crucial that each of the four terms appearing in (91) and (95)
multiplied by the string is a Gaussian. Similarly to the single block case,
the information on the effects of the transformations $S_{\tilde{A}}^{z}$ and
$P_{j}$, for $j=1,2,\ldots,6$, appearing in (54), on Gaussian operators
$\rho_{\gamma}$, for $\gamma=0,x,y,z$, are encoded in the effects on the
correlation functions. Namely, for a normalized Gaussian $\rho(\Gamma)$ the
operator $P_{j}\rho(\Gamma)P_{j}$ is for $j=0,1,\ldots,6$ a normalized
Gaussian $\rho(V_{j}\Gamma V_{j})$, where $V_{j}$ are diagonal matrices of
size $2|\tilde{A}\tilde{C}|$ defined by

$\displaystyle\left(V_{0}\right)_{\ell\ell}=\begin{cases}-1,&\ell\leq
2|\tilde{A}|\\\
1,&\mathrm{otherwise}\end{cases}\;,\qquad\left(V_{3}\right)_{\ell\ell}=\begin{cases}(-1)^{\lfloor\frac{\ell}{2}\rfloor},&\ell\leq
2|\tilde{A}|\\\
(-1)^{|\tilde{A}|},&\mathrm{otherwise}\end{cases}\;,\qquad\left(V_{6}\right)_{\ell\ell}=\begin{cases}1,&\ell\leq
2|\tilde{A}|\\\
(-1)^{\lfloor\frac{\ell}{2}\rfloor},&\mathrm{otherwise}\end{cases}\;,$ (113)
$\displaystyle\left(V_{j}\right)_{\ell\ell}=\begin{cases}-1,&\ell=2r_{j}-1,2r_{j}\\\
1,&\mathrm{otherwise}\end{cases}\qquad\textrm{for $j=1,2,4,5$, where
$r_{1}=1,\ r_{2}=|\tilde{A}|,\ r_{4}=|\tilde{A}|+1,\
r_{5}=|\tilde{A}\tilde{C}|$}\;.$ (114)

To write the final formula for the second Rényi entropy let us introduce the
signs $s_{j}^{\gamma}$, for $j=0,1,\ldots,6$, which are positive or negative
depending on whether $P_{j}$ commutes or anticommutes with
$\mathcal{S}^{\gamma}$ respectively, i.e.

$s_{j}^{\gamma}\equiv\begin{cases}1,&[P_{j},\mathcal{S}^{\gamma}]=0\\\
-1,&\\{P_{j},\mathcal{S}^{\gamma}\\}=0\end{cases}$ (115)

We remind that the string operators $\mathcal{S^{\gamma}}$ are chosen so that
$\braket{\mathcal{S}^{\gamma}}\neq 0$ (see appendix E for the exact choice).
Note that if $\mathcal{S}^{\gamma}=S^{\gamma}$ all the signs $s_{j}^{\gamma}$
are positive. Let us also introduce the compact notation

$s^{\gamma}_{j_{1},j_{2},\ldots,j_{6}}=(s_{1}^{\gamma})^{j_{1}}(s_{2}^{\gamma})^{j_{2}}\ldots(s_{6}^{\gamma})^{j_{6}}\;,$
(116)

and, similarly to the single block case,

$\mathcal{V}[j_{0},j_{1},\ldots,j_{6}](\Gamma)=\left(V_{0}^{j_{0}}V_{1}^{j_{1}}\ldots
V_{6}^{j_{6}}\right)^{\mathrm{T}}\Gamma\left(V_{0}^{j_{0}}V_{1}^{j_{1}}\ldots
V_{6}^{j_{6}}\right)\;.$ (117)

The formula for the second Rényi entropy then reads

$\begin{split}&S_{2}(AC)=5\log
2-\log\bigg{\\{}\sum_{j_{1,2,\ldots,6}=0}^{1}\bigg{(}\\\
&\mathcal{G}\Big{(}\Gamma_{0},\mathcal{V}[0,j_{1},\ldots,j_{6}](\Gamma_{0})\Big{)}+\mathcal{G}\Big{(}\Gamma_{0},\mathcal{V}[1,j_{1},\ldots,j_{6}](\Gamma_{0})\Big{)}\\\
&+(-1)^{j_{4}+j_{5}}s^{x}_{j_{1},j_{2},\ldots,j_{6}}|\braket{\mathcal{S}^{x}}|^{2}\left[\mathcal{G}\Big{(}\Gamma_{x},\mathcal{V}[0,j_{1},\ldots,j_{6}](\Gamma_{x}^{\dagger})\Big{)}+s_{0}^{x}(-1)^{|\tilde{B}|}\mathcal{G}\Big{(}\Gamma_{x},\mathcal{V}[1,j_{1},\ldots,j_{6}](\Gamma_{x}^{\dagger})\Big{)}\right]\\\
&+(-1)^{j_{4}+j_{5}}s^{y}_{j_{1},j_{2},\ldots,j_{6}}|\braket{\mathcal{S}^{y}}|^{2}\left[\mathcal{G}\Big{(}\Gamma_{y},\mathcal{V}[0,j_{1},\ldots,j_{6}](\Gamma_{y}^{\dagger})\Big{)}-s_{0}^{y}(-1)^{|\tilde{B}|}\mathcal{G}\Big{(}\Gamma_{y},\mathcal{V}[1,j_{1},\ldots,j_{6}](\Gamma_{y}^{\dagger})\Big{)}\right]\\\
&+s^{z}_{j_{1},j_{2},\ldots,j_{6}}|\braket{\mathcal{S}^{z}}|^{2}\left[\mathcal{G}\Big{(}\Gamma_{z},\mathcal{V}[0,j_{1},\ldots,j_{6}](\Gamma_{z}^{\dagger})\Big{)}-s_{0}^{z}\mathcal{G}\Big{(}\Gamma_{z},\mathcal{V}[1,j_{1},\ldots,j_{6}](\Gamma_{z}^{\dagger})\Big{)}\right]\;\bigg{)}\bigg{\\}},\end{split}$
(118)

Note that the summand $5\log 2$ comes from second term in (47), the factor
$1/2^{6}$ in (54), the factor $1/2$ in (91) and (95) and by taking into
account that in applying the formula (105) a factor $2^{|\tilde{B}|}$ has to
be included to account for the degeneracy in the enlarged space we work with.

The final expression (118) is rather lengthy. Evaluating the second Rényi
entropy exactly requires the evaluation of 512 traces of products of
Gaussians. In practical implementations we will notice that some terms become
negligible with respect to the others, beyond the machine precision, very fast
as we increase the size of the subsystem. To speed up our numerical
implementation, once some terms become negligible as we increase the subsystem
size, we simply drop them for larger subsystems.

##  Appendix E Choice of the strings for the reduced density matrix of
disjoint blocks

In appendices C and D we have not specified the choice of the string operators
$\mathcal{S}^{\gamma}$, $\gamma=x,y,z$, introduced in (96) and appearing in
the final formula (118) for the second Rényi entropy. As discussed after eq.
(96), the string operators are given by

$\mathcal{S}^{\gamma}=F_{\tilde{A}}^{\prime}F_{\tilde{C}}^{\prime}S^{\gamma}$
(119)

where $S^{\gamma}=\prod_{\ell\in\tilde{B}}\tau^{\gamma}_{\ell}$ and
$F_{\tilde{A}}^{\prime},F_{\tilde{C}}^{\prime}$ are products of fermions whose
indices are restricted to block $\tilde{A},\tilde{C}$ respectively. The
operators $F_{\tilde{A}}^{\prime},F_{\tilde{C}}^{\prime}$ have to be chosen in
such a way that $\braket{\mathcal{S}^{\gamma}}\neq 0$. We do this task by
exploiting the properties of the Majorana correlators given in appendix A. In
the stationary state of the flipped-spin quench protocol the string operators

$\displaystyle\mathcal{S}^{\gamma}=\left(\prod_{\ell=|\tilde{A}|-|\tilde{B}|}^{|\tilde{A}|-1}\sigma_{\ell}^{\gamma}\right)S^{\gamma},\qquad\gamma=x,y$
(120) $\displaystyle\mathcal{S}^{z}=\begin{cases}S^{z},&|\tilde{B}|\textrm{
even}\\\ \tau^{z}_{|\tilde{A}|}S^{z},&|\tilde{B}|\textrm{ odd}\end{cases}$
(121)

do the job when $|\tilde{A}|>|\tilde{B}|$. For $|\tilde{A}|\leq|\tilde{B}|$ in
the studied examples (the uppermost curve in figure 5) the same
$\mathcal{S}^{z}$ works, while we have been able to show directly that all
correlation functions of the form
$\braket{F_{\tilde{A}}F_{\tilde{C}}S^{\gamma}}$ vanish for $\gamma=x,y$.
Accordingly, the terms in (118) corresponding to $\gamma=x,y$ were dropped. In
producing the data for the all-spin-up quench protocol in figures 2 and 4 we
have used the same $\mathcal{S}^{z}$ as for the flipped-spin protocol, while
the terms for $\gamma=x,y$ do not contribute and have been dropped (apart from
the exception $|A|=|B|=|C|=3$ for which we have defined the strings properly,
but this particular point is not important anyway).

##  Appendix F Correlation matrix with strings

Here we derive the expression (101) for the correlation matrix, that is
suitable for numerical implementations. We suppose that the string is given by
$\mathcal{S}^{\gamma}=a_{i_{1}}a_{i_{2}}\ldots a_{i_{m}}$ for some positive
even integer $m$. We use Wick theorem (see e.g. appendix G) to write the
elements in terms of pfaffians,

$\left(\Gamma_{\gamma}\right)_{j,\ell}=\delta_{j,\ell}-\mathrm{pf}\begin{pmatrix}M&Q\\\
-Q^{\mathrm{T}}&N\end{pmatrix}/\mathrm{pf}(N)\;,$ (122)

where $M$ is a $2\times 2$ antisymmetric matrix

$M=\begin{pmatrix}0&\braket{a_{f(j)}a_{f(\ell)}}\\\
-\braket{a_{f(j)}a_{f(\ell)}}&0\end{pmatrix},$ (123)

$N$ is defined in (103) and $Q$ is a $2\times m$ matrix

$Q=\begin{pmatrix}\braket{a_{f(j)}a_{i_{1}}}&\braket{a_{f(j)}a_{i_{2}}}&\ldots&\braket{a_{f(j)}a_{i_{m}}}\\\
\braket{a_{f(\ell)}a_{i_{1}}}&\braket{a_{f(\ell)}a_{i_{2}}}&\ldots&\braket{a_{f(\ell)}a_{i_{m}}}\end{pmatrix}\;.$
(124)

Using the pfaffian identity (see e.g. [114])

$\mathrm{pf}\begin{pmatrix}M&Q\\\
-Q^{\mathrm{T}}&N\end{pmatrix}=\mathrm{pf}(N)\mathrm{pf}\left(M+QN^{-1}Q^{T}\right)$
(125)

and the property
$\mathrm{pf}\left(M+QN^{-1}Q^{T}\right)=\left(M+QN^{-1}Q^{T}\right)_{1,2}$,
that is just the definition of the pfaffian for a $2\times 2$ matrix, we get
(101).

##  Appendix G Wick theorem for fermionic Gaussians constructed with complex
antisymmetric matrices

Here we state and prove the Wick theorem, for Gaussians constructed with
complex antisymmetric matrices, directly in the formalism with Majorana
fermions. The proof takes ingredients from ref. [115, 116].

We consider a $2^{d}$ dimensional Fock space associated with Majorana fermions
$a_{1},a_{2},\ldots,a_{2d}$, that satisfy the algebra
$\\{a_{\ell},a_{n}\\}=2\delta_{\ell,n}$. We consider the Gaussian state (84)
and we are interested in the “expectation values”
$\mathrm{tr}[\rho(\Gamma)a_{i_{1}}a_{i_{2}}\ldots a_{i_{n}}]$, where
$i_{1},i_{2},\ldots,i_{n}\in\\{1,2,\ldots,2d\\}$ are not necessarily distinct
indices. For odd $n$ the expectation value is zero. For even $n$, the Wick
theorem, given in the following, expresses the expectation value as a sum over
all contractions, conveniently expressed as a pfaffian. In the proof we will
use the relations (86) (see [68]) and (87) (see [68, 113]), proven
independently of the Wick theorem given here.

##### Wick theorem.

For a Gaussian state
$\rho=\exp(\frac{1}{4}\vec{a}^{\dagger}W\vec{a})/\mathcal{Z}$, with $W$ a
complex antisymmetric matrix (not necessarily Hermitian) such that
$\mathcal{Z}=\mathrm{tr}\exp(\frac{1}{4}\vec{a}^{\dagger}W\vec{a})\neq 0$, and
Majorana fermions $a_{i_{1}},a_{i_{2}},\ldots,a_{i_{n}}$ (not necessarily
distinct), with $n$ an even number, we have

$\braket{a_{i_{1}}a_{i_{2}}\ldots
a_{i_{n}}}=\mathrm{pf}\begin{pmatrix}0&\braket{a_{i_{1}}a_{i_{2}}}&\braket{a_{i_{1}}a_{i_{3}}}&\ldots&\braket{a_{i_{1}}a_{i_{n}}}\\\
-\braket{a_{i_{1}}a_{i_{2}}}&0&\braket{a_{i_{2}}a_{i_{3}}}&\ldots&\braket{a_{i_{2}}a_{i_{n}}}\\\
-\braket{a_{i_{1}}a_{i_{3}}}&-\braket{a_{i_{2}}a_{i_{3}}}&0&\ldots&\braket{a_{i_{3}}a_{i_{n}}}\\\
\vdots&\vdots&\vdots&\ddots&\vdots\\\
-\braket{a_{i_{1}}a_{i_{n}}}&-\braket{a_{i_{2}}a_{i_{n}}}&-\braket{a_{i_{3}}a_{i_{n}}}&\ldots&0\par\end{pmatrix}\;,$
(126)

where $\braket{O}$ stands for $\mathrm{tr}(O\rho)$.

##### Proof.

For some operators $O_{1},O_{2},\ldots,O_{n}$ whose anticommutators are
$c$-numbers we have

$\begin{split}&\mathrm{tr}[O_{1}O_{2}\ldots
O_{n}\rho]=\sum_{k=2}^{n}(-1)^{k}\\{O_{1},O_{k}\\}\mathrm{tr}[(O_{2}O_{3}\ldots
O_{k-1})(O_{k+1}O_{k+2}\ldots O_{n})\rho]\\\ &-\mathrm{tr}[O_{2}O_{3}\ldots
O_{n}O_{1}\rho],\end{split}$ (127)

that is obtained by commuting $O_{1}$ with $O_{k}$ for $k=2,3,\ldots,n$ and
using the property that each anticommutator is a $c$-number to bring it
outside the trace. Here we abuse the notation slightly, as we replace the
operators proportional to the identity by numbers. Relation (127) holds, in
particular, for $O_{k}=a_{i_{k}}$, $k=2,3,\ldots,n$ and $O_{1}=a_{i}$ for any
$i\in\\{1,2,\ldots,2d\\}$. Now the idea is to commute $O_{1}$ with $\rho$ in
the last term in (127), use the cyclic property of the trace and move the term
to the left hand side of the equation.

We use the standard nested commutators identity for complex square matrices
(see e.g. Proposition 3.35 in [117])

$e^{A}Be^{-A}=B+[A,B]+\frac{1}{2}[A,[A,B]]+\ldots+\frac{1}{n!}[A,[A,\ldots[A,B]\ldots]]+\ldots,$
(128)

where in the $n$-th term $n$ commutators appear. It is easy to show

$\left[\frac{\vec{a}^{\dagger}W\vec{a}}{4},a_{i}\right]=-\sum_{j=1}^{2d}W_{ij}a_{j},$
(129)

from which we obtain recursively for $n$ nested commutators

$\left[\frac{\vec{a}^{\dagger}W\vec{a}}{4},\left[\frac{\vec{a}^{\dagger}W\vec{a}}{4},\ldots\left[\frac{\vec{a}^{\dagger}W\vec{a}}{4},a_{i}\right]\ldots\right]\right]=(-1)^{n}\sum_{j=1}^{d}\left(W^{n}\right)_{ij}a_{j}.$
(130)

From identity (128) we now get

$a_{i}\rho=\rho\sum_{j=1}^{2d}\left(e^{W}\right)_{ij}a_{j}.$ (131)

Using (131) and the cyclic property of the trace, from (127) we get

$\sum_{j=1}^{2d}\left(1+e^{W}\right)_{ij}\mathrm{tr}[a_{j}a_{i_{2}}a_{i_{3}}\ldots
a_{i_{n}}\rho]=\sum_{k=2}^{n}(-1)^{k}\\{a_{i},a_{i_{k}}\\}\mathrm{tr}[a_{i_{2}}a_{i_{3}}\ldots
a_{i_{k-1}}a_{i_{k+1}}\ldots a_{i_{n}}\rho]$ (132)

for $j=1,2,\ldots,2d$. From the assumption of the theorem we have
$\mathcal{Z}\neq 0$, and since
$\mathcal{Z}^{2}=\det(e^{\frac{W}{2}}+e^{-\frac{W}{2}})$ it follows that
$1+e^{W}$ is invertible. We can thus multiply (132) by $(1+e^{W})^{-1}_{\ell
i}$ and sum over $i$ to get

$\mathrm{tr}[a_{\ell}a_{i_{2}}\ldots
a_{i_{n}}\rho]=\sum_{k=2}^{n}(-1)^{k}C_{\ell
i_{k}}\mathrm{tr}[a_{i_{2}}a_{i_{3}}\ldots a_{i_{k-1}}a_{i_{k+1}}\ldots
a_{i_{n}}\rho],$ (133)

where $C_{\ell j}=\sum_{i=1}^{2d}(1+e^{W})^{-1}_{\ell i}\\{a_{i},a_{j}\\}$.
Using Majorana anticommutation relations $\\{a_{i},a_{j}\\}=2\delta_{ij}$
explicitly we get a simplification $C_{\ell j}=2(1+e^{W})^{-1}_{\ell j}$.
Since the correlation matrix is given by $\Gamma=\tanh(W/2)$ it follows
$C_{\ell j}=\braket{a_{\ell}a_{j}}$. We have thus shown

$\mathrm{tr}[a_{i_{1}}a_{i_{2}}\ldots
a_{i_{n}}\rho]=\sum_{k=2}^{n}(-1)^{k}\braket{a_{i_{1}}a_{i_{k}}}\mathrm{tr}[a_{i_{2}}a_{i_{3}}\ldots
a_{i_{k-1}}a_{i_{k+1}}\ldots a_{i_{n}}\rho]\;.$ (134)

The theorem follows by applying the recursive definition of the pfaffian for
an $n\times n$ antisymmetric matrix $M$ with $n$ even,

$\mathrm{pf}(M)=\sum_{k=2}^{n}(-1)^{k}M_{1k}\mathrm{pf}(M_{\hat{1},\hat{k}}),$
(135)

where $M_{\hat{1},\hat{k}}$ stands for the matrix with both the first and the
$k$-th rows and columns removed.

## References

  * Marić and Fagotti [2023] V. Marić and M. Fagotti, Universality in the tripartite information after global quenches, Phys. Rev. B 108, L161116 (2023).
  * Marić and Fagotti [2023] V. Marić and M. Fagotti, Universality in the tripartite information after global quenches: (generalised) quantum xy models, Journal of High Energy Physics 2023, 140 (2023).
  * 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).
  * Eisert  _et al._ [2015] J. Eisert, M. Friesdorf, and C. Gogolin, Quantum many-body systems out of equilibrium, Nature Phys 11, 124 (2015).
  * Gogolin and Eisert [2016] C. Gogolin and J. Eisert, Equilibration, thermalisation, and the emergence of statistical mechanics in closed quantum systems, Rep. Prog. Phys. 79, 056001 (2016).
  * Rigol  _et al._ [2008] M. Rigol, V. Dunjko, and M. Olshanii, Thermalization and its mechanism for generic isolated quantum systems, Nature 452, 854 (2008).
  * Deutsch [1991] J. M. Deutsch, Quantum statistical mechanics in a closed system, Phys. Rev. A 43, 2046 (1991).
  * Srednicki [1994] M. Srednicki, Chaos and quantum thermalization, Phys. Rev. E 50, 888 (1994).
  * Deutsch [2018] J. M. Deutsch, Eigenstate thermalization hypothesis, Reports on Progress in Physics 81, 082001 (2018).
  * Kinoshita  _et al._ [2006] T. Kinoshita, T. Wenger, and D. S. Weiss, A quantum newton’s cradle, Nature 440, 900 (2006).
  * Rigol  _et al._ [2007] M. Rigol, V. Dunjko, V. Yurovsky, and M. Olshanii, Relaxation in a Completely Integrable Many-Body Quantum System: An Ab Initio Study of the Dynamics of the Highly Excited States of 1D Lattice Hard-Core Bosons, Phys. Rev. Lett. 98, 050405 (2007).
  * Doyon [2017] B. Doyon, Thermalization and Pseudolocality in Extended Quantum Systems, Commun. Math. Phys. 351, 155 (2017).
  * Vidmar and Rigol [2016] L. Vidmar and M. Rigol, Generalized Gibbs ensemble in integrable lattice models, J. Stat. Mech. 2016, 064007 (2016).
  * Essler and Fagotti [2016] F. H. L. Essler and M. Fagotti, Quench dynamics and relaxation in isolated integrable quantum spin chains, J. Stat. Mech. 2016, 064002 (2016).
  * Alba  _et al._ [2021] V. Alba, B. Bertini, M. Fagotti, L. Piroli, and P. Ruggiero, Generalized-hydrodynamic approach to inhomogeneous quenches: correlations, entanglement and quantum effects, Journal of Statistical Mechanics: Theory and Experiment 2021, 114004 (2021).
  * Bertini  _et al._ [2016] B. Bertini, M. Collura, J. De Nardis, and M. Fagotti, Transport in Out-of-Equilibrium ${XXZ}$ Chains: Exact Profiles of Charges and Currents, Phys. Rev. Lett. 117, 207201 (2016).
  * Castro-Alvaredo  _et al._ [2016] O. A. Castro-Alvaredo, B. Doyon, and T. Yoshimura, Emergent Hydrodynamics in Integrable Quantum Systems Out of Equilibrium, Phys. Rev. X 6, 041065 (2016).
  * Bertini and Fagotti [2016] B. Bertini and M. Fagotti, Determination of the Nonequilibrium Steady State Emerging from a Defect, Phys. Rev. Lett. 117, 130402 (2016).
  * Holzhey  _et al._ [1994] C. Holzhey, F. Larsen, and F. Wilczek, Geometric and renormalized entropy in conformal field theory, Nuclear Physics B 424, 443 (1994).
  * Calabrese and Cardy [2004] P. Calabrese and J. Cardy, Entanglement entropy and quantum field theory, Journal of Statistical Mechanics: Theory and Experiment 2004, P06002 (2004).
  * Korepin [2004] V. E. Korepin, Universality of entropy scaling in one dimensional gapless models, Phys. Rev. Lett. 92, 096402 (2004).
  * Vidal  _et al._ [2003] G. Vidal, J. I. Latorre, E. Rico, and A. Kitaev, Entanglement in Quantum Critical Phenomena, Phys. Rev. Lett. 90, 227902 (2003).
  * Calabrese and Cardy [2009] P. Calabrese and J. Cardy, Entanglement entropy and conformal field theory, Journal of Physics A: Mathematical and Theoretical 42, 504005 (2009).
  * Calabrese and Cardy [2005] P. Calabrese and J. Cardy, Evolution of entanglement entropy in one-dimensional systems, J. Stat. Mech. 2005, P04010 (2005).
  * Sotiriadis and Cardy [2008] S. Sotiriadis and J. Cardy, Inhomogeneous quantum quenches, Journal of Statistical Mechanics: Theory and Experiment 2008, P11003 (2008).
  * Bastianello and Calabrese [2018] A. Bastianello and P. Calabrese, Spreading of entanglement and correlations after a quench with intertwined quasiparticles, SciPost Phys. 5, 033 (2018).
  * Bertini  _et al._ [2022] B. Bertini, K. Klobas, V. Alba, G. Lagnese, and P. Calabrese, Growth of rényi entropies in interacting integrable models and the breakdown of the quasiparticle picture, Phys. Rev. X 12, 031016 (2022).
  * Alba and Calabrese [2017] V. Alba and P. Calabrese, Entanglement and thermodynamics after a quantum quench in integrable systems, Proceedings of the National Academy of Sciences 114, 7947 (2017), https://www.pnas.org/doi/pdf/10.1073/pnas.1703516114 .
  * Alba and Calabrese [2018] V. Alba and P. Calabrese, Entanglement dynamics after quantum quenches in generic integrable systems, SciPost Phys. 4, 017 (2018).
  * Casini  _et al._ [2016] H. Casini, H. Liu, and M. Mezei, Spread of entanglement and causality, Journal of High Energy Physics 2016, 77 (2016).
  * Liu and Suh [2014] H. Liu and S. J. Suh, Entanglement growth during thermalization in holographic systems, Phys. Rev. D 89, 066012 (2014).
  * Skinner  _et al._ [2019] B. Skinner, J. Ruhman, and A. Nahum, Measurement-induced phase transitions in the dynamics of entanglement, Physical Review X 9, 10.1103/physrevx.9.031009 (2019).
  * Eisler and Zimborás [2014] V. Eisler and Z. Zimborás, Area-law violation for the mutual information in a nonequilibrium steady state, Phys. Rev. A 89, 032321 (2014).
  * Fraenkel and Goldstein [2021] S. Fraenkel and M. Goldstein, Entanglement measures in a nonequilibrium steady state: Exact results in one dimension, SciPost Phys. 11, 085 (2021).
  * Fagotti  _et al._ [2022] M. Fagotti, V. Marić, and L. Zadnik, Nonequilibrium symmetry-protected topological order: emergence of semilocal gibbs ensembles, arXiv:2205.02221 (2022).
  * Cerf and Adami [1998] N. J. Cerf and C. Adami, Information theory of quantum entanglement and measurement, Physica D: Nonlinear Phenomena 120, 62 (1998), proceedings of the Fourth Workshop on Physics and Consumption.
  * Calabrese  _et al._ [2009] P. Calabrese, J. Cardy, and E. Tonni, Entanglement entropy of two disjoint intervals in conformal field theory, Journal of Statistical Mechanics: Theory and Experiment 2009, P11001 (2009).
  * Gluza  _et al._ [2019] M. Gluza, J. Eisert, and T. Farrelly, Equilibration towards generalized Gibbs ensembles in non-interacting theories, SciPost Phys. 7, 38 (2019).
  * Zauner  _et al._ [2015] V. Zauner, M. Ganahl, H. G. Evertz, and T. Nishino, Time evolution within a comoving window: scaling of signal fronts and magnetization plateaus after a local quench in quantum spin chains, Journal of Physics: Condensed Matter 27, 425602 (2015).
  * Eisler and Maislinger [2020] V. Eisler and F. Maislinger, Front dynamics in the XY chain after local excitations, SciPost Phys. 8, 37 (2020).
  * Eisler  _et al._ [2016] V. Eisler, F. Maislinger, and H. G. Evertz, Universal front propagation in the quantum Ising chain with domain-wall initial states, SciPost Phys. 1, 014 (2016).
  * Gruber and Eisler [2021] M. Gruber and V. Eisler, Entanglement spreading after local fermionic excitations in the XXZ chain, SciPost Phys. 10, 5 (2021).
  * Eisler and Maislinger [2018] V. Eisler and F. Maislinger, Hydrodynamical phase transition for domain-wall melting in the XY chain, Phys. Rev. B 98, 161117 (2018).
  * Bidzhiev  _et al._ [2022] K. Bidzhiev, M. Fagotti, and L. Zadnik, Macroscopic Effects of Localized Measurements in Jammed States of Quantum Spin Chains, Phys. Rev. Lett. 128, 130603 (2022).
  * Zadnik  _et al._ [2022] L. Zadnik, S. Bocini, K. Bidzhiev, and M. Fagotti, Measurement catastrophe and ballistic spread of charge density with vanishing current, Journal of Physics A: Mathematical and Theoretical 55, 474001 (2022).
  * Fagotti [2022] M. Fagotti, Global Quenches after Localized Perturbations, Phys. Rev. Lett. 128, 110602 (2022).
  * Bocini and Fagotti [2023] S. Bocini and M. Fagotti, Growing schrödinger’s cat states by local unitary time evolution of product states, arXiv:2210.15585 [quant-ph] (2023).
  * Alaeian and Buča [2022] H. Alaeian and B. Buča, Exact multistability and dissipative time crystals in interacting fermionic lattices, Communications Physics 5, 318 (2022).
  * Groisman  _et al._ [2005] B. Groisman, S. Popescu, and A. Winter, Quantum, classical, and total amount of correlations in a quantum state, Phys. Rev. A 72, 032317 (2005).
  * Wolf  _et al._ [2008] M. M. Wolf, F. Verstraete, M. B. Hastings, and J. I. Cirac, Area laws in quantum systems: Mutual information and correlations, Phys. Rev. Lett. 100, 070502 (2008).
  * Kuwahara  _et al._ [2021] T. Kuwahara, A. M. Alhambra, and A. Anshu, Improved thermal area law and quasilinear time algorithm for quantum gibbs states, Phys. Rev. X 11, 011047 (2021).
  * Eisert  _et al._ [2010] J. Eisert, M. Cramer, and M. B. Plenio, Colloquium: Area laws for the entanglement entropy, Rev. Mod. Phys. 82, 277 (2010).
  * Bernigau  _et al._ [2015] H. Bernigau, M. J. Kastoryano, and J. Eisert, Mutual information area laws for thermal free fermions, Journal of Statistical Mechanics: Theory and Experiment 2015, P02008 (2015).
  * Lemm and Siebert [2023] M. Lemm and O. Siebert, Thermal Area Law for Lattice Bosons, Quantum 7, 1083 (2023).
  * Alhambra [2022] A. M. Alhambra, Quantum many-body systems in thermal equilibrium, arXiv:2204.08349 (2022).
  * Casini and Huerta [2009] H. Casini and M. Huerta, Remarks on the entanglement entropy for disconnected regions, Journal of High Energy Physics 2009, 048 (2009).
  * Fagotti [2012] M. Fagotti, New insights into the entanglement of disjoint blocks, EPL (Europhysics Letters) 97, 17007 (2012).
  * Kitaev and Preskill [2006] A. Kitaev and J. Preskill, Topological Entanglement Entropy, Phys. Rev. Lett. 96, 110404 (2006).
  * Caraglio and Gliozzi [2008] M. Caraglio and F. Gliozzi, Entanglement entropy and twist fields, Journal of High Energy Physics 2008, 076 (2008).
  * Furukawa  _et al._ [2009] S. Furukawa, V. Pasquier, and J. Shiraishi, Mutual Information and Boson Radius in a $c=1$ Critical System in One Dimension, Phys. Rev. Lett. 102, 170602 (2009).
  * Rajabpour and Gliozzi [2012] M. A. Rajabpour and F. Gliozzi, Entanglement entropy of two disjoint intervals from fusion algebra of twist fields, Journal of Statistical Mechanics: Theory and Experiment 2012, P02016 (2012).
  * Calabrese  _et al._ [2011a] P. Calabrese, J. Cardy, and E. Tonni, Entanglement entropy of two disjoint intervals in conformal field theory: II, Journal of Statistical Mechanics: Theory and Experiment 2011, P01021 (2011a).
  * Coser  _et al._ [2014] A. Coser, L. Tagliacozzo, and E. Tonni, On rényi entropies of disjoint intervals in conformal field theory, Journal of Statistical Mechanics: Theory and Experiment 2014, P01008 (2014).
  * Ruggiero  _et al._ [2018] P. Ruggiero, E. Tonni, and P. Calabrese, Entanglement entropy of two disjoint intervals and the recursion formula for conformal blocks, Journal of Statistical Mechanics: Theory and Experiment 2018, 113101 (2018).
  * Alba  _et al._ [2010] V. Alba, L. Tagliacozzo, and P. Calabrese, Entanglement entropy of two disjoint blocks in critical ising models, Phys. Rev. B 81, 060411 (2010).
  * Blanco and Casini [2011] D. D. Blanco and H. Casini, Entanglement entropy for non-coplanar regions in quantum field theory, Classical and Quantum Gravity 28, 215015 (2011).
  * Alba  _et al._ [2011] V. Alba, L. Tagliacozzo, and P. Calabrese, Entanglement entropy of two disjoint intervals in c=1 theories, Journal of Statistical Mechanics: Theory and Experiment 2011, P06012 (2011).
  * Fagotti and Calabrese [2010] M. Fagotti and P. Calabrese, Entanglement entropy of two disjoint blocks in XY chains, Journal of Statistical Mechanics: Theory and Experiment 2010, P04016 (2010).
  * Fries and Reyes [2019] P. Fries and I. A. Reyes, Entanglement and relative entropy of a chiral fermion on the torus, Phys. Rev. D 100, 105015 (2019).
  * Balasubramanian  _et al._ [2011] V. Balasubramanian, A. Bernamonti, N. Copland, B. Craps, and F. Galli, Thermalization of mutual and tripartite information in strongly coupled two dimensional conformal field theories, Phys. Rev. D 84, 105017 (2011).
  * Grava  _et al._ [2021] T. Grava, A. P. Kels, and E. Tonni, Entanglement of two disjoint intervals in conformal field theory and the 2d coulomb gas on a lattice, Phys. Rev. Lett. 127, 141605 (2021).
  * Ares  _et al._ [2021] F. Ares, R. Santachiara, and J. Viti, Crossing-symmetric twist field correlators and entanglement negativity in minimal CFTs, Journal of High Energy Physics 2021, 175 (2021).
  * Agón and Faulkner [2016] C. A. Agón and T. Faulkner, Quantum corrections to holographic mutual information, Journal of High Energy Physics 2016, 118 (2016).
  * Agón  _et al._ [2022] C. A. Agón, P. Bueno, and H. Casini, Tripartite information at long distances, SciPost Phys. 12, 153 (2022).
  * Ali-Akbari  _et al._ [2021] M. Ali-Akbari, M. Rahimi, and M. Asadi, Holographic mutual and tripartite information in a non-conformal background, Nuclear Physics B 964, 115329 (2021).
  * Hayden  _et al._ [2013] P. Hayden, M. Headrick, and A. Maloney, Holographic mutual information is monogamous, Phys. Rev. D 87, 046003 (2013).
  * Carollo and Alba [2022] F. Carollo and V. Alba, Entangled multiplets and spreading of quantum correlations in a continuously monitored tight-binding chain, Phys. Rev. B 106, L220304 (2022).
  * Parez  _et al._ [2023] G. Parez, P.-A. Bernard, N. Crampé, and L. Vinet, Multipartite information of free fermions on hamming graphs, Nuclear Physics B 990, 116157 (2023).
  * Hosur  _et al._ [2016] P. Hosur, X.-L. Qi, D. A. Roberts, and B. Yoshida, Chaos in quantum channels, Journal of High Energy Physics 2016, 4 (2016).
  * Schnaack  _et al._ [2019] O. Schnaack, N. Bölter, S. Paeckel, S. R. Manmana, S. Kehrein, and M. Schmitt, Tripartite information, scrambling, and the role of Hilbert space partitioning in quantum lattice models, Phys. Rev. B 100, 224302 (2019).
  * Sünderhauf  _et al._ [2019] C. Sünderhauf, L. Piroli, X.-L. Qi, N. Schuch, and J. I. Cirac, Quantum chaos in the Brownian SYK model with large finite N : OTOCs and tripartite information, Journal of High Energy Physics 2019, 38 (2019).
  * Kuno  _et al._ [2022] Y. Kuno, T. Orito, and I. Ichinose, Purification and scrambling in a chaotic hamiltonian dynamics with measurements, Phys. Rev. B 106, 214304 (2022).
  * Parez and Bonsignori [2023] G. Parez and R. Bonsignori, Analytical results for the entanglement dynamics of disjoint blocks in the xy spin chain, Journal of Physics A: Mathematical and Theoretical 55, 505005 (2023).
  * Caceffo and Alba [2023] F. Caceffo and V. Alba, Negative tripartite mutual information after quantum quenches in integrable systems, Phys. Rev. B 108, 134434 (2023).
  * Bluhm  _et al._ [2022] A. Bluhm, Á. Capel, and A. Pérez-Hernández, Exponential decay of mutual information for Gibbs states of local Hamiltonians, Quantum 6, 650 (2022).
  * Hastings and Koma [2006] M. B. Hastings and T. Koma, Spectral gap and exponential decay of correlations, Communications in Mathematical Physics 265, 781 (2006).
  * Zadnik and Fagotti [2021] L. Zadnik and M. Fagotti, The Folded Spin-1/2 XXZ Model: I. Diagonalisation, Jamming, and Ground State Properties, SciPost Phys. Core 4, 10 (2021).
  * Zadnik and Garrahan [2023] L. Zadnik and J. P. Garrahan, Slow heterogeneous relaxation due to constraints in dual xxz models, Phys. Rev. B 108, L100304 (2023).
  * Pai and Pretko [2020] S. Pai and M. Pretko, Fractons from confinement in one dimension, Phys. Rev. Res. 2, 013094 (2020).
  * Eck and Fendley [2023] L. Eck and P. Fendley, From the xxz chain to the integrable rydberg-blockade ladder via non-invertible duality defects, arXiv:2302.14081 [cond-mat.stat-mech] (2023).
  * Lieb  _et al._ [1961] E. Lieb, T. Schultz, and D. Mattis, Two soluble models of an antiferromagnetic chain, Annals of Physics 16, 407 (1961).
  * Calabrese  _et al._ [2011b] P. Calabrese, F. H. L. Essler, and M. Fagotti, Quantum Quench in the Transverse-Field Ising Chain, Phys. Rev. Lett. 106, 227203 (2011b).
  * Calabrese  _et al._ [2012a] P. Calabrese, F. H. L. Essler, and M. Fagotti, Quantum quench in the transverse field ising chain: I. time evolution of order parameter correlators, Journal of Statistical Mechanics: Theory and Experiment 2012, P07016 (2012a).
  * Calabrese  _et al._ [2012b] P. Calabrese, F. H. L. Essler, and M. Fagotti, Quantum quenches in the transverse field ising chain: Ii. stationary state properties, Journal of Statistical Mechanics: Theory and Experiment 2012, P07022 (2012b).
  * Fagotti and Essler [2013] M. Fagotti and F. H. L. Essler, Reduced density matrix after a quantum quench, Phys. Rev. B 87, 245107 (2013).
  * Fagotti [2020] M. Fagotti, Locally quasi-stationary states in noninteracting spin chains, SciPost Phys. 8, 48 (2020).
  * Kennedy and Tasaki [1992] T. Kennedy and H. Tasaki, Hidden ${{\mathbb{Z}}_{2}\times{\mathbb{Z}}_{2}}$ symmetry breaking in Haldane-gap antiferromagnets, Phys. Rev. B 45, 304 (1992).
  * Else  _et al._ [2013] D. V. Else, S. D. Bartlett, and A. C. Doherty, Hidden symmetry-breaking picture of symmetry-protected topological order, Phys. Rev. B 88, 085114 (2013).
  * Duivenvoorden and Quella [2013] K. Duivenvoorden and T. Quella, From symmetry-protected topological order to Landau order, Phys. Rev. B 88, 125115 (2013).
  * Iglói and Peschel [2010] F. Iglói and I. Peschel, On reduced density matrices for disjoint subsystems, Europhysics Letters 89, 40001 (2010).
  * Pollmann and Turner [2012] F. Pollmann and A. M. Turner, Detection of symmetry-protected topological phases in one dimension, Phys. Rev. B 86, 125441 (2012).
  * Coser  _et al._ [2016] A. Coser, E. Tonni, and P. Calabrese, Spin structures and entanglement of two disjoint intervals in conformal field theories, Journal of Statistical Mechanics: Theory and Experiment 2016, 053109 (2016).
  * Fagotti [2014] M. Fagotti, On conservation laws, relaxation and pre-relaxation after a quantum quench, J. Stat. Mech. 2014, P03016 (2014).
  * Ilievski  _et al._ [2016] E. Ilievski, M. Medenjak, T. Prosen, and L. Zadnik, Quasilocal charges in integrable lattice systems, J. Stat. Mech. 2016, 064008 (2016).
  * Dalmonte  _et al._ [2022] M. Dalmonte, V. Eisler, M. Falconi, and B. Vermersch, Entanglement hamiltonians: From field theory to lattice models and experiments, Annalen der Physik 534, 2200064 (2022), https://onlinelibrary.wiley.com/doi/pdf/10.1002/andp.202200064 .
  * Fagotti [2016] M. Fagotti, Charges and currents in quantum spin chains: late-time dynamics and spontaneous currents, J. Phys. A: Math. Theor. 50, 034005 (2016).
  * Vidal and Werner [2002] G. Vidal and R. F. Werner, Computable measure of entanglement, Phys. Rev. A 65, 032314 (2002).
  * Plenio [2005] M. B. Plenio, Logarithmic negativity: A full entanglement monotone that is not convex, Phys. Rev. Lett. 95, 090503 (2005).
  * Calabrese  _et al._ [2012c] P. Calabrese, J. Cardy, and E. Tonni, Entanglement negativity in quantum field theory, Phys. Rev. Lett. 109, 130502 (2012c).
  * 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).