Multiple relaxation times in perturbed XXZ chain

M. Mierzejewski1, J. Pawłowski1, P. Prelovšek2,3, and J. Herbrych1$\star$

1 Department of Theoretical Physics, Faculty of Fundamental Problems of
Technology, Wrocław University of Science and Technology, 50-370 Wrocław,
Poland

2 Jožef Stefan Institute, SI-1000 Ljubljana, Slovenia

3 Faculty of Mathematics and Physics, University of Ljubljana, SI-1000
Ljubljana, Slovenia

⋆<EMAIL_ADDRESS>
## Abstract

We numerically study the relaxation of correlation functions in weakly
perturbed integrable XXZ chain. The decay of the spin-current and the energy-
current correlations at zero magnetization are well described by single, but
quite distinct, relaxation rates governed by the square of the perturbation
strength $g$. However, at finite magnetization a single correlation function
reveals multiple relaxation rates. The result can be understood in terms of
multi-scale relaxation scenario, where various relaxation times are linked
with various quantities which are conserved in the reference integrable
system. On the other hand, the correlations of non-commuting quantities, being
conserved at particular anisotropies $\Delta$, decay non-exponentially with
characteristic time scale linear in $g$.

###### Contents

  1. 1 Introduction
  2. 2 Relaxation of spin and energy currents
  3. 3 Memory-function analysis of the energy-current correlations
  4. 4 Dependence on the form of perturbation
  5. 5 Non-commuting conserved quantities and their non-exponential relaxation
  6. 6 Conclusions

##  1 Introduction

Integrable quantum many-body systems attract a lot of attention due to their
unique properties, as well as due to development of new analytical and
numerical tools to deal with them (for a recent review see Ref. [1]). The
crucial role in the behavior of such systems is played by the presence of
extensive number of local and/or quasilocal conserved quantities (CQ), which
have important consequences for the (lack of) relaxation of observables and
for the transport properties. The latter consequences are formally expressed
via the Mazur bounds which relate the long-time correlations (stiffness) of
observables with their projections on the local CQ (charges) [2, 3]. However,
the microscopic models are integrable only for fine-tuned sets of parameters,
while more realistic systems might be described in terms of nearly integrable
(NI) models which contain small but non-vanishing perturbation that breaks the
integrability [4, 5]. An important question concerns the details of the
integrability breaking, in particular whether the asymptotic dynamics becomes
consistent with the generic dissipative diffusive-type transport [6, 7, 8]. So
far, the properties of NI models are better understood at intermediate time-
scales, when the dynamics resembles that of integrable models, the phenomenon
known as prethermalization [9, 10, 11].

The problem of asymptotic dynamics of NI models at long time scales [12, 13]
appears to be more complex. It has been argued that the time evolution can be
accurately described as the generalized Gibbs enseble [14, 15] with the time-
dependent Lagrange parameters [16]. Several analytical and numerical studies
have demonstrated that breaking of integrability leads to exponential
relaxation of typical observables and that the corresponding relaxation rates
scale quadratically with the strength of the perturbation [12, 17, 18, 19,
20]. Description of NI models within the framework of generalized
hydrodynamics (GHD) [21, 22, 23, 24, 25] seems also demanding, since the
generic integrability-breaking processes involve large momenta transfers [26,
27]. Nevertheless, recent results suggest that arbitrarily weak perturbation
applied to a macroscopic integrable system restores the generic chaotic
(dissipative) dynamics [28]. Still, more detailed understanding or even
theoretical analysis of the relaxation of different quantities in NI systems
is mostly lacking so far.

In this work we numerically study relaxation of several operators in a NI XXZ
model, employing the microcanonical Lanczos method (MCLM) [29, 30, 31] which
allows to reach long-enough times even for NI model. In particular, we confirm
that both spin-current and energy-current decay exponentially in time with
relaxation rates that scale quadratically with the strength of perturbation.
However, the energy-current relaxation rate turns out to be much smaller than
the relaxation rate for the spin current. The presence of distinct relaxation
times is consistent with the predictions of GHD [26, 27]. This result can be
simply explained by single (but different for both quantities)
local/quasilocal CQ (charges) involved in the relaxation process [19]. Still
even in this case, a more detailed analysis using the memory-function
indicates a weak contribution of other CQ from the same symmetry sector.
Moreover, we can construct a single observable that clearly reveals multiple
relaxation rates or, in other words, that its relaxation cannot be described
by a single exponential function. We confirm the latter possibility studying
relaxation of current correlations in sectors with small magnetization, i.e.,
with non-zero total spin $S^{z}_{\rm tot}\neq 0$. Results in this case can be
explained with a projection on the decay of several CQ with different
relaxation times. Furthermore, we find that the relaxation explicitly depends
on the form and the symmetry of the perturbation.

Finally, we study the weakly perturbed XXZ model for the specific values
anisotropy parameters, $\Delta$, when the many-body spectra are
macroscopically degenerate. Such degeneracy in the integrable model allows for
additional local CQ which do not commute with $S^{z}_{\rm tot}$ or with other
local CQ [32]. It is the case, e.g., for the anisotropy parameter $\Delta=0.5$
which we use in the present work. However, the simplest example can be studied
for $\Delta=1$, where $S^{x}_{\rm tot}$ is conserved and but does not commute
with $S^{z}_{\rm tot}$. Upon introducing a perturbation, the correlation
functions of quantities which do not commute with $S^{z}_{\rm tot}$ show -
instead of exponential decay - an approximately Gaussian decay with the
characteristic relaxation time that scales linearly with the perturbation
strength.

##  2 Relaxation of spin and energy currents

We consider the one-dimensional XXZ chain with $L$ sites assuming periodic
boundary conditions, with a specific perturbation involving the next-nearest-
neighbor interaction

$H=H_{0}+gH^{\prime},\quad H_{0}=\sum_{i}h_{i},\quad
h_{i}=\frac{J}{2}\left(S^{+}_{i}S^{-}_{i+1}+S^{-}_{i}S^{+}_{i+1}\right)+J\Delta
S^{z}_{i}S^{z}_{i+1},\quad H^{\prime}=J\sum_{i}S^{z}_{i}S^{z}_{i+2}\,,$ (1)

where $S^{\pm,z}$ are spin-$\frac{1}{2}$ operators. The model is integrable
for $g=0$ and the integrability is broken for $g\neq 0$. In this section we
study the relaxation of the spin current, $j_{\sigma}$, as well as of the
energy current, $j_{\kappa}$, obtained for the unperturbed system

$j_{\sigma}=i\sum_{l,l^{\prime}}l[h_{l^{\prime}},S^{z}_{l}]\,,\qquad
j_{\kappa}=i\sum_{l,l^{\prime}}l[h_{l^{\prime}},h_{l}]\,,$ (2)

for which we calculate the corresponding normalized current-current
correlation functions

$C_{\sigma}(t)=\frac{\langle j_{\sigma}(t)j_{\sigma}\rangle}{\langle
j_{\sigma}j_{\sigma}\rangle}\,,\quad\quad C_{\kappa}(t)=\frac{\langle
j_{\kappa}(t)j_{\kappa}\rangle}{\langle j_{\kappa}j_{\kappa}\rangle}\,.$ (3)

We focus on the result for large (infinite) temperature $T=1/\beta\gg J$.
Here, $\langle\dots\rangle={\rm Tr}(\dots)$ denotes averaging either over the
canonical ensemble with fixed $S^{z}_{\rm tot}$ (in Secs. 2-4) or
grandcanonical ensemble (in Sec. 5), and
$j_{\alpha}(t)=e^{iHt}j_{\alpha}e^{-iHt}$. Note that $C_{\sigma}(t)$ and
$C_{\kappa}(t)$ are related via the Fourier transform to the dynamical spin
conductivity, $\sigma(\omega)=\beta\langle
j_{\sigma}j_{\sigma}\rangle\tilde{C}_{\sigma}(\omega)$, and the thermal
conductivity, $\kappa(\omega)=\beta^{2}\langle
j_{\kappa}j_{\kappa}\rangle\tilde{C}_{\kappa}(\omega)$, respectively.

For the numerical evaluation of the spectral functions
$\tilde{C}_{\alpha}(\omega)$, on systems up to $L=28$ sites, we use the MCLM
which allows for very high frequency resolution, $\delta\omega\lesssim
10^{-3}J$, and consequently enables evaluation of $C_{\alpha}(t)$ up to
$t\lesssim t_{\rm max}\sim 10^{3}/J$. Reaching large $t_{\rm max}J\gg 1$ is
crucial to follow the slow relaxation at weak perturbations, $g\ll 1$, and to
resolve possible distinct relaxation times. This is achieved within MCLM by
using large number of Lanczos steps $M_{L}\sim 5.10^{4}$. In the following, we
numerically calculate canonical $C_{\alpha}(t)$, i.e., in sectors with fixed
total spin projections $S^{z}_{\rm tot}$. Large but finite $M_{L}$ also sets
the smallest reliable value of $C_{\alpha}(t)\gtrsim 2.10^{-3}$.

Figure 1: Normalized correlation functions $C_{\alpha}(t)$ for (a) spin
current $j_{\sigma}$ and (b) energy current $j_{\kappa}$, respectively, in
$S^{z}_{\rm tot}=4$ sector, calculated for $L=28$ and $\Delta=0.5$. Insets:
the same as in main panels but with the time $t$ rescaled by the square of the
perturbation strength $g^{2}$.

As a first step, we establish the range of the perturbation strength, $g$,
which is relevant for NI systems. On the one hand, $g$ should be small to
remain a perturbation but, on the other hand, $g$ should be sufficiently large
so that the effective mean-free path is smaller then the system size ($L\leq
28$). The latter requirement means that the accessible system sizes impose
lower bounds on the accessible values of $g$. Since we are not able to
directly evaluate the mean-free path, we use the range of parameters for which
the decay of the correlation functions does not reveal any significant $L$
dependence but also exhibits a clear quadratic dependence on $g$ [12, 19, 20].
Main panels in Fig. 1(a,b) show, respectively, $C_{\sigma}(t)$ and
$C_{\kappa}(t)$ within the sector $S^{z}_{\rm tot}=4$ for various $g=0.1-0.4$,
whereas the time in the insets is rescaled by $g^{2}$. For $0.15\leq g\leq
0.3$ we observe a perfect collapse of all curves, whereas the latter collapse
is worse for weaker $g$ indicating that the effective mean-free path becomes
larger than the systems size. Therefore from now on, we set $g=0.15$.

Figure 2: Correlation functions $C_{\alpha}(t)$ in sectors with various
$S^{z}_{\rm tot}$. Dashed curves in (a) depict exponential fits for
$\alpha=j_{\sigma},j_{\kappa}$ at $s=S^{z}_{\rm tot}=0$, while in (c,d) for
$s\neq 0$ only long-time fits are presented. Calculated for $L=28$,
$\Delta=0.5$, and $g=0.15$.

In Fig. 2(a) we directly compare the relaxation of spin and energy currents
$C_{\alpha}(t)$ at zero-magnetization, $S^{z}_{\rm tot}=0$. Here, $H$ is
invariant under the $Z_{2}$ spin-flip transformation, generated by the parity
operator, $P=\prod_{J}(S^{+}_{j}+S^{-}_{j})$. Since the spin/energy current is
odd/even under this transformation, $Pj_{\sigma}P=-j_{\sigma}$ and
$Pj_{\kappa}P=j_{\kappa}$, in this sector both currents are mutually
orthogonal, i.e., $\langle j_{\sigma}j_{\kappa}\rangle=0$. Consequently,
relaxation of both correlations are evidently different. At longer times,
$tJ>100$, the decays do not show any clear deviations from simple exponential,
$C_{\alpha}(t)\propto\exp(-t/\tau^{0}_{\alpha})$, $\alpha=\sigma,\kappa$,
where the upper index in the relaxation times, $\tau^{s}_{\alpha}$, marks the
value of $s=S^{z}_{\rm tot}$. Nevertheless, it is clear that the relaxation of
$j_{\sigma}$ is much faster than that of $j_{\kappa}$,
$\tau^{0}_{\sigma}\simeq 3\tau^{0}_{\kappa}$, confirming that the breaking of
integrability leads to multiple distinct relaxation times, here due to
different symmetry sectors involved. The evident differences are also at short
times $t$. Since $j_{\sigma}$ is not CQ in the reference $H_{0}$, there is an
incoherent drop to $C^{0}_{\sigma}\sim 0.5$ at $tJ\sim O(1)$, reflecting the
nontrivial spin stiffness [3] and finite overlap with quasilocal CQ [33, 34].
On the other hand, $j_{\kappa}$ is CQ at $g=0$, so one might expect a single
exponential decay in the whole range of $t$. To good approximation this is
indeed the case, but there is still some visible deviation at $tJ<100$ about
which we comment in more detail in Sec. 3.

More challenging question is whether distinct relaxation rates can be observed
in the dynamics of a single observable. This might happen when we consider
$S^{z}_{\rm tot}\neq 0$ sectors where the above symmetry arguments do not
apply. In Figs. 2(b-d) we show that a departure from a simple exponential
relaxation becomes increasingly visible for $S^{z}_{\rm tot}\neq 0$, when also
$\langle j_{\sigma}j_{\kappa}\rangle\neq 0$. Figs. 2(c,d) for $S^{z}_{\rm
tot}=4,6$ confirm that for longest $tJ>500$ the relaxation is asymptotically
determined by the same $\tau^{s}_{\kappa}$ for both currents, while
$S^{z}_{\rm tot}=2$ case on Fig. 2(b) is marginal due to fast decay of
$C_{\sigma}(t)$ (and the limitation $C_{\sigma}(t)>2.10^{-3}$). Nevertheless,
both correlation functions, and in particular $C_{\sigma}(t)$, reveal a clear
deviation from a single-exponential decay.

The modest dependence of $\tau^{s}_{\kappa}$ on $s$, as extracted from Fig.
2(a-d), is shown in Fig. 3(a). More importantly, our analysis indicate that
one can fit the decay of correlations $C_{\alpha}(t)$ for $S^{z}_{\rm tot}>0$
by a sum of two exponential functions with distinct relaxation rates
$\tau^{s}_{\kappa},\tau^{s}_{\sigma}$. This is shown in Fig. 3(b,c) for
$C_{\sigma}(t)$, but also for $C_{\kappa}(t)$ in Fig. 3(d). It is indicative
that the fit is consistent with only two relaxation times for both
correlations, i.e., with longer relaxation rate $\tau^{s}_{\kappa}$ still
weakly dependent on $s$ as given in Fig. 3(a), and much faster
$\tau^{s}_{\sigma}\sim\tau^{0}_{\sigma}$ which we can approximate just with
the $s=0$ result for $C_{\sigma}(t)$ in Fig. 2(a).

Figure 3: (a): Relaxation time for the energy current obtained from long-time
fits marked by dashed lines in Fig. 2. (b,c): Continuous line shows spin-
current $C_{\sigma}(t)$ for $S^{z}_{\rm tot}=2,4$, respectively, while (d)
energy-current $C_{\kappa}(t)$ for $S^{z}_{\rm tot}=4$. The results in (b,d)
are fitted by a sum of two exponential functions (dashed, green curves) with
the relaxation times $\tau^{s}_{\kappa}$ as in panel (a), and
$\tau^{s}_{\sigma}=\tau^{0}_{\sigma}$. Calculated for $L=28$, $\Delta=0.5$,
and $g=0.15$.

The above results may serve as a motivation for a simple phenomenological
description. The value of the correlation functions for $t\to\infty$ in the
integrable model, $g=0$, is determined via the Mazur bound [2, 3] by the
projections of the studied operators on local/quasilocal CQ [33, 34, 35, 36,
37, 38, 39] $Q_{n}$, which should be chosen as mutually orthogonal $\langle
Q_{n}Q_{n^{\prime}}\rangle\propto\delta_{n,n^{\prime}}$ in the canonical
ensemble with fixed $S^{z}_{\rm tot}$. Assuming the completeness (saturation)
of the bound, this means for considered correlations,

$C_{\alpha}(t\to\infty)=C^{0}_{\alpha}=\frac{1}{\langle
j_{\alpha}j_{\alpha}\rangle}\sum_{n}\frac{\langle
j_{\alpha}Q_{n}\rangle^{2}}{\langle Q_{n}Q_{n}\rangle}\,,$ (4)

where the term $\langle j_{\alpha}j_{\alpha}\rangle^{-1}$ arises from
normalization in Eq. (3). Eq. (4) is invariant under the orthogonal
transformation of normalized CQ

$\frac{Q_{n}}{\sqrt{\langle
Q_{n}Q_{n}\rangle}}=\sum_{s}\left(\hat{O}\right)_{ns}\frac{Q_{s}}{\sqrt{\langle
Q_{s}Q_{s}\rangle}}\,,$ (5)

where $\hat{O}$ is arbitrary orthogonal matrix. In NI system, $Q_{n}$ are not
any more CQ and decay with the characteristic times $\propto g^{2}$. Based on
the results in Fig. 3(b-d) and previous numerical studies in Ref. [19], we
conjecture that projections on $Q_{n}$ are essential also for the asymptotic
dynamics in NI models,

$C_{\alpha}(t\gg J^{-1})\simeq\frac{1}{\langle
j_{\alpha}j_{\alpha}\rangle}\sum_{n}\frac{\langle
j_{\alpha}Q_{n}\rangle^{2}}{\langle
Q_{n}Q_{n}\rangle}\exp\left(-\frac{t}{\tau_{n}}\right)\,,$ (6)

where $\tau_{n}\propto 1/g^{2}$. It should be, however, stressed that in
contrast to Eq. (4) the appropriate set of $Q_{n}$ in Eq. (6) is not
arbitrary, but determined by the perturbation [17, 5, 27]. A numerical
algorithm for identifying such modes has been discussed in Ref. [19]. It
amounts in finding eigenvector of the matrix $M_{AB}=\langle A(t)B\rangle$
calculated in the long-time regime, $t\gg 1$, for a set of local orthogonal
operators, $A$ and $B$. Since the number of local operators supported on $M$
sites grows exponentially with $M$, this method is numerically demanding.
Additionally, it requires full diagonalization of the Hamiltonian hence it is
applicable only to small systems. Therefore in this work, we do not attempt to
determine relevant $Q_{n}$ in more detail.

However, for $S^{z}_{\rm tot}=0$ one can assume that one of $Q_{n}$ in Eq. (6)
can be well approximated by $j_{\kappa}$ (being one of CQ in the reference
system), which explains nearly perfect exponential decay of $C_{\kappa}(t)$ in
Fig. 2(a), even though a small correction might be needed (as tested in more
detail in Sec. 3). It is well known [3] that in the integrable model ($g=0$)
the spin-current stiffness, $C^{0}_{\sigma}$, is nonzero for $S^{z}_{\rm
tot}=0$ provided that $\Delta<1$. Since $\langle
j_{\sigma}j_{\kappa}\rangle=0$ for $S^{z}_{\rm tot}=0$, $Q_{n}$ that are
relevant for the decay time $\tau^{0}_{\sigma}$ have to be related to the
quasilocal CQ [33, 34, 40, 19].

On the other hand, for nonzero $S^{z}_{\rm tot}$, we clearly need at least two
distinct relaxation times $\tau^{s}_{\kappa},\tau^{s}_{\sigma}$ to fit both
decays, $C_{\alpha}(t)$. The necessity of multiple relaxation times is most
evident for $C_{\sigma}(t)$ at $S^{z}_{\rm tot}=2,4$ presented on Fig. 3(b-c).
Upon changing $S^{z}_{\rm tot}$, one has also to modify $\tau^{s}_{\kappa}$,
as suggested by the results in Fig. 3(a), while the dependence of
$\tau^{s}_{\sigma}$ on $S^{z}_{\rm tot}$ is less evident, so we can
approximate $\tau^{s}_{\sigma}\sim\tau^{0}_{\sigma}$. In Fig. 3(d) we show
that $C_{\kappa}(t)$ can also be well fitted by a sum of two exponential
functions, with the very same relaxation times which are used for fitting
$C_{\sigma}(t)$ in Fig. 3(c). This indicates that the same pair of $Q_{n}$ is
relevant for the decay of the spin and the energy currents for $S^{z}_{\rm
tot}\neq 0$.

##  3 Memory-function analysis of the energy-current correlations

It is desirable to have more explicit way to evaluate the long-time decay of
correlations, $C_{\alpha}(t)$, for chosen perturbations $H^{\prime}$, and in
particular a direct expression for relevant relaxation rates $1/\tau_{n}$.
Here, it is convenient to follow the Mori formalism [41, 42], where one can
generally express the current-correlation relaxation function
$\phi_{\alpha}(\omega)$ in terms of the corresponding memory functions (MF),
$M_{\alpha}(\omega)$,

$\displaystyle\phi_{\alpha}(\omega)$ $\displaystyle=$
$\displaystyle\frac{1}{L}\left(j_{\alpha}\left|\left[{\cal
L}-\omega\right]^{-1}\right|j_{\alpha}\right)=\frac{\chi_{\alpha}(\omega)-\chi^{0}_{\alpha}}{\omega}=\frac{-\chi^{0}_{\alpha}}{\omega+M_{\alpha}(\omega)}\,,$
$\displaystyle\chi_{\alpha}(\omega)$ $\displaystyle=$
$\displaystyle\frac{i}{L}\int_{0}^{\infty}e^{i\omega^{+}t}\langle[j_{\alpha}^{\dagger}(t),j_{\alpha}]\rangle\mathrm{d}t\,,\qquad\chi^{0}_{\alpha}=\chi_{\alpha}(0)>0\,,$
(7)

where $\left(A|B\right)=1/L\int_{0}^{\beta}{\rm d}\tau\langle
AB(i\tau)\rangle$. At high temperatures, $\beta\to 0$,
$\beta\tilde{C}_{\alpha}(\omega)=\mathrm{Im}\phi_{\alpha}(\omega)$ and
$\chi^{0}_{\alpha}=\beta\langle j_{\alpha}j_{\alpha}\rangle$. Knowing
numerical result for $\tilde{C}_{\alpha}(\omega)$ at finite $g$, one can
evaluate (via Kramers-Kronig relation) the complex $\phi_{\alpha}(\omega)$,
extract directly the corresponding complex MF, $M_{\alpha}(\omega)$, and in
particular the dynamical relaxation rate
$\Gamma_{\alpha}(\omega)=\mathrm{Im}\,M_{\alpha}(\omega)$.

On the other hand, for a NI system one can find an explicit expression in the
case where the current is a CQ in the reference system [17, 19]. This is
particularly the case for the energy current,
$[H_{0},Q_{3}]=[H_{0},j_{\kappa}]=0$. Then, within the lowest order in the
perturbation $g\ll 1$ one can approximate MF as

$M_{\kappa}(\omega)=\frac{1}{\chi^{0}_{\alpha}}N_{\kappa}(\omega),\qquad
N_{\kappa}(\omega)=g^{2}({\cal F}|[{\cal L}-\omega]^{-1}|{\cal
F})\,,\qquad{\cal F}=[H^{\prime},j_{\kappa}]\,,$ (8)

i.e., as the correlation function of the the force ${\cal F}$ in the
unperturbed - integrable system. Eq. (8) is derived under the assumption that
only a single charge $Q_{3}=j_{\kappa}$ is relevant for the decay, as well as
that ${\cal F}$ has no overlap with any of $Q_{n}$. This can be true for the
sector with $S^{z}_{\rm tot}=0$ and can be directly tested for operator ${\cal
F}$ for the particular $H^{\prime}$ in Eq. (1),

${\cal
F}=i\sum_{i}\left[T_{i-1}^{i+1}S^{z}_{i}\left(S^{z}_{i+3}-S^{z}_{i-3}\right)-\Delta
T_{i}^{i+1}\left(S^{z}_{i+2}+S^{z}_{i-1}\right)\left(S^{z}_{i+3}-S^{z}_{i+2}\right)\right]\,,$
(9)

with $T_{i}^{l}=(S^{+}_{l}S^{-}_{i}+\mathrm{H.c.})/2$. Due to the parity
symmetry at $S^{z}_{\rm tot}=0$, the above ${\cal F}$ appears orthogonal to
known local/quasilocal CQ. We can then evaluate numerically
$M_{\kappa}(\omega)$ as the correlation of ${\cal F}$ in the reference $H_{0}$
system. The MCLM result for $\Gamma_{\kappa}(\omega)/g^{2}$ at $S^{z}_{\rm
tot}=0$ is presented in Fig. 4, both as extracted directly via Eq. (7) from
$\tilde{C}_{\kappa}(\omega)$ for various finite $g=0.1-0.4$, as well as by
calculating the result from Eq. (8) for the force given in Eq. (9). The
overall agreement of perturbation result with the numerically extracted MF
$\Gamma_{\kappa}(\omega)/g^{2}$, (being essentially $g$-independent) is very
satisfactory in the whole $\omega$ regime [see the inset of Fig. 4(a)]. Still,
at low $\omega/J<0.4$ there appears a visible difference. Partly responsible
is the (apparently) singular contribution in perturbative
$\Gamma^{0}_{\kappa}=\Gamma_{\kappa}(\omega\to 0)$, which seems to indicate a
small overlap with some CQ (possibly being a finite-size effect). Apart from
that, also visible is a quantitative mismatch at $\omega\to 0$ which indicates
that the relevant $Q_{n}$ in Eq. (6) is not just $j_{\kappa}$, and
consequently also Eq. (8) is not a full description of relaxation. In other
words, even for $S^{z}_{\rm tot}=0$ more than a single mode is needed in Eq.
(6) to properly describe $C_{\kappa}(t)$.

Figure 4: Energy-current memory function (relaxation-rate),
$\Gamma_{\kappa}(\omega)$, (a) extracted directly from
$\tilde{C}_{\kappa}(\omega)$ for $S^{z}_{\mathrm{tot}}=0$ with different
perturbations $g=0.1-0.4$, compared with the perturbation result, Eqs.
(8),(9), and (b) extracted for various $S^{z}_{\mathrm{tot}}$ from
$\tilde{C}_{\kappa}(\omega)$ for fixed $g=0.15$. Calculated for $L=28$ and
$\Delta=0.5$.

In Fig. 4(b) we present $\Gamma_{\kappa}(\omega)/g^{2}$ for fixed $g=0.15$ but
various $S^{z}_{\rm tot}$. It is indicative that for $\omega/J>0.15$ the MF is
essentially independent of $S^{z}_{\rm tot}$. On the other hand, the decrease
of $\Gamma^{0}_{\kappa}=1/\tau^{s}_{\kappa}$ with $S^{z}_{\rm tot}$ reflects
the observed increase of $\tau^{s}_{\kappa}$ in Fig. 3(a). Still, it is not
straightforward to capture this in a perturbative approach, Eq. (8).

##  4 Dependence on the form of perturbation

In Sec. 2 we analysed a particular form of perturbation, i.e., the next
nearest-neighbor interaction, Eq. (1), which does not break the translational
symmetry or the spin parity $P$. This is also essential for the
phenomenological explanation, Eq. (6), in terms of decaying CQ, as well as for
the MF analysis in Sec. 3. However, using numerical MCLM we can check also
other perturbations. Let us consider as a perturbation

$H^{\prime\prime}=J\sum_{i}S^{z}_{i}S^{z}_{i+1}S^{z}_{i+2}\,,$ (10)

which breaks the parity symmetry, $P$, of the total Hamiltonian,
$H=H_{0}+gH^{\prime\prime}$. Fig. 5 shows current $\alpha=\sigma\,,\kappa$
correlation functions $C_{\alpha}(t)$ for $S^{z}_{\rm tot}=0$. We observe
exponential decay with a single relaxation time $\tau_{\alpha}\propto
1/g^{2}$, which is different from two distinct relaxation times obtained for
the $P$-preserving $H^{\prime}$, Eq. (1), i.e., the relaxation evidently
depends on the form of perturbation. One may interpret this behavior in terms
of conjecture, Eq. (6). All observables which have non-vanishing projection on
$Q_{n}$ with the longest relaxation time should asymptotically decay with the
same decay rate. In order to obtain different rate, one needs to build an
operator that is strictly orthogonal to the latter $Q_{n}$. It is highly
nontrivial task, since approximate $Q_{n}$ in Eq. (6) can be obtained
numerically only for small system. For the parity-preserving perturbations,
the total Hamiltonian is even and all $Q_{n}$ in Eq. (6) have well defined
parity being either odd or even. Then operators from one parity sectors are
strictly orthogonal to $Q_{n}$ from the other parity sector. Due to this
orthogonality we observe different relaxation times in odd and even sectors
without any fine-tuning of the studied observables. However for odd
perturbations, the total Hamiltonian contains even and odd terms, hence both
parity sectors are mixed during the time evolution. As a consequence, $Q_{n}$
in Eq. (6) may not have well defined parity.

Figure 5: Current correlation function $C_{\alpha}(t)$ for the perturbation,
Eq. (10), breaking the parity symmetry $P$ (odd perturbations), as compared
with $P$-preserving perturbation Eq. (1) (even perturbations). Calculated for
$L=28$, $\Delta=0.5$, $g=0.15$, and $S^{z}_{\mathrm{tot}}=0$.

##  5 Non-commuting conserved quantities and their non-exponential relaxation

In this section, we focus on specific values of the anisotropy parameter,
$\Delta$, where the many-body spectra exhibit additional massive degeneracies
[31]. The latter originate from eigenstates corresponding to different
$S^{z}_{\rm tot}$ with equal energies. This property allows for the presence
of non-commuting local/quasilocal CQ [43, 44, 32, 45]. As a first example we
take $\Delta=1$ in which case $H_{0}$ is the SU(2)-symmetric Heisenberg chain
for which one can study the total $S^{x}_{\rm tot}$ spin operator

$O_{1}=\sum_{j}(S^{+}_{j}+S^{-}_{j})\,.$ (11)

It is clear that $O_{1}$ is a local operator which commutes with $H_{0}$ at
$\Delta=1$, but does not commute with other local CQ, e.g., not with the
$S^{z}_{\rm tot}$.

The second nontrivial case is related with the commensurate
$\Delta=\cos(\pi/3)=1/2$, for which the non-commuting local CQ have been
derived (for $g=0$) in Ref. [32]. Here, we study the relaxation of

$O_{3}=\sum_{j}(-1)^{j}(S^{+}_{j-1}S^{+}_{j}S^{+}_{j+1}+{\rm H.c.})\,,$ (12)

which does not commute with $S^{z}_{\rm tot}$ and is not invariant under
translations by odd number of sites. It should be mentioned that analogous
local operators exist for other commensurate cases, in particular for
$\Delta=\cos(\pi/2)=0$, where
$O_{2}=\sum_{j}(-1)^{j}(S^{+}_{j}S^{+}_{j+1}+\mathrm{H.c})$ is local CQ for
$H_{0}$[45].

In similarity to Eq. (3), we calculate normalized correlations functions
$C_{1}(t)$ and $C_{3}(t)$ for the operators $O_{1}$ and $O_{3}$, respectively.
In contrast to the preceding section, now the averaging $\langle\dots\rangle$
is carried out over grand canonical ensemble. We note that $O_{3}$ is similar
to the previously studied spin-current $j_{\sigma}$ in the sense, that it does
not commute with $H_{0}$, but has large projection on the corresponding
quasilocal CQ [32]. In order to confirm this, in the inset in Fig. 6(c) we
show the finite-size scaling of the relevant stiffness
$\tilde{C}_{3}(\omega\to 0^{+},g=0)$ (with value $\simeq 0.4$ thermodynamic
limit $L\to\infty$). Without imposing the translational symmetry, we studied
all local non-commuting operators $\sum_{i}(A_{i}+A^{\dagger}_{i})$, where
$A_{i}$ are supported on up to $3$ sites and do not commute with $S^{z}_{\rm
tot}$. Utilizing the alhorithm from Ref. [40], we have found (for
$\Delta=1/2$, $g=0$ and $L\leq 14$) that $O_{3}$ has the largest stiffness out
of all these operators (not shown).

Next, we focus on the asymptotic decay of $C_{\alpha}(t)$, $\alpha=1,3$, and
their dependence on $g\neq 0$. To this end, we first calculate the Fourier
transform, $C_{\alpha}(\omega)$, and then its cumulative spectral function

$\tilde{C}_{\alpha}(\omega,g)=\int_{-\omega}^{\omega}\mathrm{d}\omega^{\prime}\;C_{\alpha}(\omega^{\prime})=\frac{\sum_{mn}\theta(\omega-|E_{m}-E_{n}|)\langle
m|O_{\alpha}|n\rangle^{2}}{\sum_{mn}\langle m|O_{\alpha}|n\rangle^{2}}\,,$
(13)

expressed in terms of eigenstates, $H|n\rangle=E_{n}|n\rangle$, obtained using
the exact diagonalization of $H$ on up to $L=16$ sites. We note that
$\tilde{C}_{1}(\omega\to 0^{+},g=0)=1$ (because $O_{1}$ commutes with $H_{0}$
at $\Delta=1$), whereas $\tilde{C}_{3}(\omega^{+}\to 0,g=0)\simeq 0.4$. Since
we are interested in the low-$\omega$ part of $C_{\alpha}(\omega)$, it is
convenient to normalize the spectral function,

$R_{\alpha}(\omega)=\frac{\tilde{C}_{\alpha}(\omega,g)}{\tilde{C}_{\alpha}(\omega\to
0^{+},g=0)}\,,$ (14)

so that one may directly compare relaxations of both considered observables.
The numerical results are shown in figures Figs. 6(a,b) and 6(d,e) for
$R_{3}(\omega)$ at $\Delta=1/2$ and $R_{1}(\omega)$ at $\Delta=1$,
respectively. In contrast to the spin- and energy-currents discussed in the
preceding sections, curves for various $g$ do not collapse when one rescales
the frequency by $g^{2}$, as it is shown in the insets in Fig. 6(b,e).
However, a convincing collapse may be obtained for the scaling $\omega/g$.
Moreover, the rescaled curves can be quite accurately fitted by the error
function, shown as the dashed curves in Fig. 6, implying
$\tilde{C}_{\alpha}(\omega)\propto\exp[-a\;(\omega/g)^{2}]$ where the
coefficient $a$ does not depend on $g$. Consequently, the decay is not
exponential but Gaussian, $C_{\alpha}(t)\propto\exp[-(t/\tau_{\alpha})^{2}]$,
where the characteristic relaxation rate, $1/\tau_{\alpha}\propto g$.

Figure 6: Cumulative spectral functions $R_{\alpha}(\omega)$, Eq. (14), as
calculated for $L=16$ system. We show the results for (a) $O_{3}$ with
$\Delta=1/2$ and (d) for $O_{1}$ with $\Delta=1$. Panels (b) and (e) show
respectively the same data but with rescaled frequency $\omega/g$ (main
panels) and $\omega/g^{2}$ (insets). Panels (c) and (f) show the results for
integrable system ($g=0$) but for shifted anisotropy parameters
$\Delta=0.5+\delta$ and $\Delta=1+\delta$, respectively. Dashed curves show
fittings with the error function. Inset in (c) shows finite-size scaling of
the stiffness $C_{3}(t\to\infty)$ for $\Delta=0.5$ and $g=0$.

The Gaussian relaxation occurs in the vicinity of $\Delta$ characterized by
additional degeneracies [45, 31], originating from that eigenstate with
different $S^{z}_{\rm tot}$ having the same energy. The perturbation breaks
integrability but also lifts the latter degeneracy. In order to disentangle
these two mechanisms, we also study the correlation functions,
$\tilde{C}_{\alpha}(\omega,g)$ for $g=0$ but with shifted $\Delta=1/2+\delta$
(for $\alpha=3$) and $\Delta=1+\delta$ (for $\alpha=1$). Then, the degeneracy
is lifted without destroying the integrability. The results are shown in Fig.
6(c,f) for both $\alpha=1,3$. One observes the same behavior as for the NI
system above, i.e., $C_{\alpha}(t)\propto\exp[-(t/\tau_{\alpha})^{2}]$, with
$1/\tau_{\alpha}\propto\delta$. This result explains additionally the origin
of the anomalous scaling of characteristic $\tau_{\alpha}$. The nonzero
stiffnesses, $\tilde{C}_{\alpha}(\omega\to 0^{+},g=0)$, emerges from states
$|m\rangle$ and $|n\rangle$ in Eq. (13) with different $S^{z}_{\rm tot}$, but
with equal energies $E_{m}=E_{n}$. The latter degeneracy is lifted either by
$g\neq 0$ or $\delta\neq 0$ already in the first order perturbation theory,
$|E_{m}=E_{n}|\propto g,\delta$ [45], leading to the scaling
$\omega/(g,\delta)$, shown in Fig. 6.

##  6 Conclusions

We numerically analyzed the decay of normalized correlation function
$C_{\alpha}(t)$ of different local quantities in the nearly integrable XXZ
model. Correlations generally reveal a fast drop at short times $tJ\sim O(1)$,
consistent with finite stiffnesses of studied quantities in the integrable
model $H_{0}$. A weak integrability breaking $g\ll 1$ then leads to further
slow - exponential-like - decay which is characterized with a single or
multiple relaxation times, all scaling with the perturbation strength
$\tau_{n}\propto 1/g^{2}$. The simplest cases appear to be the energy-current
$j_{\kappa}$ and spin-current $j_{\sigma}$ correlations at zero magnetization
$S^{z}_{\rm tot}=0$, where - due to symmetry - different CQ are involved in
the relaxation of both quantities, and consequently relevant $\tau_{n}$ are
quite distinct. Since $j_{\kappa}$ is by itself CQ within $H_{0}$, one can go
a step further and give an explicit memory-function analysis and perturbative
expression for the relaxation-rate spectral function
$\Gamma_{\kappa}(\omega)$. In this case the extracted and perturbative
$\Gamma_{\kappa}(\omega)$ match quite well in the whole $\omega$ range. Still,
some deviations in $\Gamma_{\kappa}(\omega\sim 0)=1/\tau^{0}_{\kappa}$ as well
as in $C_{\kappa}(t)$ at shorter $t/g^{2}$ indicate on possible multiple
relaxation times and more than one relevant $Q_{n}$ even in this case.

The existence of multiple relaxation times is becoming evident for $S^{z}_{\rm
tot}\neq 0$. The conjecture Eq. (6) concerning the presence of multiple
relaxation times $\tau_{n}$, which are linked with (appropriately rotated)
various CQ of the parent integrable model, is not in conflict with previous
results [20] which report a simple exponential relaxation. In order to
demonstrate at least two relaxation times, we have carefully selected
observables so that they have a large projection on quickly decaying CQ for
$S^{z}_{\rm tot}\neq 0$ ($Q_{n}$ relevant for $j_{\sigma}$ at $S^{z}_{\rm
tot}=0$), as well as much smaller projection on slowly decaying one
($Q_{3}=j_{\kappa}$). Namely, we made us of the result from Ref. [3] that the
projection $\langle j_{\sigma}j_{\kappa}\rangle$ is proportional to
$S^{z}_{\rm tot}$ thus the projection may be easily tuned via changing the
magnetization sector. Typically, the opposite holds true: one studies
observables which have largest projection on simplest CQ, which are supported
only on the few sites and the related $\tau_{n}$ are the longest relaxation
times [19]. As a consequences, small and fast decaying projections on more
complicated CQ in Eq. (6) may not be visible in the numerical results. Such
generic scenario of the multi-scale relaxation is consistent with numerical
results obtained in the present work as well as in Ref. [19]. Nevertheless, it
should be considered as conjecture that requires further studies and should be
verified for other nearly integrable systems.

Furthermore, our analysis indicates that the form, in particular the symmetry,
of the perturbation is relevant for the decay of $C_{\alpha}(t)$. In contrast
to quantities considered above at general anisotropy $\Delta$ and their decay,
the correlations $C_{\alpha}(t)$ of particular quantities $O_{l}$, being
conserved by $H_{0}$ only at commensurate $\Delta_{0}=\cos(\pi/m)$, but not
commuting with $S^{z}_{\mathrm{tot}}$, behave qualitatively different. Under
finite perturbation, which here can be introduced either by $g\neq 0$ or by
deviation $\Delta=\Delta_{0}+\delta$, we observe effectively a (non-
exponential) Gaussian-like decay of $C_{\alpha}(t)$ with different scaling of
characteristic decay time, i.e. $\tau_{\alpha}\propto 1/g$. The origin here is
the lifting of macroscopic degeneracy of many-body states at these particular
$\Delta_{0}$.

### Funding information

J.H. acknowledges the support by the Polish National Agency of Academic
Exchange (NAWA) under contract PPN/PPO/2018/1/00035. M.M. acknowledges the
support by the National Science Centre, Poland via projects
2020/37/B/ST3/00020. P.P. acknowledges the support by the project N1-0088 of
the Slovenian Research Agency. The numerical calculation were partly carried
out at the facilities of the Wroclaw Centre for Networking and Supercomputing.

## References

  * [1] B. Bertini, F. Heidrich-Meisner, C. Karrasch, T. Prosen, R. Steinigeweg and M. Žnidarič,  _Finite-temperature transport in one-dimensional quantum lattice models_ ,  Rev. Mod. Phys. 93, 025003 (2021),  10.1103/RevModPhys.93.025003.
  * [2] P. Mazur,  _Non-ergodicity of phase functions in certain systems_ ,  Physica 43, 533 (1969),  http://dx.doi.org/10.1016/0031-8914(69)90185-2.
  * [3] X. Zotos, F. Naef and P. Prelovšek,  _Transport and conservation laws_ ,  Phys. Rev. B 55, 11029 (1997),  10.1103/PhysRevB.55.11029.
  * [4] X. Zotos,  _High temperature thermal conductivity of two-leg spin- $1/2$ ladders_,  Phys. Rev. Lett. 92, 067202 (2004),  10.1103/PhysRevLett.92.067202.
  * [5] G. P. Brandino, J.-S. Caux and R. M. Konik,  _Glimmers of a quantum kam theorem: Insights from quantum quenches in one-dimensional bose gases_ ,  Phys. Rev. X 5, 041043 (2015),  10.1103/PhysRevX.5.041043.
  * [6] T. Prosen,  _Time evolution of a quantum many-body system: Transition from integrability to ergodicity in the thermodynamic limit_ ,  Phys. Rev. Lett. 80, 1808 (1998),  10.1103/PhysRevLett.80.1808.
  * [7] Y. Huang, C. Karrasch and J. E. Moore,  _Scaling of electrical and thermal conductivities in an almost integrable chain_ ,  Phys. Rev. B 88, 115126 (2013),  10.1103/PhysRevB.88.115126.
  * [8] F. H. L. Essler, S. Kehrein, S. R. Manmana and N. J. Robinson,  _Quench dynamics in a model with tuneable integrability breaking_ ,  Phys. Rev. B 89, 165104 (2014),  10.1103/PhysRevB.89.165104.
  * [9] M. Kollar, F. A. Wolf and M. Eckstein,  _Generalized gibbs ensemble prediction of prethermalization plateaus and their relation to nonthermal steady states in integrable systems_ ,  Phys. Rev. B 84, 054304 (2011),  10.1103/PhysRevB.84.054304.
  * [10] B. Bertini, F. H. L. Essler, S. Groha and N. J. Robinson,  _Prethermalization and thermalization in models with weak integrability breaking_ ,  Phys. Rev. Lett. 115, 180601 (2015),  10.1103/PhysRevLett.115.180601.
  * [11] K. Mallayya, M. Rigol and W. De Roeck,  _Prethermalization and thermalization in isolated quantum systems_ ,  Phys. Rev. X 9, 021027 (2019),  10.1103/PhysRevX.9.021027.
  * [12] P. Jung, R. W. Helmes and A. Rosch,  _Transport in almost integrable models: Perturbed heisenberg chains_ ,  Phys. Rev. Lett. 96, 067202 (2006),  10.1103/PhysRevLett.96.067202.
  * [13] B. Bertini, F. H. L. Essler, S. Groha and N. J. Robinson,  _Thermalization and light cones in a model with weak integrability breaking_ ,  Phys. Rev. B 94, 245117 (2016),  10.1103/PhysRevB.94.245117.
  * [14] 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),  10.1103/PhysRevLett.98.050405.
  * [15] A. C. Cassidy, C. W. Clark and M. Rigol,  _Generalized thermalization in an integrable lattice system_ ,  Phys. Rev. Lett. 106, 140405 (2011),  10.1103/PhysRevLett.106.140405.
  * [16] F. Lange, Z. Lenarčič and A. Rosch,  _Time-dependent generalized gibbs ensembles in open quantum systems_ ,  Phys. Rev. B 97, 165138 (2018),  10.1103/PhysRevB.97.165138.
  * [17] P. Jung and A. Rosch,  _Lower bounds for the conductivities of correlated quantum systems_ ,  Phys. Rev. B 75, 245104 (2007),  10.1103/PhysRevB.75.245104.
  * [18] P. Jung and A. Rosch,  _Spin conductivity in almost integrable spin chains_ ,  Phys. Rev. B 76, 245108 (2007),  10.1103/PhysRevB.76.245108.
  * [19] M. Mierzejewski, T. Prosen and P. Prelovšek,  _Approximate conservation laws in perturbed integrable lattice models_ ,  Phys. Rev. B 92, 195121 (2015),  10.1103/PhysRevB.92.195121.
  * [20] K. Mallayya and M. Rigol,  _Quantum quenches and relaxation dynamics in the thermodynamic limit_ ,  Phys. Rev. Lett. 120, 070603 (2018),  10.1103/PhysRevLett.120.070603.
  * [21] J. De Nardis, D. Bernard and B. Doyon,  _Hydrodynamic diffusion in integrable systems_ ,  Phys. Rev. Lett. 121, 160603 (2018),  10.1103/PhysRevLett.121.160603.
  * [22] E. Ilievski, J. De Nardis, M. Medenjak and T. Prosen,  _Superdiffusion in one-dimensional quantum lattice models_ ,  Phys. Rev. Lett. 121, 230602 (2018),  10.1103/PhysRevLett.121.230602.
  * [23] S. Gopalakrishnan, D. A. Huse, V. Khemani and R. Vasseur,  _Hydrodynamics of operator spreading and quasiparticle diffusion in interacting integrable systems_ ,  Phys. Rev. B 98, 220303(R) (2018),  10.1103/PhysRevB.98.220303.
  * [24] U. Agrawal, S. Gopalakrishnan, R. Vasseur and B. Ware,  _Anomalous low-frequency conductivity in easy-plane XXZ spin chains_ ,  Phys. Rev. B 101, 224415 (2020),  10.1103/PhysRevB.101.224415.
  * [25] V. B. Bulchandani, S. Gopalakrishnan and E. Ilievski,  _Superdiffusion in spin chains_ ,  J. Stat. Mech. 2021, 084001 (2021),  10.1088/1742-5468/ac12c7.
  * [26] A. J. Friedman, S. Gopalakrishnan and R. Vasseur,  _Diffusive hydrodynamics from integrability breaking_ ,  Phys. Rev. B 101, 180302 (2020),  10.1103/PhysRevB.101.180302.
  * [27] A. Bastianello, A. D. Luca and R. Vasseur,  _Hydrodynamics of weak integrability breaking_ ,  J. Stat. Mech. 2021, 114003 (2021),  10.1088/1742-5468/ac26b2.
  * [28] T. LeBlond, D. Sels, A. Polkovnikov and M. Rigol,  _Universality in the Onset of Quantum Chaos in Many-Body Systems_ ,  arXiv e-prints 2012.07849 (2020),  2012.07849.
  * [29] M. W. Long, P. Prelovšek, S. El Shawish, J. Karadamoglou and X. Zotos,  _Finite-temperature dynamical correlations using the microcanonical ensemble and the lanczos algorithm_ ,  Phys. Rev. B 68, 235106 (2003),  10.1103/PhysRevB.68.235106.
  * [30] P. Prelovšek and J. Bonča,  _Ground state and finite temperature lanczos methods_ ,  In A. Avella and F. Mancini, eds., _Strongly Correlated Systems \- Numerical Methods_. Springer, Berlin,  10.1007/978-3-642-35106-8 (2013).
  * [31] P. Prelovšek, M. Mierzejewski and J. Herbrych,  _Coexistence of diffusive and ballistic transport in integrable quantum lattice models_ ,  Phys. Rev. B 104, 115163 (2021),  10.1103/PhysRevB.104.115163.
  * [32] L. Zadnik, M. Medenjak and T. Prosen,  _Quasilocal conservation laws from semicyclic irreducible representations of $u_{q}(sl_{2})$ in xxz spin-1/2 chains_,  Nucl. Phys. B. 902, 339 (2016),  10.1016/J.NuclPhysB.2015.11.023.
  * [33] T. Prosen,  _Open $\mathrm{XXZ}$ spin chain: Nonequilibrium steady state and a strict bound on ballistic transport_,  Phys. Rev. Lett. 106, 217206 (2011),  10.1103/PhysRevLett.106.217206.
  * [34] T. Prosen and E. Ilievski,  _Families of quasilocal conservation laws and quantum spin transport_ ,  Phys. Rev. Lett. 111, 057203 (2013),  10.1103/PhysRevLett.111.057203.
  * [35] M. Mierzejewski, P. Prelovšek and T. Prosen,  _Breakdown of the generalized $\mathrm{Gibbs}$ ensemble for current-generating quenches_,  Phys. Rev. Lett. 113, 020602 (2014),  10.1103/PhysRevLett.113.020602.
  * [36] B. Pozsgay, M. Mestyán, M. A. Werner, M. Kormos, G. Zaránd and G. Takács,  _Correlations after quantum quenches in the $\mathrm{XXZ}$ spin chain: Failure of the generalized $\mathrm{Gibbs}$ ensemble_,  Phys. Rev. Lett. 113, 117203 (2014),  10.1103/PhysRevLett.113.117203.
  * [37] G. Goldstein and N. Andrei,  _Failure of the local generalized $\mathrm{Gibbs}$ ensemble for integrable models with bound states_,  Phys. Rev. A 90, 043625 (2014),  10.1103/PhysRevA.90.043625.
  * [38] R. G. Pereira, V. Pasquier, J. Sirker and I. Affleck,  _Exactly conserved quasilocal operators for the $\mathrm{XXZ}$ spin chain_,  J. Stat. Mech. 2014, P09037 (2014),  10.1088/1742-5468/2014/09/P09037.
  * [39] E. Ilievski, J. De Nardis, B. Wouters, J.-S. Caux, F. H. L. Essler and T. Prosen,  _Complete generalized gibbs ensembles in an interacting theory_ ,  Phys. Rev. Lett. 115, 157201 (2015),  10.1103/PhysRevLett.115.157201.
  * [40] M. Mierzejewski, P. Prelovšek and T. Prosen,  _Identifying local and quasilocal conserved quantities in integrable systems_ ,  Phys. Rev. Lett. 114, 140601 (2015),  10.1103/PhysRevLett.114.140601.
  * [41] H. Mori,  _Transport, collective motion and Brownian motion_ ,  Prog. Theor. Phys. 33, 423 (1965),  10.1143/PTP.33.423.
  * [42] P. Prelovšek and J. Herbrych,  _Self-consistent approach to many-body localization and subdiffusion_ ,  Phys. Rev. B 96, 035130 (2017),  10.1103/PhysRevB.96.035130.
  * [43] M. Fagotti,  _On conservation laws, relaxation and pre-relaxation after a quantum quench_ ,  J. Stat. Mech. 2014, P03016 (2014),  10.1088/1742-5468/2014/03/p03016.
  * [44] F. H. L. Essler and M. Fagotti,  _Quench dynamics and relaxation in isolated integrable quantum spin chains_ ,  J. Stat. Mech. 2016, 064002 (2016),  10.1088/1742-5468/2016/06/064002.
  * [45] M. Mierzejewski, J. Herbrych and P. Prelovšek,  _Ballistic transport in integrable quantum lattice models with degenerate spectra_ ,  Phys. Rev. B 103, 235115 (2021),  10.1103/PhysRevB.103.235115.