Anomalous relaxation and multiply time scales in the quantum $XY$ model with boundary dissipation

TL;DR

This paper exactly analyzes the relaxation dynamics of the quantum XY model with boundary dissipation, revealing multiple time scales, counter-intuitive dissipation effects, and odd-even phenomena, providing new insights into topological qubits in environments.

Contribution

It provides an exact calculation of the relaxation dynamics in the quantum XY model with boundary dissipation, uncovering multiple time scales and counter-intuitive dissipation effects.

Findings

Relaxation characterized by multiple time scales.

Weak dissipation induces relaxation, strong dissipation can suppress it.

Odd-even effects in long chain relaxation phenomena.

Abstract

The relaxation of many-body system is still a challenging problem that has not well been understood. In this work we exactly calculate the dynamics of the quantum $XY$ model with boundary dissipation, in which the density matrix in terms of Majorana operators can be decoupled into independent subspaces represented by different number of Majorana fermions. The relaxation is characterized by multiply time scales, and in the long-time limit it is determined by the single particle relaxation process in a typical time scale $T^*$. For the bulk bands, we find $T^* \propto N^3/\gamma n^2 $ in the weak dissipation limit; and $T^* \propto \gamma N^3/ n^2$ in the strong dissipation limit, where $N$ is the chain length, $\gamma$ is the dissipation rate and $n$ is the band index. For the edge modes $T^* \propto 1/\gamma$, indicating of most vulnerable to dissipation in the long chain limit. These results are counter-intuitive because it means any weak dissipation can induce relaxation, while strong dissipation can induce weak relaxation. We find that these two limits correspond to two different physics, which are explained based on the first and second-order perturbation theory in an equivalent non-Hermitian model.Furthermore, we show that even in the long chain limit the relaxation may exhibit strong odd-even effect. These results shade new insight into the dynamics of topological qubits in environment.

Figures

Equations

$$\dot{\rho }=\mathcal{L}(\rho )=-i[H_{\text{XY}},\rho ]+\mathcal{D}(\rho ).$$

$$H=i{g}_1\sum _{j=1}^{N-1}{\alpha }_{2j}{\alpha }_{2j+1}-i{g}_2\sum _{j=1}^{N-1}{\alpha }_{2j-1}{\alpha }_{2j+2}.$$

$$\mathcal{D}(\rho )=-\gamma \sum _{j=1,N}({\alpha }_{2j-1}{\alpha }_{2j}\rho {\alpha }_{2j-1}{\alpha }_{2j}+\rho ).$$

$$E=\sum _i{n}_i{ϵ}_i,n_i=\{0,1\},$$

$$H=\frac{i}{2}\sum _{k=1}^{N-1}{ϵ}_k{b}_k^{\prime }{b}_k^{\prime \prime }+\frac{i\delta E_c}{2}{\alpha }_{\text{L}}{\alpha }_{\text{R}},$$

$$|1⟩=\frac{1+{\alpha }_{\text{L}}}{\sqrt{2}}|0⟩,|2⟩=\frac{1-{\alpha }_{\text{L}}}{\sqrt{2}}|0⟩,$$

$$\bar{\rho }=\underset{t\to \infty }{\lim }\rho (t)=(I+cP)/2^N,$$

$$\rho =\frac{1}{2^N}\sum c_{a_1,a_2,\cdots ,a_{2N}}{\alpha }_1^{a_1}{\alpha }_2^{a_2}\cdots {\alpha }_{2N}^{a_{2N}},$$

$$\mathcal{K}={\mathcal{K}}_0\oplus {\mathcal{K}}_1\oplus {\mathcal{K}}_2\oplus \cdots \oplus {\mathcal{K}}_{2N},$$

$$\text{Tr}(\rho (t)|n⟩⟨n|)$$

$$i{\partial }_t\Psi =2(H_0+i\Gamma )\Psi ,$$

$$H_{0}=\left(\begin{array}[]{ccccc}0&ig_{2}&0&\cdots&0\\ -ig_{2}&0&ig_{1}&\cdots&0\\ 0&-ig_{1}&0&\ddots&0\\ \vdots&\vdots&\ddots&\ddots&ig_{2}\\ 0&0&0&-ig_{2}&0.\end{array}\right).$$

$${\Delta }_N=\frac{(g_1{g}_2{)}^{m-1}}{\sin \theta }[g_1{g}_2(2i\gamma +ϵ)\sin (m+1)\theta$$

$$+(-{\gamma }^2ϵ+i{g}_1^2\gamma +i{g}_2^2\gamma )\sin (m\theta )]=0,$$

$${\Delta }_N=\frac{(g_1{g}_2{)}^{m-1}}{\sin \theta }[(-{\gamma }^2+g_1^2+i2\gamma ϵ)\sin (m\theta )$$

$$+g_1{g}_2\sin (m+1)\theta -{\gamma }^2\frac{g_2}{g_1}\sin (m-1)\theta ]=0,$$

$${ϵ}^2=g_1^2+g_2^2+2g_1{g}_2\cos \theta .$$

$${ϵ}_n={ϵ}_{n,\text{r}}-i{ϵ}_{n,\text{i}},{\theta }_n=\frac{n\pi }{m+1}+z_{n,\text{r}}+i{z}_{n,\text{i}},$$

$${ϵ}_{n,\text{i}}^{\text{o}}=\frac{(g_1^2+g_2^2)n^2{\pi }^2\gamma }{(g_1+g_2{)}^2{m}^3},{ϵ}_{n,\text{i}}^{\text{e}}=\frac{2g_2^2{n}^2{\pi }^2\gamma }{(g_1+g_2{)}^2{m}^3}.$$

$$T^{\ast }=\text{max}(1/(2{ϵ}_{n,\text{i}})).$$

$${ϵ}_{n,\text{i}}=\frac{2g_1^2{g}_2^2{n}^2{\pi }^2}{(g_1+g_2{)}^2{m}^3\gamma }.$$

$${ϵ}_{\text{L}}={ϵ}_{\text{R}}=-i\gamma (1-{\lambda }^2),$$

$${ϵ}_{\text{L}}=-i\gamma +i{g}_2^2/\gamma ,{ϵ}_{\text{R}}=-i\gamma +i{g}_1^2/\gamma ,$$

$${ϵ}_{\text{L}}={ϵ}_{\text{R}}=-i\gamma +i{g}_2^2/\gamma ,$$

$$V=\left(\begin{array}[]{ccccc}0&ig_{2}&0&\cdots&0\\ -ig_{2}&0&0&\cdots&0\\ 0&0&0&\ddots&0\\ \vdots&\vdots&\ddots&\ddots&ig_{2}\\ 0&0&0&-ig_{2}&0.\end{array}\right),$$

$${ϵ}_{n,\text{i}}\propto -\sum _{j=a,b}\frac{⟨{\tilde{\psi }}_N^n|V|j⟩⟨j|V|{\tilde{\psi }}_N^n⟩}{(E_n-E_j)},$$

$$\ln T^{\ast }=\alpha \ln N+{\alpha }^{\prime },$$

Full text

While the dynamics of qubits in environment has been well studied Gardiner and Haken (1991); Johansson et al. (2012), the same issue in the many-body systems is still one major challenge in theory Weimer (2015); Le Boité et al. (2013); Karevski et al. (2013); Prosen (2014, 2011); Rowlands and Lamacraft (2018); Banchi et al. (2017); Medvedyeva et al. (2016); Ribeiro and Prosen (2019); Cai and Barthel (2013); Žnidarič (2015, 2011, 2010) due to more expensive computation cost com . However, this is an important question at least from two diverse aspects. The many-body systems may possess some features that are totally different from the single particle systems, such as ergodicity and thermalization Eisert et al. (2015); Nandkishore and Huse (2015), which are fundamental concepts in statistics. In the trapped ions, it may exhibit different dynamics depending strongly on the initial states, which are explained based on quantum many-body scar Bernien et al. (2017); Turner et al. (2018a, b). Moreover, it is also an important issue in topological quantum computation Nayak et al. (2008); Freedman et al. (2003); Sau et al. (2010), in which the two ground states are separated from the excited bands by a finite energy gap Zhang et al. (2018); Nadj-Perge et al. (2014); Wang et al. (2018); Sun et al. (2016); Xu et al. (2015); Liu et al. (2018); Mourik et al. (2012); Deng et al. (2012); Das et al. (2012). Thus if the temperature is much lower than the excitation gap, the occupation of the excited states are exponentially small. This picture is not necessarily true in the presence of dissipation, which can induce direct coupling between ground states and excited states.

Here we explore the roles of edge modes and bulk bands in the dynamics of the quantum $XY$ model with boundary dissipation, in which the two ground states are protected by a finite gap. (I) In terms of Majorana operators, the density matrix is decoupled into different subspaces represented by different number of Majorana fermions. In time evolution the density matrix exhibits multiply relaxation scales, in which the slowest decay is given by the single particle relaxation. This time scale $T^{*}$ is used to determine the relaxation time of the many-body ground state with dissipation. (II) In the weak dissipation limit, $T^{*}\propto N^{3}/\gamma$, where $N$ is the total chain length and $\gamma$ is the boundary dissipation rate. However, in the strong dissipation limit, $T^{*}\propto\gamma N^{3}$. The edge modes are shown to be most vulnerable to dissipation due to $T^{*}\propto 1/\gamma$. These results are counter-intuitive because it means that weak dissipation can induce fast relaxation, while strong dissipation can induce weak relaxation. We understand these results by mapping the single particle dynamics to a non-Hermitian model. (III) This dynamics exhibits strong odd-even effect in the weak dissipation limit, which reduces to a unified form in the strong dissipation limit. These results shade new insight into the dissipation and relaxation of the topological qubits in environment.

Model and Master Equation. We consider a quantum $XY$ model with boundary dissipation, which reads as

[TABLE]

Here $H_{\text{XY}}=-\sum_{i=1}^{N-1}(g_{1}\sigma_{i}^{x}\sigma_{i+1}^{x}+g_{2}\sigma_{i}^{y}\sigma_{i+1}^{y})$ Osterloh et al. (2002); Zhu (2006); Eisert and Osborne (2006); Osborne and Nielsen (2002); Cheng et al. (2017); Osborne and Nielsen (2002) and the dissipation by the Lindblad operator is Hein et al. (2005); Cai et al. (2006) $\mathcal{D}(\rho)=\frac{\gamma}{2}\sum_{j=1,N}(2\sigma_{j}^{z}\rho\sigma_{j}^{z}-2\rho)$. This model is fermionized using Majorana operators via Jordan-Wigner transformation as $\alpha_{2j-1}=(\prod_{k}^{j-1}\sigma_{k}^{z})\sigma_{j}^{x},\alpha_{2j}=(\prod_{k}^{j-1}\sigma_{k}^{z})\sigma_{j}^{y}$, after which Pascual Jordan (1928); Sachdev (2011); Kitaev (2001)

[TABLE]

In this case, the Lindblad operator is still made by local dissipation dis

[TABLE]

The spectra between the $XY$ model and $p$-wave superconducting model are related by

[TABLE]

thus the $N$ eigenvalues from the single particle Hamiltonian can be used to construct all the $2^{N}$ eigenvalues in the quantum XY model (see Fig. 1 (a) and (c)). In the fermion representation the two localized zero modes at the open ends give rise to the two-fold degeneracy in the $XY$ model. We focus on $g_{1}>0$ and $g_{2}>0$, and the phase diagram for the $XY$ model and its corresponding single fermion phases are given in Fig. 1 (b), with boundary at $\lambda=g_{2}/g_{1}=1$. In fermion representation it can be regarded as two separate Majorana chains $A$($\alpha_{1},\alpha_{4},\alpha_{5},\cdots$) and $B$($\alpha_{2},\alpha_{3},\alpha_{6},\cdots$), as shown in Fig. 1 (d). We can bring Eq. 2 to a form of paired Majorana operators Kitaev (2001); Kitaev and Laumann (2009)

[TABLE]

where $\delta E_{c}\sim e^{-N/\xi}[(-1)^{N}+1]$. After a special orthogonal transformation the Hamiltonian in Eq. 2 can be brought into small Jordan blocks, with $b_{k}^{\prime}$ and $b_{k}^{\prime\prime}$ are some new Majorana operators, following the prominent work by Kitaev Kitaev (2001). By a transformation from Majorana fermion to canonical fermion $c_{k}=(b_{k}^{\prime}+ib_{k}^{\prime\prime})/2$ and $c_{k}^{\dagger}=(b_{k}^{\prime}-ib_{k}^{\prime\prime})/2$, the spectra of the bulk bands can be obtained as $\pm\epsilon_{k}$. In the last term $\alpha_{\text{L}}$ and $\alpha_{\text{R}}$ are edge modes at the left and right edges, which can be written as $\alpha_{\text{L}}=\alpha_{1}-\lambda\alpha_{5}+\lambda^{2}\alpha_{9}-\cdots$ and $\alpha_{\text{R}}=\alpha_{2N}-\lambda\alpha_{2N-4}+\lambda^{2}\alpha_{2N-8}-\cdots$ (for $\lambda<1$). In $\delta E_{c}$, $\xi\propto 1/|\ln g_{1}-\ln g_{2}|$ defines the correlation length Kitaev (2001). This odd-even effect is a typical feature of coupling between two distant zero modes, which may happen even in continuous space Cheng et al. (2009). With these operators, the two ground states in Fig. 1 (a) can be written as

[TABLE]

where $|0\rangle$ is the ground state satisfying $c_{k}|0\rangle=0$ for all $k$ (thus $|0\rangle=\prod_{k}c_{k}|\text{vac}\rangle$, where $|\text{vac}\rangle$ is the vaccum state).

This dynamics respects the parity symmetry $[\mathcal{L},P]=0$ with $P=\prod_{j}^{N}\sigma_{j}^{z}=i^{N}\prod_{j}^{2N}\alpha_{j}$. In the long time limit,

[TABLE]

where $c=\langle\Psi_{0}|P|\Psi_{0}\rangle$. This corresponds to the maximal mixed state at infinite high temperature.

Evolution of the density matrix. The density matrix can be written as Kronecker product of Pauli matrices. However, it is more convenient to write this matrix in terms of Majorana operators as following Calabrese and Cardy (2005); Vidal et al. (2003); Eisert et al. (2010),

[TABLE]

where $a_{j}=\{0,1\}$. Thus the dynamics of $\rho$ is decoupled into different subspaces denoted as $\mathcal{K}_{n}$ for $n=0-2N$, where $n$ is the number of Majorana operators given by $n=\sum_{i}a_{i}$. Formally, we have

[TABLE]

where the dimension of $\mathcal{K}_{i}$ is $C_{2N}^{i}$. We can readily check that the dimension of $\mathcal{K}$ (the whole space of the Hamiltonian) is $C_{2N}^{0}+C_{2N}^{1}+C_{2N}^{2}+\cdots+C_{2N}^{2N}=4^{N}$. This decoupling is essentially the same as the probability distribution in the classical Ising model by Glauber Glauber (1963). By a direct comparison with Eq. 7, we can find that while the first term ($a_{j}\equiv 0$) and the last term ($a_{j}\equiv 1$) are unchanged, which give Eq. 7; while all the other terms will disappear in the long time limit.

In principle, the dynamics of $\rho$ can be calculated using hierarchy equations Glauber (1963). We may define the following variables $\psi_{i}=\text{Tr}(\rho\alpha_{i})$, $\psi_{ij}=\text{Tr}(\rho\alpha_{i}\alpha_{j})$, $\psi_{ijk}=\text{Tr}(\rho\alpha_{i}\alpha_{j}\alpha_{k})$ (with $i\neq j\neq k$) and then calculate their time evolution based on the Heisenberg equations. We find that for the boundary dissipation in Eq. 1, the dynamics of these variables are restricted to their own subspaces $\mathcal{K}_{1,2,3}$. This limit will show up the multiply time scales during dynamics much more clearly. Obviously, this approach can be generalized to models with more complicated dissipation and many-body interaction.

We first compare the full calculation of $\rho$ against the dynamics in subspace $\mathcal{K}_{1}$. In Fig. 2 (a), we present the real part of the eigenvalues of superoperator $\mathcal{L}$ of Eq. 1, in which the two zero eigenvalues correspond to the unchanged state in $\mathcal{K}_{0}$ and $\mathcal{K}_{2N}$ subspaces. We also calculate the corresponding eigenvalues of this superoperator in the subspace constructed by $\mathcal{K}_{1}$, which is given in Fig. 2 (b). We find that the smallest eigenvalue Re($\lambda_{1\sigma}$) of the $\mathcal{L}$ in the $XY$ model is the same as the spectra in $\mathcal{K}_{1}$ subspace, indicating that in the long time limit, all the higher-order terms in $\mathcal{K}_{k\geq 2}$ subspaces decay much faster than that in $\mathcal{K}_{1}$, leaving $\mathcal{K}_{1}$ to be the dominated relaxation channel for the quantum $XY$ model. In Fig. 2 (c) and (d), we show the projection of $\rho(t)$ to the ground states of the $XY$ model and the edge modes. In (c), it will approaches $1/2^{N}$ ($c=0$); while in (d), it will approaches zero, as expected.

To reinforce this conclusion, we also calculate the eigenvalues of the superoperator $\mathcal{L}$ in subspaces $\mathcal{K}_{2}$ and $\mathcal{K}_{3}$. Let us denote the eigenvalues as $\lambda$ in each subspace. We find that, roughly, the slowest decay rate in $\mathcal{K}_{2}$ is two times faster than that in $\mathcal{K}_{1}$. Similarly, the decay rate in $\mathcal{K}_{3}$ is much faster than that in $\mathcal{K}_{2}$. The similar relations can be found in Ref. Glauber (1963) by a finite truncation of the hierarchy equations. Thus in the long time and long chain limits, we can fully characterize the relaxation time of the many-body system in terms of single particle decay rate.

With this density matrix, we can understand the dynamics in the many-body state for any given initial wave function. For example, for the result in Fig. 2 (c) and (e), we can express the dynamics of $|n\rangle\langle n|$ as

[TABLE]

where the coefficients can be determined by calculation the dynamics of the density matrix in each subspace $\mathcal{K}_{i}$. Initially, we may find that the coefficients $c_{a_{1},a_{2},\cdots,a_{2N}}$ and the overlap $\langle n|\alpha_{1}^{a_{1}}\alpha_{2}^{a_{2}}\cdots\alpha_{2N}^{a_{2N}}|n\rangle$ are in the order of unity, which after a long time relaxation will finally approach the steady solution in Eq. 7. For this reason not only the single particle terms, but also all the higher order terms, which correspond to the many-body relaxation, will contribute to the relaxation process. As a result the dynamics of the many-body state will exhibit multiply time scales during relaxation, in which some of these lifetime scales are shown in Fig. 3. Due to the presence of these multi-particle relaxations, we find that in general the decay rate of $|n\rangle\langle n|$ is faster than the single particle ones (see Fig. 2). Noticed that in Fig. 2 (e) and (f), the decay processes have the same oscillation period due to the finite coupling between the two edge modes, which in the many-body case will have the same energy splitting in the two ground state energy in the XY model (see Fig. 1 (a) and Eq. 4). However, in an odd chain, this oscillation will disappear due to the absence of coupling.

Relaxation in the long chain limit. The above results have established a connection between many-body dynamics and single particle dynamics in the long time limit. Some more issues need to be explained. (1) Why Re($\lambda_{2\sigma}$) exhibits an inflexion point at $\gamma_{c}\sim 1$; (2) What will happen in the long chain limit? and (3) What are the different roles played by the edge modes and bulk modes during relaxation? We focus on subspace $\mathcal{K}_{1}$, in which the dynamics of $\psi_{i}$ is given by the following non-Hermitian schrödinger equation (see Fig. 5 (a))

[TABLE]

where $\Gamma=\text{diag}(-\gamma,0,\cdots,0,-\gamma)$ and

[TABLE]

Here we have defined $\Psi=(\psi_{1},\psi_{4},\psi_{5},\cdots)$ in chain A, while its treatment for chain B is similar. The eigenvalues in this bordered matrix are determined by Kouachi (2006); Yueh (2005); Willms (2008)

[TABLE]

when $N=2m+1$ is odd; and

[TABLE]

when $N=2m$ is even, via $\text{Det}(H-\epsilon)$. In the above equations, $\epsilon$ is the eigenvalue and its relation to $\theta$ is determined by

[TABLE]

When $\gamma=0$ for odd chain, we have $\epsilon_{n}=\pm\sqrt{g_{1}^{2}+g_{2}^{2}+2g_{1}g_{2}\cos\theta_{n}}$ and $\theta_{n}=n\pi/(m+1)$ for $n=1,\cdots,m$ and $\epsilon_{0}=0$. Thus the energy gap in Fig. 1 (a) is given by $\varepsilon_{g}=|g_{1}-g_{2}|$.

We find that for the extended bands, the eigenvalues and the phase can be written as

[TABLE]

where $n\ll m$ and $\epsilon_{n,\text{i}}$, $z_{n,\text{r}}$, $z_{n,\text{i}}$ are small numbers in the sense that $\lim_{m\rightarrow\infty}mz_{n,i/r}=0$. These solutions can be obtained by linearizing the above nonlinear equations.

(a) In the weak dissipation limit ($\gamma\ll g_{1},g_{2}$) and in the odd (o) and even (e) chains, we have

[TABLE]

One can easily check our previous approximation that $\lim_{m\rightarrow\infty}mz_{n,i/r}=0$. The imaginary part of $\epsilon$ is responsible for the characteristic relaxation time $T^{*}$ as

[TABLE]

This result accounts for the multiply time scales during relaxation. When $\gamma\rightarrow 0$, $T^{*}\rightarrow\infty$, indicating of persistent coherent dynamics. Since $\epsilon_{n,\text{i}}\propto\gamma$, it means that in the weak dissipation limit, relaxation is still important and can happen in a finite system. Moreover, we find that the odd-even effect is still visible in the long chain limit.

(b) In the strong dissipation limit, the odd-even effect will vanish, and we find

[TABLE]

We are surprised to find that in the strong dissipation limit, the relaxation time $T^{*}\propto\gamma N^{3}$, thus it will be prolonged by the dissipation. The crossover between these two cases are determined by $\epsilon_{n,\text{i}}=\epsilon_{n,\text{i}}^{i}$ with $i=$e, o, which yields $\gamma_{c}=g_{1}$ in an even chain; and $\gamma_{c}=\sqrt{2}g_{1}g_{2}/\sqrt{g_{1}^{2}+g_{2}^{2}}$ in an odd chain. Thus the strong dissipation regime can be assigned by $\gamma>\gamma_{c}$. These inflexion points are also numerical verified, which are presented in Fig. 4 (b). The above solutions are also approximately correct even in a short chain in regarding of the fast decay of the imaginary energy according to $\epsilon_{n,\text{i}}\propto 1/N^{3}$.

For the localized edge modes, new decomposition is required. In the weak dissipation limit, the edge modes are not changed by the dissipation and we can make perturbation around $\cos(\theta_{c})=-{(g_{1}^{2}+g_{2}^{2})\over 2g_{1}g_{2}}$, which yields

[TABLE]

for the two modes at the left and right ends. In the strong dissipation limit, the edge modes will be fully localized at the two ends, thus we can assume $\epsilon\simeq-i\gamma$, and by perturbation about $\cos(\tilde{\theta}_{c})=-{(\gamma^{2}+g_{1}^{2}+g_{2}^{2})\over 2g_{1}g_{2}}$, we have

[TABLE]

in odd chian and

[TABLE]

in even chain. These results are independent of chain length, indicating that in the long chain limit the relaxation time is fully determined by the bulk bands, while the edge modes are most vulnerable to dissipation.

Understand of these results. Let us try to understand these anomalous results using perturbation theory (see Fig. 5). Firstly, in the weak dissipation limit we can treat $\Gamma$ as perturbation. Based on the first-order perturbation theory, we find that for the localized edge modes, we have $\epsilon_{\text{L}}=\imaginary\langle\psi_{\text{e}}|i\Gamma|\psi_{\text{e}}\rangle=-i\gamma|\psi_{\text{e}}(1)|^{2}=-i\gamma(1-\lambda^{2})$, with wave function at the left end as $\psi_{\text{e}}\sim(1,0,\lambda,0,\lambda^{2},\cdots)$ (see Fig. 5 (b)). In this case, the dissipation can be independent of total chain length. For the extended bands in an even chain, the wave function at each site $x$ is $\psi_{N\in\text{e}}^{n}(x\in\text{e})\propto\sin({n\pi x\over N+{2\lambda\over 1+\lambda}})$ and $\psi_{N\in\text{e}}^{n}(x\in\text{o})\propto\sin({n\pi(N+1-x)\over N+{2\lambda\over 1+\lambda}})$, where $n$ is the band index. In an odd chain, $\psi_{N\in\text{o}}^{n}(x\in\text{e})\propto\sqrt{\lambda^{2}+1+2\lambda\cos\theta_{n}}\sin({(N-x+1)\over 2}\theta_{n})$ and $\psi_{N\in\text{o}}^{n}(x\in\text{o})\propto\lambda\sin(({N-x\over 2}+1)\theta_{n})+\sin({L-x\over 2}\theta_{n})$, where $\theta_{n}={2n\pi\over(N+1)}$. Then we find $\epsilon_{n,\text{i}}=\langle\psi_{N}^{n}|\Gamma|\psi_{N}^{n}\rangle$, which will recover the expression in Eq. 17. In the strong dissipation limit, we need a different decomposition $H_{0}+i\Gamma=\mathcal{H}_{0}+i\Gamma+V$, where $\mathcal{H}_{0}=\mathcal{P}_{1}H_{0}\mathcal{P}_{1}$ and $V=\mathcal{P}_{2}H_{0}\mathcal{P}_{2}$ ($\mathcal{P}_{1}$ is the projector into the zero subspace of $i\Gamma$, while $\mathcal{P}_{2}$ is orthogonal to $\mathcal{P}_{1}$). After a few algebra, we find

[TABLE]

which contains only the coupling at sites $2$ and $L-1$, while all the couplings equal to zero exactly. $\mathcal{H}_{0}$ can be obtained by extracting $V$ from $H_{0}$ (see Eq. 12). We can treat $V$ as a perturbation, and via the second-order perturbation theory, we have

[TABLE]

where the extended bands $\tilde{\psi}_{N}^{n}(1)=\tilde{\psi}_{N}^{n}(N)=0$, $\tilde{\psi}_{N}^{n}(x)\propto\psi_{N-2}^{n}(x-1)$, with $x=2,\dots,N-1$ and the localized modes $|a\rangle=(1,0,0\cdots),|b\rangle=(\cdots 0,0,1)$. This expression yields Eq. 19 (see the overlap between the bulk modes and the edge dissipation in Fig. 5 (b)).

We present the dynamics of $\text{tr}(\rho\alpha_{\text{L}})$ in Fig. 4 (c), where $\alpha_{\text{L}}$ is the localized edge modes without dissipation. Since strong dissipation can influence the profile of the edge modes, $\alpha_{\text{L}}$ is no longer the eigenvector of $H_{0}+i\Gamma$. By a linear fitting we find

[TABLE]

where $\alpha\sim 2.3$. The exponent of $\alpha=-3$ is not reached due to the finite time window in simulation. This is because we only consider the dissipation of the extended bands $\tilde{\psi}_{N}^{n}$ from its overlap with the edge modes ($\langle\psi_{\text{e}}|\tilde{\psi}_{N}^{n}\rangle\neq 0$), thus the dynamics of $\text{tr}(\rho\alpha_{\text{L}})$ is dominated by the extended bands in the long time limit. However, this will not influence the scaling of $T^{*}$ with respect to $\gamma$ (see Fig. 4 (d)), since all the bulk bands have the same scaling law $T^{*}\propto\gamma$.

Conclusion. Dissipation in the many-body system is a fundamental problem that up to date has not yet been well understood. We explore the dissipation induced relaxation in a quantum $XY$ model with boundary dissipation, in which the relaxation is characterized by a characteristic time $T^{*}$. In the long-time limit it is fully determined by the dynamics of single particle physics. We explore the roles played by edge modes and bulk bands, and their scaling laws with respect to chain length and dissipation rate. An intuitive picture based on an equivalent non-Hermitian model is proposed. These results may also suggest that the lowest two states protected by a finite gap width can not serve as quantum memory under dissipation, in consistent with the conclusions in the previous literature Budich et al. (2012); Pedrocchi and DiVincenzo (2015); Rainis and Loss (2012); Mazza et al. (2013); Goldstein and Chamon (2011); Carmele et al. (2015).