Integrable quenches in nested spin chains II: fusion of boundary transfer matrices

TL;DR

This paper analyzes the spectrum of integrable boundary transfer matrices in the $SU(3)$-invariant spin chain, deriving functional relations and applying them to compute the Loschmidt echo after quantum quenches.

Contribution

It provides a detailed spectral analysis of fused boundary transfer matrices using functional relations, advancing understanding of integrable quenches in nested spin chains.

Findings

Derived functional relations for eigenvalues of fused operators.

Computed Loschmidt echo for various quench protocols.

Connected boundary transfer matrix spectrum to physical observables.

Abstract

We consider quantum quenches in the integrable $SU(3)$-invariant spin chain (Lai-Sutherland model), and focus on the family of integrable initial states. By means of a Quantum Transfer Matrix approach, these can be related to "soliton-non-preserving" boundary transfer matrices in an appropriate transverse direction. In this work, we provide a technical analysis of such integrable transfer matrices. In particular, we address the computation of their spectrum: this is achieved by deriving a set of functional relations between the eigenvalues of certain "fused operators" that are constructed starting from the soliton-non-preserving boundary transfer matrices (namely the $T$- and $Y$-systems). As a direct physical application of our analysis, we compute the Loschmidt echo for imaginary and real times after a quench from the integrable states. Our results are also relevant for the study of the spectrum of $SU(3)$-invariant Hamiltonians with open boundary conditions.

Equations

$$H_L=\sum _{j=1}^L\left[{\mathbf{s}}_j\cdot {\mathbf{s}}_{j+1}+{\left({\mathbf{s}}_j\cdot {\mathbf{s}}_{j+1}\right)}^2\right]-2L,$$

$${\hat{R}}_{12}(\lambda )=\frac{1}{\lambda +i}\left(\lambda +i{\mathcal{P}}_{12}\right),$$

$${\mathcal{P}}_{12}|a{⟩}_1\otimes |b{⟩}_2=|b{⟩}_1\otimes |a{⟩}_2,$$

$${\hat{R}}_{12}(0)={\mathcal{P}}_{12},$$

$${\hat{R}}_{12}(u){\hat{R}}_{21}(-u)=1,$$

$${\hat{R}}_{21}(u)={\mathcal{P}}_{12}{\hat{R}}_{12}(u){\mathcal{P}}_{12}.$$

$$R_{12}(\lambda )=(\lambda +i){\hat{R}}_{12}(\lambda )=\lambda +i{\mathcal{P}}_{12},$$

$$R_{12}(\lambda )R_{21}(-\lambda )=\zeta (\lambda ),$$

$$R_{12}^{t_1}(\lambda )R_{21}^{t_1}(-\lambda -3i)=\bar{\zeta }(\lambda +3i/2),$$

$$\zeta (\lambda )=(\lambda +i)(-\lambda +i),\bar{\zeta }(\lambda )=(\lambda +3i/2)(-\lambda +3i/2).$$

$$A^t=V^{-1}{A}^{T}V,$$

$$V=\left(\begin{array}[]{ccc}0&0&1\\ 0&1&0\\ 1&0&0\end{array}\right)\,.$$

$${\bar{R}}_{12}(\lambda ):=R_{12}^{t_1}\left(-\lambda -\frac{3i}{2}\right).$$

$${\bar{R}}_{12}(\lambda ){\bar{R}}_{21}(-\lambda )=\bar{\zeta }(\lambda ),$$

$${\bar{R}}_{12}^{t_1}(\lambda ){\bar{R}}_{12}^{t_1}(-\lambda -3i)=\zeta (\lambda +3i/2).$$

$$R_{12}(\lambda )R_{13}(\lambda +\mu )R_{23}(\mu )=R_{23}(\mu )R_{13}(\lambda +\mu )R_{12}(\lambda ),$$

$${\bar{R}}_{12}(\lambda ){\bar{R}}_{13}(\lambda +\mu )R_{23}(\mu )=R_{23}(\mu ){\bar{R}}_{13}(\lambda +\mu ){\bar{R}}_{12}(\lambda ).$$

$$t(\lambda )$$

$$\bar{t}(\lambda )$$

$$\left[t(\lambda ),t(\mu )\right]=\left[\bar{t}(\lambda ),t(\mu )\right]=\left[\bar{t}(\lambda ),\bar{t}(\mu )\right]=0.$$

$$H_{L}=i\frac{\partial}{\partial\lambda}\ln t(\lambda)\Bigr{|}_{\lambda=0}\,.$$

$${\left(\frac{k_j+i/2}{k_j-i/2}\right)}^L=\prod _{\genfrac{}{}{0ex}{}{p=1}{p\ne j}}^N\frac{k_j-k_p+i}{k_j-k_p-i}\prod _{\ell =1}^M\frac{{\lambda }_{\ell }-k_j+i/2}{{\lambda }_{\ell }-k_j-i/2},j=1,\dots ,N,$$

$$1=\prod _{j=1}^N\frac{k_j-{\lambda }_{\ell }-i/2}{k_j-{\lambda }_{\ell }+i/2}\prod _{\genfrac{}{}{0ex}{}{m=1}{m\ne \ell }}^M\frac{{\lambda }_{\ell }-{\lambda }_m-i}{{\lambda }_{\ell }-{\lambda }_m+i},\ell =1,\dots ,M.$$

$$t(\lambda )|\{k_j\},\{{\lambda }_j\}⟩=\nu (\{k_j\},\{{\lambda }_j\},\lambda )|\{k_j\},\{{\lambda }_j\}⟩,$$

$$\nu \left(\{k_j\},\{{\lambda }_j\},\lambda \right)={\left[a(\lambda )\right]}^L\prod _{j=1}^N\frac{1}{a(\lambda -k_j+i/2)}{\nu }_1\left(\{k_j\},\{{\lambda }_j\},\lambda \right)$$

$$+\prod _{j=1}^N\frac{1}{a(k_j-i/2-\lambda )},$$

$${\nu }_1(\{k_j\},\{{\lambda }_j\},\lambda )=\prod _{j=1}^{N}a(\lambda -k_j+i/2)\prod _{r=1}^M\frac{1}{a(\lambda -{\lambda }_r+i/2)}+\prod _{r=1}^M\frac{1}{a({\lambda }_r-i/2-\lambda )},$$

$$a(\lambda )=\frac{\lambda }{\lambda +i}.$$

$$\tau (\lambda )$$

$$\bar{\tau }(\lambda )$$

$$T_a(\lambda )=R_{aN}(\lambda -{\xi }_N){\bar{R}}_{a(N-1)}(\lambda -{\xi }_{N-1})\dots R_{a2}(\lambda -{\xi }_2){\bar{R}}_{a1}(\lambda -{\xi }_1),$$

$${\hat{T}}_{\bar{a}}(\lambda )=R_{1a}(\lambda +{\xi }_1){\bar{R}}_{2a}(\lambda +{\xi }_2)\dots R_{(N-1)a}(\lambda +{\xi }_{N-1}){\bar{R}}_{Na}(\lambda +{\xi }_N),$$

$$T_{\bar{a}}(\lambda )={\bar{R}}_{aN}(\lambda -{\xi }_N)R_{a(N-1)}(\lambda -{\xi }_{N-1})\dots {\bar{R}}_{a2}(\lambda -{\xi }_2)R_{a1}(\lambda -{\xi }_1),$$

$${\hat{T}}_a(\lambda )={\bar{R}}_{1a}(\lambda +{\xi }_1)R_{2a}(\lambda +{\xi }_2)\dots {\bar{R}}_{(N-1)a}(\lambda +{\xi }_{N-1})R_{Na}(\lambda +{\xi }_N).$$

$$R_{ab}(\lambda -\mu )K_a(\lambda ){\bar{R}}_{ba}(\lambda +\mu )K_b(\mu )=K_b(\mu ){\bar{R}}_{ab}(\lambda +\mu )K_a(\lambda )R_{ba}(\lambda -\mu ),$$

$${\bar{R}}_{ab}(\lambda -\mu )K_{\bar{a}}(\lambda )R_{ba}(\lambda +\mu )K_b(\mu )=K_b(\mu )R_{ab}(\lambda +\mu )K_{\bar{a}}(\lambda ){\bar{R}}_{ba}(\lambda -\mu ).$$

Full text

1. 1 Introduction

2. 2 The general setting

3. 2.1 The periodic algebraic Bethe ansatz

4. 2.2 The boundary algebraic Bethe ansatz

5. 2.3 The integrable states

6. 3 Fusion relations of boundary transfer matrices

7. 3.1 $T$- and $Y$-systems for periodic transfer matrices

8. 3.2 The $T$-system for open transfer matrices

9. 3.3 Fusion of $R$-matrices

10. 3.4 Fusion of $K$-matrices

11. 3.5 Product of transfer matrices

12. 3.6 The final result

13. 4 The $Y$-system and the Loschmidt echo

14. 4.1 The non-linear integral equations for the delta-state

15. 4.2 The analytic solution to the $Y$-system at vanishing times

16. 4.3 The Loschmidt echo

17. 5 Conclusions

18. A Inversion relations of transfer matrices

19. References

1 Introduction

This is the second of two works devoted to the study of quantum quenches in the $SU(3)$-invariant Lai-Sutherland chain [2, 3]

[TABLE]

where the spin-$1$ operators $s^{a}_{j}$ are given by the standard three-dimensional representation of the $SU(2)$ generators. This integrable model is interesting from the physical point of view, as it displays different quasi-particle species, each one forming an infinite number of bound states. While many properties of the system at zero and finite temperature are by now well understood [6, 7, 4, 5, 8, 9], including the knowledge of its correlation functions [10, 11], until recently the analytical study of quench problems in this model has been out of our reach. The main reason for this lies in the fact that the solution to the Hamiltonian (1) involves a complicated nested Bethe ansatz [12, 13], for which generalizations of recent analytic advances in integrability out of equilibrium [14] are not straightforward, including the string-charge duality [15, 16, 17, 18, 19] or the Quench Action method [20, 21, 22, 23, 24]. Indeed, until very recently no initial state was known for which an explicit characterization of the corresponding post-quench steady state could be achieved. We note that, while it is established that for non-nested integrable systems the latter always corresponds to a generalized Gibbs ensemble (GGE) [25, 26, 27, 28, 29, 31, 30, 32], the picture is less clear in nested systems, where even a satisfactory understanding of the structure of quasi-local conservation laws is missing.

In the first of our works [33], it was shown that one can successfully tackle the study of quantum quenches in the nested chain (1) by means of an approach based on the Quantum Transfer Matrix (QTM) formalism [34, 35]. The latter is complementary to the string-charge duality and the Quench Action method, and was initially introduced for non-nested $XXZ$ Heisenberg chains in [36] (see also [37]). The starting point of this approach is to identify a family of initial integrable states [38] which, in close analogy with classical results in the context of quantum field theory [39, 40, 41], can be related to integrable boundary conditions in an appropriate transfer direction (or channel). More precisely, to each initial integrable state one can associate a given transfer matrix with open boundary conditions; in the case of the $SU(3)$-invariant Hamiltonian (1) these were shown to correspond to the “soliton-non-preserving” case studied in [42, 43, 44].

Importantly, the steady state reached at large times after a quench from integrable states could be characterized in [33] in terms of the corresponding quasi-particle distribution functions, which were computed analytically. This result is non-trivial: indeed, except for the special matrix product states found in [45, 46] , whose quench dynamics was considered in [47], no other example is known where the post-quench steady state can be determined, even approximately, for nested Hamiltonians. The derivation of [33] was based on a comparison between the QTM and Quench Action methods, and exploited as a technical ingredient certain “fusion relations” of the $SU(3)$-invariant boundary transfer matrices associated with the integrable states. These relations generalize to the nested case the ones derived in [48, 49] for the $SU(2)$-invariant case (and the corresponding $q$-deformations).

Albeit important, the derivation of such fusion relations is very technical, and was beyond the scope of Ref. [33]. This task is thus carried out in this more specialized work, where the fusion relations (and the corresponding $Y$-system [50]) are derived and analyzed in detail. In turn, these findings in principle give us access to the full spectrum of the boundary transfer matrix [51], and are thus relevant also for the study of $SU(3)$-invariant Hamiltonians with open boundary conditions. As an immediate physical application in the context of non-equilibrium physics, we address the computation of the Loschmidt echo [52] after a quench from integrable states, which is both interesting at real times and imaginary times, due to its connection to the study of dynamical quantum phase transitions [53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78] and the statistics of the work performed by the quench [79, 80], respectively. We work it out exactly at imaginary time and for small real times. Altogether, the calculations presented in [33] and in this work, provide a comprehensive analysis of integrable quenches in the $SU(3)$-invariant spin chain (1). Our findings are likely to have ramifications in the study of other nested integrable systems such as, for instance, multi-component Fermi and Bose ultra-cold gases, which are of special relevance for cold-atom experimental realizations [81, 83, 82].

The organization of this article is as follows. In Sec. 2 we introduce the algebraic Bethe ansatz, and recall the construction of the integrable states carried out in [33]. We introduce in particular the corresponding soliton-non-preserving boundary transfer matrices, whose fusion relations are derived explicitly in Sec. 3. In Sec. 4 we use them to compute the Loschmidt echo, both at imaginary and real times. Our conclusions are consigned to Sec. 5.

2 The general setting

2.1 The periodic algebraic Bethe ansatz

We begin by recalling the technical tools which are needed in order to analyze the Hamiltonian (1). In [33], we reviewed the main aspects of the coordinate nested Bethe ansatz, with which one is able to directly write down the wave functions corresponding to the eigenstates of the Hamiltonian (1). While this approach provides a clear physical picture, a quantitative analysis is better carried out by means of its algebraic version, which we briefly review. The main object of the theory is the $R$-matrix $\widehat{R}_{12}(\lambda)$ which acts on the tensor product $h_{1}\otimes h_{2}$, with $h_{j}\simeq\mathbb{C}^{3}$. Explicitly, it reads

[TABLE]

where $\mathcal{P}_{12}$ is the permutation matrix exchanging the spaces $h_{1}$ and $h_{2}$

[TABLE]

and which can be written in terms of the spin-1 operators as $\mathcal{P}_{12}=-1+\vec{S}_{1}\cdot\vec{S}_{2}+(\vec{S}_{1}\cdot\vec{S}_{2})^{2}$. The $R$-matrix satisfies the following properties of regularity and unitarity

[TABLE]

where

[TABLE]

Note that one has simply $\widehat{R}_{12}(u)=\widehat{R}_{21}(u)$. We also introduce the $R$-matrix with a different normalization

[TABLE]

which satisfies

[TABLE]

where

[TABLE]

Here we introduced the transposition

[TABLE]

where $A^{T}$ is the usual transposed of the matrix $A$, while

[TABLE]

Moreover, the notation $A^{t_{1}}$ in (9) signals the partial transposition, which is performed only on the indices refering to $h_{1}$ in the tensor product $h_{1}\otimes h_{2}$. In the following, it will be also necessary to introduce the $R$-matrix involving the conjugate representation of $SU(3)$ [84, 85, 86]. The latter is given as follows

[TABLE]

The $R$-matrix $\bar{R}_{12}(\lambda)$ can be interpreted physically as the scattering matrix describing the interaction between a soliton and an anti-soliton, and fulfills

[TABLE]

The $R$-matrices introduced above satisfy the Yang-Baxter equations, which take the form

[TABLE]

These fundamental relations have several immediate consequences. Indeed, defining the transfer matrices acting on the global Hilbert space $\mathcal{H}=h_{1}\otimes\ldots\otimes h_{N}$ as

[TABLE]

it follows that

[TABLE]

The traces in (18), (19) are taken over the auxiliary space $h_{j}\simeq\mathbb{C}^{3}$. The importance of this relation becomes manifest when complemented with the following trace formula

[TABLE]

From (20) and (21) it is clear that in order to diagonalize the Hamiltonian $H_{L}$ it is sufficient to find the eigenspectrum of $t(\lambda)$. The algebraic Bethe ansatz directly provides us with the tools to do exactly this. By a standard procedure, it is possible to prove [12] that the eigenvalues of $t(\lambda)$ are labeled by two sets of rapidities $\{k_{j}\}_{j=1}^{N}$, $\{\lambda_{j}\}_{j=1}^{M}$ satisfying the Bethe equations

[TABLE]

In particular, by denoting with $\ket{\{k_{j}\},\{\lambda_{j}\}}$ the corresponding eigenvectors we have

[TABLE]

where

[TABLE]

with the additional definitions

[TABLE]

and

[TABLE]

2.2 The boundary algebraic Bethe ansatz

The transfer matrices $t(\lambda)$ and $\bar{t}(\lambda)$ defined in (18) and (19) are translationally invariant operators corresponding to periodic boundary conditions. It is possible to define analogous operators in the case where open boundaries are assumed, via the boundary algebraic Bethe ansatz [87]. As we have already seen in [33], in the case of higher rank algebras $SU(\mathcal{N})$ with $\mathcal{N}\geq 3$, there exist two inequivalent boundary transfer matrices, which correspond to the so-called soliton-preserving [8] and soliton-non-preserving [42] boundary conditions [43]. In this work we will be interested in the latter, which are reviewed in the following.

We recall that the soliton-non-preserving boundary transfer matrices acting on a chain of length $N$, composed of an alternating product of fundamental and conjugate representations, are defined as

[TABLE]

where the following definitions are used

[TABLE]

Here we introduced the free parameters $\xi_{j}$ (called inhomogeneities) while the trace in (28) and (29) is taken over the auxiliary space $h_{a}\simeq\mathbb{C}^{3}$. Finally, $K_{a}^{\pm}(\lambda)$, $K_{\bar{a}}^{\pm}(\lambda)$ are $3\times 3$ matrices.

Importantly, the $K$-matrices $K_{a}^{\pm}(\lambda)$, $K_{\bar{a}}^{\pm}(\lambda)$ have to be chosen in such a way that transfer matrices with different spectral parameter commute. This can be done by setting $K_{a}^{-}(\lambda)=K_{a}(\lambda)$, $K_{\bar{a}}^{-}(\lambda)=K_{\bar{a}}(\lambda)$, $K_{a}^{+}(\lambda)=K^{t}_{a}(-\lambda-i3/2)$ and $K_{\bar{a}}^{+}(\lambda)=K^{t}_{\bar{a}}(-\lambda-i3/2)$, where $K_{a}(\lambda)$ and $K_{\bar{a}}(\lambda)$ are a solution to the following twisted boundary Yang-Baxter (or reflection) equations

[TABLE]

A classification of all invertible $K$-matrices satisfying the reflection equations was performed in [43]. In particular, it was shown that for the $SU(3)$-case the only invertible solutions to the latter equations are scalar matrices $K^{-}(\lambda)=K^{-}$ such that $(K^{-})^{t}=K^{-}$, namely such that $\tilde{K}^{-}=K^{-}V$ is symmetric. The most general matrix of this form is written as

[TABLE]

where $K_{s}$ a symmetric numerical matrix

[TABLE]

with ${\rm det}K_{s}\neq 0$.

The classification performed in [43] does not treat non-invertible matrices. It follows from continuity that any symmetric matrix with ${\rm det}K_{s}=0$ also satisfies the reflection equations, but here we discard these solutions. On the one hand, the fusion relations to be treated below depend crucially on a non-vanishing determinant. On the other hand, solutions with ${\rm det}K_{s}=0$ can be rotated to lie in an $SU(2)$ sector, and thus the earlier results in the literature can be used to solve these cases (see also Ref. [33]).

2.3 The integrable states

As we have seen in [33], there is a close relation between the soliton-non-preserving boundary transfer matrices introduced above and a special class of product states, namely the integrable states. We recall that the latter were initially introduced in [38] as those matrix product states [88] (with a finite bond-dimension) which are annihilated by all the parity-odd conserved charges of the Hamiltonian. In the case of the $XXZ$ Heisenberg chain, it was shown how to derive an infinite family of such states from integrable boundary conditions in an appropriate rotated channel. The constructions of [38] were generalized to the $SU(3)$-invariant spin chain (1) in [33], and a class of integrable product states was explicitly derived.

The family of product states found in [33] reads

[TABLE]

where the two-site state has to be chosen as

[TABLE]

Here we denoted with $|1\rangle$, $|2\rangle$ and $|3\rangle$ the basis vectors of the locals spaces $h_{j}\simeq\mathbb{C}^{3}$ and used the convention

[TABLE]

The coefficients $\kappa_{ij}$ have to be such that the matrix (37) has a non-vanishing determinant. It was shown in [33] that the product states with the building block (39) are integrable according to the definition of [38]. Furthermore, it was argued that, up to global $SU(3)$ transformations, one can always restrict themselves to the simpler case of diagonal boundary conditions

[TABLE]

with

[TABLE]

which includes as a particular case the following $SO(3)$-invariant “delta-state”

[TABLE]

As explained in [33], the relation between integrable states and transfer matrices with soliton-non-preserving boundary conditions can be derived by studying the following partition function

[TABLE]

Indeed, using the following Suzuki-Trotter decomposition [34, 10]

[TABLE]

one can show that, for a given integrable initial state with two-site block (39), the partition function (44) can be expressed as

[TABLE]

where

[TABLE]

Here $\tau^{s}_{N}(\lambda)$ is a soliton-non-preserving boundary transfer matrix acting on a chain of $N$ sites, which is related to (28) by a shift in the spectral parameters and rapidities

[TABLE]

Namely

[TABLE]

where

[TABLE]

with inhomogeneities

[TABLE]

Note that with the same shift (48), the transfer (29) is rewritten as

[TABLE]

where now

[TABLE]

The $K$-matrices $K^{\pm}$ in (49) and (54) are given by $K^{-}=K$, $K^{+}=K^{t}$, where $K$ is defined in (36).

Eqs. (46) and (47) provide the bridge between quantum quenches from integrable states and soliton-non-preserving boundary transfer matrices. In fact, as first observed in [36], several properties of the quench can be deduced from the study of such transfer matrices. Most prominently, we will show in Sec. 4.3 how the leading eigenvalue of $\tau^{s}_{N}(\lambda)$ directly yields the Loschmidt echo after the quench, which is closely related to the partition function $\mathcal{Z}(\beta)$. In the next section we will derive a very important set of functional relations for the eigenvalues of the transfer matrix $\tau^{s}_{N}(\lambda)$. On the one hand, these can be used to compute the spectrum of the latter; on the other hand, as we showed in Ref. [33], they can also be used to derive the quasi-particle distribution functions of the steady state reached at long times.

3 Fusion relations of boundary transfer matrices

The computation of the spectrum of $\tau^{s}_{N}(\lambda)$ for generic boundary conditions is a non-trivial problem. In fact, even for the well-known case of $XXZ$ Heisenberg chains, this has been solved only recently for finite system sizes [89, 90, 91, 92, 93, 94, 95]. Furthermore, non-trivial challenges arise for the study of the large-$N$ limit. A convenient strategy to address the computation of the spectrum of generic boundary transfer matrices, both for finite and infinite system sizes, was pursued in [36], by exploiting a set of functional relations which go under the name of $T$-system [50] and which is obtained by an appropriate “fusion procedure” [96]. In practice, one exploits the fact that the boundary transfer matrix can be embedded into a family of commuting operators, such that functional equations between their eigenvalues can be established, and eventually solved.

In the case of $XXZ$ spin-$1/2$ chains, the $T$-system for boundary transfer matrices was worked out in [48, 49], while to our knowledge it has not been explicitly obtained yet for $SU(3)$-invariant soliton-non-preserving boundary transfer matrices. This task is carried out in this section. Applications of this result will be presented in the subsequent sections.

3.1 $T$- and $Y$-systems for periodic transfer matrices

It is useful to start our discussion from the case of periodic transfer matrices. It is an established result [97, 98] that the operators $t(\lambda)$ defined in (18) can be embedded into a family of commuting transfer matrices $\{t_{m}^{(a)}(\lambda)\}$, with $m=1,2,\ldots+\infty$ and $a=1,2$, such that

[TABLE]

These satisfy a set of functional relations which go under the name of $T$-system [9]

[TABLE]

with the convention

[TABLE]

and where $\bar{t}(\lambda)$ is defined in (19). The operators $t^{(a)}_{m}(\lambda)$ can be explicitly obtained from $t(\lambda)$ and $\bar{t}(\lambda)$ by means of a geometrical construction called fusion [96]. For this reason, the operators $t^{(a)}_{m}(\lambda)$ are called fused transfer matrices.

Each transfer matrix $t^{(a)}_{m}(\lambda)$, can be written as (18), but with a different irreducible representation of $SU(3)$ in the auxiliary space. We recall that each irreducible representation is labeled by a pair $(m_{1},m_{2})$ of two non-negative integers [99]. The $m$-fold symmetric tensor of the fundamental representation is associated with $(m,0)$, while the $m$-fold symmetric tensor of its conjugate representation is associated with the pair $(0,m)$. The dimension of a representation corresponding to the pair $(m_{1},m_{2})$ is given by the formula [99]

[TABLE]

For each irreducible representation $(n,m)$, we label as $t_{(n,m)}(\lambda)$ the transfer matrix which has that representation in the auxiliary space. Note that, in the $SU(2)$ case, irreducible representations are instead labeled by a single index, and so all the possible transfer matrices are labeled by a single index. The operators $t^{(a)}_{m}(\lambda)$ entering the $T$-system are identified as follows

[TABLE]

The validity of the $T$-system (58), (59) was proven in [98]. The functions $\Phi_{n}^{(1)}(u)$, $\Phi_{n}^{(2)}(u)$ are scalar and depend on the physical situation. In the case of thermal quantum transfer matrices, their explicit form can be found in [9]. In any case, they can be seen to satisfy the relations

[TABLE]

with

[TABLE]

From the $T$-system encoded in (58) and (59), it is possible to derive another crucial set of functional relations which is usually called the $Y$-system [50]. The latter is important in particular to obtain a description in terms of integral equations for functions defined on the real line. Following e.g. [9] we define the $Y$-functions

[TABLE]

which are easily seen to satisfy the $Y$-system

[TABLE]

with the convention

[TABLE]

In the next sections, we will see how to transform these kinds of functional relations into integral equations.

3.2 The $T$-system for open transfer matrices

Analogously to the periodic case, also transfer matrices with open boundary conditions can be embedded into a family of fused operators $\tau_{(m_{1},m_{2})}(\lambda)$, characterized by a $SU(3)$ representation in the auxiliary space, corresponding to the pair $(m_{1},m_{2})$. As in the previous section we define

[TABLE]

with the convention

[TABLE]

and where $\tau(\lambda)$ and $\bar{\tau}(\lambda)$ are defined in (28) and (29).

In the case of $XXZ$ spin-$1/2$ Heisenberg chains, fusion relations for open transfer matrices were derived in [49], where it was found that the latter satisfy the same $T$-system of periodic transfer matrices. A direct generalization of the work [49] for the $SU(3)$-invariant spin chain is missing in the literature, even though some particular fused transfer matrices have been constructed, both for the soliton-preserving and soliton-non-preserving cases. For example, in the paper [43] the transfer matrix $\tau_{(1,1)}(u)$ was constructed in the soliton-non-preserving case from fusion of $\tau_{(1,0)}(u)$ and $\tau_{(0,1)}(u)$.

The purpose of this section is to write down the $T$-system in the case of boundary transfer matrices in the soliton-non-preserving case. In analogy with the spin-$1/2$ case, we assume that the latter is of the same form of the periodic case encoded in (58) and (59). Hence, we only need to determine the auxiliary functions $\Phi_{n}^{(a)}(\lambda)$. In fact, given the recursive relations (66) and (67), the whole calculation boils down to determining the functions $\Phi_{1}^{(1)}(\lambda)$ and $\Phi_{1}^{(2)}(\lambda)$. In order to obtain the latter, we compute directly the products $\tau(\lambda)\tau(\mu)$ and $\bar{\tau}(\lambda)\bar{\tau}(\mu)$. As we will see in the following, by means of the fusion of $R$- and $K$-matrices introduced in [42, 43], one can explicitly see that these products decompose exactly as in the functional relations (58) and (59) for $n=1$.

In the following calculations we will restrict to the case of invertible $K$-matrices, which are solution to the reflection equations (34) and (35). In particular, according to the results of [43], we know that the latter are of the form (36), which will be explicitly exploited.

3.3 Fusion of $R$-matrices

As a first ingredient, we need to review the fusion properties of the $R$-matrix [96, 43]. We note that $R_{23}(-i)=P_{\mathcal{H}_{3}}$, where $P_{\mathcal{H}_{3}}$ is the projector onto the three-dimensional space in the decomposition

[TABLE]

where we denote with $d$ and $\bar{d}$ the irreducible fundamental and conjugate representation of dimension $d$ respectively. Consider now the Yang-Baxter equations (16), (17). Choosing $\mu=-i$ it is immediate to see from (17) that the operators

[TABLE]

leave stable the symmetric six-dimensional representation in the decomposition (79), which implies

[TABLE]

where $P_{ab}^{-}=P_{\mathcal{H}_{3}}$ while $P_{ab}^{+}$ is the projector on the six-dimensional symmetric representation in the decomposition (79). Analogously, one can see that the operators

[TABLE]

leave stable the anti-symmetric three-dimensional representation in the decomposition (79). Accordingly, one can define the fused $R$-matrices

[TABLE]

where two indices between brackets stand for one single index labeling the six-dimensional representation. These fused $R$-matrices also satisfy the Yang-Baxter relations, as it can be easily seen for example through a graphical proof.

3.4 Fusion of $K$-matrices

Analogously to the case of $R$-matrices, one can also address the fusion of boundary reflection matrices [100]. We consider the soliton-non-preserving case, and start with the reflection equation (34). Setting $\mu=\lambda+i$ we find

[TABLE]

from which we obtain straightforwardly that the operator

[TABLE]

leaves stable the symmetric six-dimensional representation in the decomposition (79). Analogously, the operator

[TABLE]

leaves stable the anti-symmetric three dimensional representation in the decomposition (79). Note that in the definition of the boundary transfer matrix, one has the matrices

[TABLE]

Accordingly, using the considerations above, one can see that

[TABLE]

Hence, as in the previous section, one can define the fused $K$-matrices

[TABLE]

which are seen, once again, to satisfy an appropriate fused version of the Yang-Baxter reflection equations.

3.5 Product of transfer matrices

We can now make a step forward and compute the product of two transfer matrices $\tau(\lambda)$, as defined in (28). We start with the identity

[TABLE]

where we used the crossing unitarity

[TABLE]

We now make use of the following property of the partial transpose. Let $M$ be an operator acting on the product $\mathcal{H}_{a}\otimes\mathcal{H}_{b}\otimes\mathcal{H}_{p}$, and $\mathcal{B}$ an operator acting on $\mathcal{H}_{b}\otimes\mathcal{H}_{p}$. Then

[TABLE]

Applying (98) to Eq. (96) we get

[TABLE]

Next, define

[TABLE]

so that

[TABLE]

where we omitted the dependence of $\mathcal{K}_{ab}^{\pm}$ on the spectral parameters. Finally, we compute

[TABLE]

where in the second and third line we have introduced contracted indices $\alpha,\beta,\varphi,a_{2},b_{2},p_{2}$ spanning the basis of the auxilliary space. The final result of this computation is extremely simple and reads

[TABLE]

As a second step of the derivation, we show that the product (104) decomposes into the sum of two terms, provided that $\lambda$ and $\mu$ are chosen appropriately. Our goal is to extract the function $\Phi^{(1)}_{1}(\lambda)$ which is needed to write down the $T$-system for general $n$. First, we define the fused transfer matrices $\tau^{(1)}_{2}(\lambda)$ and $\tau^{(2)}_{2}(\lambda)$ as

[TABLE]

where

[TABLE]

Here we used the fused $R$- and $K$-matrices defined (84),(85) and (94), (95) respectively. Note that the traces in (105) and (106) are now over a $d$-dimensional space with $d=6$.

Next, we need to specify the $K$-matrices $K^{\pm}_{a}$ and $K^{\pm}_{\bar{a}}$ entering the definition of the boundary transfer matrices involved in the $T$-system which we intend to derive. First, we restrict to invertible solutions to the reflection equations (34) and (35). Next, in order to construct a good set of fused transfer matrices, we need to start with two operators $\tau(\lambda)$ and $\bar{\tau}(\mu)$ such that $[\tau(\lambda),\bar{\tau}(\mu)]=0$. This enforces a relation between the $K$-matrices $K^{\pm}_{a}$ and $K^{\pm}_{\bar{a}}$ of $\tau(\lambda)$ and $\bar{\tau}(\mu)$ respectively. In fact, it was shown in [43] that, in order for the transfer matrices $\tau(\lambda)$ and $\bar{\tau}(\mu)$ to commute, the $K$-matrix $K_{\bar{a}}$ has to be related to $K_{a}$ as

[TABLE]

where $\alpha$ is a non-vanishing constant number. In the following, we will choose $\alpha=1$, and derive the corresponding $T$-system. Explicitly, we will consider boundary transfer matrices $\tau(\lambda)$ and $\bar{\tau(\lambda)}$ as in (28) and (29) with the $K$-matrices

[TABLE]

where $K_{s}$ is any invertible symmetric numerical matrix. We can now plug the resolution of the identity

[TABLE]

into (104), and setting $\mu=\nu-i/2$, $\lambda=\nu+i/2$ we obtain

[TABLE]

From this expression one can see that the product $\tau_{b}(\nu-i/2)\tau_{a}(\nu+i/2)$ decomposes into the sum of two terms. In particular, it is straightforward to see

[TABLE]

In order to conclude our calculation, we show

[TABLE]

Making use of (81) and (82), we get immediately

[TABLE]

where

[TABLE]

and

[TABLE]

One can straightforwardly compute

[TABLE]

where the matrices $R^{(2)}_{\langle ab\rangle,c}(\lambda)$ and $\bar{R}^{(2)}_{\langle ab\rangle,c}(\lambda)$ are defined in (84), (85), while

[TABLE]

Note that the following basis has been chosen for the fused spaced obtained out of $h_{a}{\otimes}h_{a}$

[TABLE]

Analogously, for any $K$-matrix which is constant, invertible and symmetric, one can compute

[TABLE]

Recalling Eq. (111) we can collect all the prefactors and obtain the final result of the calculation. Performing then analogous steps for the product $\bar{\tau}(\lambda-i/2)\bar{\tau}(\lambda+i/2)$ we get

[TABLE]

where

[TABLE]

From this result, we can finally exploit the functional relations (66) and (67) to obtain higher functions $\Phi^{(1)}_{n}(\lambda)$ and $\Phi^{(2)}_{n}(\lambda)$.

3.6 The final result

We can now summarize the final result of our calculations, namely the explicit $T$-system for soliton-non-preserving open boundary conditions. It is convenient to write down the latter for the transfer matrices $\tau^{s}_{N}(\lambda)$ with shifted spectral parameters defined in (49). After performing the shift (48) in the formulas of the previous section, we find that the fused transfer matrices obtained from $\tau_{N}(\lambda)$ satisfy the $T$-system

[TABLE]

with the convention

[TABLE]

and where $\tau^{s}_{N}(u)$ and $\bar{\tau}^{s}_{N}(u)$ are defined in (49) and (54). Here the functions $\tilde{\Phi}_{m}^{(a)}(\lambda)$ are given by

[TABLE]

where

[TABLE]

Note that the general expressions (144) and (145) are obtained exploiting the relations (66) and (67), and using the convention $\tilde{\Phi}^{(1)}_{0}\equiv\tilde{\Phi}^{(2)}_{0}\equiv 1$. This result generalizes the findings of [49] to $SU(3)$-invariant boundary transfer matrices in the case of soliton-non-preserving boundary conditions, and will be used throughout the rest of this work.

4 The $Y$-system and the Loschmidt echo

As we have anticipated, the usefulness of the $T$-system in (139), (140) relies on the possibility of solving it for any length $N$, and also in the limit $N\to\infty$. In turn, this provides an explicit solution to the problem of computing the spectrum of the boundary transfer matrix $\tau^{s}_{N}(\lambda)$. In practice, this is done as follows. First, one needs to write down the $Y$-system corresponding to the $T$-system; then the latter is cast, within a standard procedure, into a set of non-linear integral equations, which can be finally solved numerically. Since the $T$-system in the case of boundary transfer matrices is the same of the periodic case, also the $Y$-system [50] will have the exact same form, encoded in Eqs. (71) and (72). Then, we are left with the task of explicitly writing down the corresponding non-linear integral equations. This procedure is well-known in the literature, and has been repeatedly used also in previous applications of the QTM formalism to quantum quenches [36, 51]. For this reason, here we only briefly summarize it, referring to the work [36] for a detailed explanation of the method.

In summary, one needs to take the logarithmic derivative of both sides of (71) and (72), and then perform a Fourier transform to obtain integrals along segments with constant imaginary parts. The integration contours can then be shifted back to the real axis, by carefully taking into account all the singularities of the logarithmic derivatives within the physical strip, which is the portion of the complex plane defined as

[TABLE]

where ${\rm Im}(\lambda)$ denotes the imaginary part of the complex number $\lambda$. From the above discussion, it is obvious that the final form of the integral equations will essentially depend on the analytic structure of the $Y$-functions inside the physical strip. Different boundary matrices $K(\lambda)$ and different inhomogeneities $\xi_{j}$ will give rise to different singularities, and ultimately to different integral equations (in fact, differences will also arise in the asymptotic behavior of $y^{(r)}_{n}(\lambda)$ as $n\to\infty$, cf. Ref. [36]).

In order to illustrate the above procedure in detail, we consider the problem of computing the leading eigenvalue of the boundary transfer matrix associated to the delta-state (43) (with inhomogeneities (52) and (53)). As we will see in Sec. 4.3, this will immediately give us access to the Loschmidt echo at imaginary and real times after the quench from the delta-state. Once again, we stress that, even though we focus on a particular boundary transfer matrix corresponding to a given boundary state, our treatment is quite general, and analogous calculations can be carried out in more general cases.

4.1 The non-linear integral equations for the delta-state

The delta-state (43) is associated with the boundary transfer matrix (49) with $K$-matrices $K^{-}V=K^{+}V=\mathbf{1}$. The inhomogeneties are chosen as (52) and (53), which corresponds to imaginary time evolution. At the end of the calculation, we will address the real-time case.

First, we note that in this case the boundary transfer matrix commutes with the global spin operator

[TABLE]

and that the leading eigenvalue of $\tau^{s}_{N}(0)$ lies in the $S^{z}=0$ sector. In order to get an intuition of the corresponding $Y$-functions $y^{(1)}_{1}(\lambda)$ and $y^{(2)}_{1}(\lambda)$, we compute it numerically in the physical strip (148) for small values of $N$. Remarkably, this can be done even without the knowledge of the Bethe roots as follows. First, one computes the leading eigenvector of $\tau^{s}_{N}(0)$ using exact diagonalization calculations for a finite $N$. Next, since all of the transfer matrices introduced in the previous section commute with one another, one can act on the latter with $\tau^{s}_{N}(\lambda)$ and $\bar{\tau}^{s}_{N}(\lambda)$ to get the associated eigenvalues as a function of $\lambda$. Then, the $Y$-functions $y^{(1)}_{1}(\lambda)$ and $y^{(2)}_{1}(\lambda)$ are obtained straightforwardly from (69), (70). Note that a crucial ingredient for this is the knowledge of the functions $\Phi_{1}^{(1)}(\lambda)$ and $\Phi_{1}^{(2)}(\lambda)$ obtained from the explicit form of the $T$-system.

Having constructed the $Y$-functions $y^{(1)}_{1}(\lambda)$ and $y^{(2)}_{1}(\lambda)$ for finite $N$, we can obtain numerically all the higher $Y$-functions through the $Y$-system, and study the corresponding analytic behavior inside the physical strip (148). Analytical inspection for finite values of $N$ reveals the following behavior

• •

$y^{(1)}_{1}(\lambda)$ and $y^{(2)}_{1}(\lambda)$ have a double pole at $\lambda=0$ for $j$ odd, and a double zero at $\lambda=0$ for $j$ even;

• •

in addition, $y^{(1)}_{1}(\lambda)$ has a pair of order $N$ zeros at $\lambda=\pm i\left(\frac{1}{2}+\frac{\beta}{N}\right)$ and a pair of order $N$ poles at $\lambda=\pm i\left(\frac{1}{2}-\frac{\beta}{N}\right)$, where we recall that $\beta$ enters the definition of the transfer matrices through the inhomogeneities (52),(53).

The above analytical structure has been verified numerically up to $N=8$. In the following we assume it holds for all the values of $N$, as can be checked a posteriori from the explicit knowledge of the $Y$ functions, determined as a solution of a set of non-linear integral equations. There, the poles and zeroes are determined in a self-consistent manner, and we indeed check that the analytic structure inside and in the vicinity of the physical strip is stable as $N$ increases.

Next, we note that the following asymptotic behavior of the $Y$-functions can be deduced

[TABLE]

where $C_{m}$ are determined from the recursive relations:

[TABLE]

with the additional constraint $C_{0}=0$, and $C_{1}=2$. The recursion relation (151) can be solved explicitly, yielding

[TABLE]

Using the properties of the $Y$-functions stated above, it is straightforward to perform the procedure sketched in the previous section (see [36] for more details). In particular, we arrive at

[TABLE]

and, for $m>1$,

[TABLE]

where

[TABLE]

Here we introduced the convolution between two functions

[TABLE]

The above integral equations can be solved numerically to high precision. In order to do this, the infinite system has to be truncated to a maximum number of equations $n_{\rm M}$. The truncation can be performed efficiently by imposing the boundary condition

[TABLE]

which automatically enforces the asymptotic behavior (152).

Once the $Y$-functions are computed, the eigenvalues of $\tau^{s}_{N}(\lambda)$ and $\bar{\tau}^{s}_{N}(\lambda)$ can be obtained from the definition of the $Y$-functions

[TABLE]

where

[TABLE]

In fact, in the following we will be interested in the leading eigenvalue of $\mathcal{T}$, defined in (47). Hence, it is convenient to consider the normalized eigenvalues $\Lambda(\lambda)$ and $\bar{\Lambda}(\lambda)$

[TABLE]

The function $\Lambda(\lambda)$, by construction, is such that $(1/3)\Lambda(0)$ is the leading eigenvalue of $\mathcal{T}$ (the prefactor $1/3$ takes into account the normalization of the initial delta-state). It is now straightforward to rewrite

[TABLE]

We now make the assumption that ${\Lambda}(\lambda)$ and $\bar{\Lambda}(\lambda)$ have no poles inside the physical strip, which we verified numerically for sizes up to $N=8$ and can once again be checked by examining the numerical solutions the non-linear integral equations. Introducing

[TABLE]

we finally arrive at the non-linear integral equation

[TABLE]

As a final step, we modify the integral equations above in such a way that $\ln\Lambda(\lambda)$ and $\ln\bar{\Lambda}(\lambda)$ are explicitly expressed in terms of $y^{(r)}_{1}$ with $r=1,2$. Following [101], we define

[TABLE]

so that we can rewrite the equations as

[TABLE]

where $r=1,2$ while

[TABLE]

and

[TABLE]

Here we introduced the operator $s^{-1}$ which is understood to act in Fourier space as multiplication by $\hat{s}^{-1}(k)$ [where $s(\lambda)$ is defined in (157), and where $\hat{s}(k)$ denotes the Fourier transform of $s(\lambda)$]. A formal solution to Eq. (173) can be written down explicitly (see [47]). In particular, one has

[TABLE]

where

[TABLE]

Writing down everything explicitly we obtain the final result

[TABLE]

The results presented so far give us access to the leading eigenvalue of the boundary transfer matrix $\tau^{s}_{N}(\lambda)$ for arbitrary length $N$ of the chain. However, it is also possible to compute the limit $N\to\infty$: this can be simply done by taking the limit of the $N$-dependent driving terms, which can be easily performed. Explicitly, the relevant equations for the leading eigenvalue of $\Lambda(\lambda)$ in the limit $N\to\infty$ can be summarized as follows

[TABLE]

with the convention $y^{(r)}_{0}(\lambda)\equiv 0$ and

[TABLE]

where

[TABLE]

In Sec. 4.3 we will present an explicit numerical solution to the above equations. Indeed, as we will see, the value of $\Lambda(0)$ is closely related to the Loschmidt echo in the thermodynamic limit. In the next subsection, instead, we will show how the solution of the integral equations above can be obtained analytically in the limit $\beta\to 0$.

4.2 The analytic solution to the $Y$-system at vanishing times

From the discussion of the previous section, it follows that in order to determine the solution to the $Y$-system, one needs in general to resort to a numerical evaluation of the corresponding integral equations. This is true, for example, for boundary transfer matrices associated to a given integrable state, in the case of generic values of the parameter $\beta$ parametrizing the inhomogeneities [cf. Eq. (53)]. However, it is possible to derive an analytic solution to the $Y$-system in the limit $\beta\to 0$. This result is particularly important due to its connection to the problem of deriving the quasi-momentum distribution functions of the quasi-particles in the post-quench steady state [33]. In this section, we present such an analytical solution for arbitrary integrable states.

Following [36], we start by introducing the operators

[TABLE]

where $T_{a}(\lambda)$ and $T_{\bar{a}}(\lambda)$ are defined in (50) and (55) respectively, while $|\psi_{0}\rangle$ is given in (39). Here we also defined

[TABLE]

where $\tilde{k}_{ij}$ are the entries of the inverse (symmetric) matrix $K_{s}^{-1}$, namely

[TABLE]

Note that one has the identification

[TABLE]

where $\mathcal{N}=\langle\psi_{0}|\psi_{0}\rangle$ and $\widetilde{\mathcal{N}}=\langle\tilde{\psi}_{0}|\tilde{\psi}_{0}\rangle$, and where $\tau^{s}_{N}(\lambda)$, $\bar{\tau}^{s}_{N}(\lambda)$ are defined in (49), (54). Now, using simple manipulations, it is not difficult to show

[TABLE]

where we defined

[TABLE]

and analogously for $|\tilde{\Psi}_{0}^{\pm}\rangle$ and $|\tilde{\psi}\rangle_{\mathcal{N}}$. Here we introduced the following transfer matrices $t_{\lambda,-\lambda}(w)$ acting in the original time direction

[TABLE]

It follows that

[TABLE]

Let $\Lambda_{0}(\lambda)$ and $\widetilde{\Lambda}_{0}(\lambda)$ be the eigenvalues of $\mathcal{U}(\lambda)$, $\widetilde{\mathcal{U}}(\lambda)$ such that they are the leading one for $\lambda=0$. Then, for small values of $\lambda$ and large system sizes $L$, we have

[TABLE]

On the other hand, the following inversion relations hold for small values of $\lambda$

[TABLE]

The rigorous meaning of (201) and (202) and their proof is reported in A. Putting everything together, we obtain

[TABLE]

We now have all the ingredients to compute the analytic solution to the $y$-system in the limit $\beta\to 0$. First, it follows from the definitions (69), (70) that

[TABLE]

Exploiting the relations (189), (190), between $\mathcal{U}(\lambda)$, $\widetilde{\mathcal{U}}(\lambda)$ and $\tau^{s}_{N}(\lambda)$, $\bar{\tau}^{s}_{N}(\lambda)$, using the results (144)–(147), and collecting all the factors, we arrive at the final result, which is extremely simple and reads

[TABLE]

Note that these expressions do not depend on $N$, and that higher $Y$-functions can be obtained using the $Y$-system relations. In the case of diagonal reflection matrices $K_{s}$, the above expressions can be explicitly written as

[TABLE]

4.3 The Loschmidt echo

In this section we finally present an immediate physical application of the results derived above, namely the computation of the Loschmidt echo after a quench from integrable states. While we present numerical evaluation for the case of the delta-state (43), we stress that the treatment is general, and analogous calculations hold for arbitrary integrable states.

We recall that the Loschmidt echo is defined by

[TABLE]

and gives us information about the probability of finding the system at the initial configuration after a time $t$ after the quench. For global quenches, the latter vanishes exponentially with the system size, so it is convenient to define the return probability

[TABLE]

The Loschmidt echo is simply related to the partition function (44) by $\mathcal{L}(t)=|\mathcal{Z}(it)|^{2}$. Furthermore, from (46) one can see that the return probability (211) can be obtained from the leading eigenvalue $\Lambda(0)$ of the transfer matrix $\mathcal{T}$ defined in (47). In the case of the delta-state (43), we can use directly the results of the previous section and obtain

[TABLE]

where the term $\log(3)$ comes from the normalization of the delta-state. The leading eigenvalue $\Lambda(0)$ can be computed using the formula (183), where the $Y$-functions are obtained by solving the integral equations (181), (182). We evaluated numerically these formulas, and reported our results in Fig.1. In particular, the return probability (212) at imaginary times is displayed in subfigure $(a)$ of Fig. 1. We see that the behavior is the expected one: the curve is monotonic and for large values of $\beta$ we have a linear growth. The slope of the curve for large $\beta$ can in principle be computed along the lines of [37], and is expected to be expressed in terms of both the ground-state energy of the Hamiltonian (1) and the overlap between the ground state and the initial state.

Next, following the same steps outlined in the previous section, it is also possible to compute the return probability for real times. This problem was tackled for the spin-$1/2$ XXZ chain in [51], where it was shown that several complications arise. In particular, in this case the spectrum of the boundary transfer matrix $\tau^{s}_{N}(\lambda)$ might display level crossings, which are responsible for the emergence of non-analytic points in the real-time dynamics of the return probability. In order to locate these crossings, one may examine the finite-$N$ spectrum of the QTM from exact diagonalization. In the present case, we checked that the first crossing occurs at a time $t\sim 1.5$, which guarantees in particular that the return probability remains analytic in the time window considered in Figure 1 (b). This approach has to be compared with that of [63], where the QTM is diagonalized numerically directly in the $N\to\infty$ limit. Furthermore, following the dynamics of a single eigenvalue, one finds that the analytic properties of the corresponding $Y$-functions might change in time, so that one is forced to employ an “excited-state” thermodynamic Bethe ansatz approach [102, 103, 104, 105]. We mention that while non-analytic points in the Loschmidt echo have been extensively studied in the context of dynamical quantum phase transitions [53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78] the QTM approach provided one of the few analytical studies of the latter in the presence of interactions [51].

It follows from the above discussion that for real times the expression (212) (after substitution $i\beta\to t$) will be correct only up to a time $t^{\ast}$, after which the correct formula must be modified using some “excited-state” thermodynamic Bethe ansatz techniques [102, 103, 104, 105]. In the $SU(3)$-invariant spin chain, this can be done in the same way explained in [51] for the $XXZ$ spin-$1/2$ model. The implementation of the method goes beyond the scope of the present paper; here we simply evaluated the real-time return probability for small times, for which Eq. (212) can be used. The result is plotted in subfigure $(b)$ of Fig. 1. We explicitly checked the latter against numerical calculations based on exact diagonalization for small system sizes $L$, finding perfect agreement.

5 Conclusions

This is the second one of two papers devoted to the study of quantum quenches in the nested $SU(3)$-invariant spin chain (1). In this and in the previous work [33], we have generalized the QTM approach, initially developed for XXZ Heisenberg chains [36, 38], to quantum quenches in nested systems. In particular, we have shown how it is possible to define a class of product integrable states, and how to compute exactly several quantities after the quench.

In this technical work we have focused on the derivation of some fusion relations of integrable boundary transfer matrices, which were required in order to derive the results of [33], and which are an important ingredient in the application of the QTM approach to quantum quenches. In particular, in this work we have explicitly exploited them to obtain exact formulas for the real- and imaginary-time Loschmidt echo after a quench from integrable states. We believe that the technical tools derived in this paper are interesting per se, and might find applications in other problems, possibly related to equilibrium physics in the presence of open boundary conditions.

This work, together with [33], provides a comprehensive picture for quantum quenches from integrable states, complementing the works [36, 38] where this formalism was initially developed. Furthermore, the successful application of these methods to the complicated nested models shows that integrable states are indeed ideal probes to explore quite generally the quench-dynamics in integrable systems. Finally, these works also open several directions to be investigated. Among these, one of the most fascinating questions pertains the study of quantum quenches in nested, multi-component Fermi and Bose gases, which are of great relevance for cold-atom realizations [81, 82, 83].

Acknowledgments

PC acknowledges support from ERC under Consolidator grant number 771536 (NEMO). EV acknowledges support by the EPSRC under grant EP/N01930X. BP was supported by the BME-Nanotechnology FIKP grant of EMMI (BME FIKP-NAT), by the National Research Development and Innovation Office (NKFIH) (K-2016 grant no. 119204 and KH-17 grant no. 125567), and by the “Premium” Postdoctoral Program of the Hungarian Academy of Sciences. Part of this work has been carried out during a visit of BP and EV to SISSA whose hospitality is kindly acknowledged.

Appendix A Inversion relations of transfer matrices

In this section we will prove

[TABLE]

where we introduced the Hilbert-Schmidt norm of an operator as

[TABLE]

and $\alpha,\beta>0$ are real numbers. We follow analogous calculations carried out for example in [106]. We derive in particular (213), as (214) can be established completely analogously.

We define $L(x)$, $M(x)$ via

[TABLE]

where we introduced the matrices $E^{a,b}$ whose matrix elements are given by

[TABLE]

From the definitions, we have

[TABLE]

so that it is easy to compute

[TABLE]

where

[TABLE]

The leading eigenvalue $\Gamma_{0}(\lambda)$ of $\left[\sum_{d}X_{dd}(\lambda)\right]\left[\sum_{d}X_{dd}(-\lambda)\right]$ can be easily found to be

[TABLE]

Furthermore, it can be argued that $\Gamma_{0}(\lambda)$ has a finite gap with respect to the next to leading eigenvalues in a neighborhood of $\lambda=0$. Finally, we can compute

[TABLE]

and so

[TABLE]

where

[TABLE]

Similar to the previous case, we find that the leading eigenvalue $\Omega_{0}(\lambda)$ of

[TABLE]

has a finite gap with respect to the next to leading ones, in a neighborhood of $\lambda=0$ and it reads

[TABLE]

Putting everything together, it is straightforward to compute

[TABLE]

with $\alpha>0$ in a neighborhood of $\lambda=0$. This proves (213).

The other relation (214) can be proven in a similar way, by making use of the local inversion relation (9).