Solution of the Lindblad equation for spin helix states

TL;DR

This paper presents an exact solution to the Lindblad equation for quantum spin systems with boundary-driven dissipation, revealing pure states with helical magnetization and superdiffusive spin currents in various dimensions.

Contribution

It introduces a method to solve the many-body Lindblad equation exactly for boundary-driven spin systems in any dimension, including non-integrable models, with explicit solutions for the stationary states.

Findings

Exact pure stationary states with helical magnetization profiles

Superdiffusive ballistic spin currents independent of system size

Applicability to higher-dimensional and non-integrable spin models

Abstract

Using Lindblad dynamics we study quantum spin systems with dissipative boundary dynamics that generate a stationary nonequilibrium state with a non-vanishing spin current that is locally conserved except at the boundaries. We demonstrate that with suitably chosen boundary target states one can solve the many-body Lindblad equation exactly in any dimension. As solution we obtain pure states at any finite value of the dissipation strength and any system size. They are characterized by a helical stationary magnetization profile and a superdiffusive ballistic current of order one, independent of system size even when the quantum spin system is not integrable. These results are derived in explicit form for the one-dimensional spin-1/2 Heisenberg chain and its higher-spin generalizations (which include for spin-1 the integrable Zamolodchikov-Fateev model and the bi-quadratic Heisenberg chain). The extension of the results to higher dimensions is straightforward.

Equations

$$H=\sum _{k=1}^{N-1}{h}_k$$

$$h_k$$

$$=$$

$${\hat{n}}_k=\frac{1}{2}\left(1-{\sigma }_k^z\right),{\hat{v}}_k=\frac{1}{2}\left(1+{\sigma }_k^z\right)$$

$$0=\frac{\mathrm{d}}{\mathrm{d}t}\rho=-i\mbox{$[{H},\,{\rho}]$}+\mathcal{D}_{L}(\rho)+\mathcal{D}_{R}(\rho)$$

$${\mathcal{D}}_j(\rho )=D_j\rho D_j^{†}-\frac{1}{2}\{D_j^{†}{D}_j,\rho \}.$$

$$\rho=\mbox{$|\,{\Phi}\,\rangle$}\mbox{$\langle\,{\Phi}\,|$}$$

$$\mbox{$|\,{\Phi}\,\rangle$}=\mbox{$|\,{\phi_{1}}\,\rangle$}\otimes\dots\otimes\mbox{$|\,{\phi_{N}}\,\rangle$}.$$

$$\rho=\mbox{$|\,{\phi_{1}}\,\rangle$}\mbox{$\langle\,{\phi_{1}}\,|$}\otimes\dots\otimes\mbox{$|\,{\phi_{N}}\,\rangle$}\mbox{$\langle\,{\phi_{N}}\,|$}.$$

$$\mbox{$|\,{\phi_{k}}\,\rangle$}=\frac{1}{\sqrt{|a|^{2}+|b|^{2}}}\left(\begin{array}[]{c}a\,\mathrm{e}^{-i\frac{1}{2}\phi_{k}}\\ b\,\mathrm{e}^{i\frac{1}{2}\phi_{k}}\end{array}\right)$$

$${\varphi }_k=\phi k$$

$$m_k^x=\frac{r}{1+r^2}\cos (\phi k-{\phi }_B),m_k^y=\frac{r}{1+r^2}\sin (\phi k-{\phi }_B),m_k^z=\frac{1}{2}\frac{1-r^2}{1+r^2}.$$

$$\phi =\frac{\Phi +2\pi K}{N-1}$$

$$\epsilon =J\left[{\left(\frac{2r}{1+r^2}\right)}^2\cos \phi +\Delta \left({\left(\frac{1-r^2}{1+r^2}\right)}^2-1\right)\right].$$

$$\mathcal{L}(\rho)=-i\mbox{$[{H},\,{\rho}]$}+\sum_{j}\mathcal{D}_{j}(\rho)=0$$

$$D_{j}\mbox{$|\,{\Psi}\,\rangle$}=\lambda_{j}\mbox{$|\,{\Psi}\,\rangle$},\quad\tilde{H}\mbox{$|\,{\Psi}\,\rangle$}=\mu\mbox{$|\,{\Psi}\,\rangle$}$$

$$\tilde{H}=H+\sum _j\frac{i}{2}\left({\bar{\lambda }}_j{D}_j-{\lambda }_j{D}_j^{†}\right).$$

$$\sum_{j}\left(\mbox{$\langle\,{\Psi}\,|$}D_{j}\mbox{$|\,{\Psi}\,\rangle$}\mbox{$\langle\,{\Psi}\,|$}D_{j}^{\dagger}\mbox{$|\,{\Psi}\,\rangle$}-\mbox{$\langle\,{\Psi}\,|$}D_{j}^{\dagger}D_{j}\mbox{$|\,{\Psi}\,\rangle$}\right)=0$$

$$D_{j}\mbox{$|\,{\Psi}\,\rangle$}=\lambda_{j}\mbox{$|\,{\Psi}\,\rangle$}$$

$$\mathcal{D}_{j}(\rho)=\frac{1}{2}\lambda_{j}\mbox{$[{\rho},\,{D_{j}^{\dagger}}]$}+\frac{1}{2}\bar{\lambda}_{j}\mbox{$[{D_{j}},\,{\rho}]$}=\mbox{$[{\frac{1}{2}(\bar{\lambda}_{j}D_{j}-\lambda_{j}D_{j}^{\dagger})},\,{\rho}]$}$$

$$\mbox{$[{H+\sum_{j}\frac{i}{2}(\bar{\lambda}_{j}D_{j}-\lambda_{j}D_{j}^{\dagger})},\,{\rho}]$}=0.$$

$$\mbox{$\langle\,{k}\,|$}A\mbox{$|\,{\Psi}\,\rangle$}\mbox{$\langle\,{\Psi^{\prime}}\,|\,{l}\,\rangle$}=\mbox{$\langle\,{k}\,|\,{\Psi}\,\rangle$}\mbox{$\langle\,{\Psi^{\prime}}\,|$}A\mbox{$|\,{l}\,\rangle$}$$

$$\frac{\mbox{$\langle\,{k}\,|$}A\mbox{$|\,{\Psi}\,\rangle$}}{\mbox{$\langle\,{k}\,|\,{\Psi}\,\rangle$}}=\frac{\mbox{$\langle\,{\Psi^{\prime}}\,|$}A\mbox{$|\,{l}\,\rangle$}}{\mbox{$\langle\,{\Psi^{\prime}}\,|\,{l}\,\rangle$}}\quad\forall k,l.$$

$$\mbox{$\langle\,{k}\,|$}A\mbox{$|\,{\Psi}\,\rangle$}=\mu\mbox{$\langle\,{k}\,|\,{\Psi}\,\rangle$},\quad\mbox{$\langle\,{\Psi^{\prime}}\,|$}A\mbox{$|\,{k}\,\rangle$}=\mu\mbox{$\langle\,{\Psi^{\prime}}\,|\,{k}\,\rangle$}\quad\forall k$$

$$A\mbox{$|\,{\Psi}\,\rangle$}=\mu\mbox{$|\,{\Psi}\,\rangle$},\quad\mbox{$\langle\,{\Psi^{\prime}}\,|$}A=\mu\mbox{$\langle\,{\Psi^{\prime}}\,|$}.$$

$$h_k=e_k(\eta )+i\sin \eta ({\sigma }_{k+1}^z-{\sigma }_k^z)=e_k(-\eta )-i\sin \eta ({\sigma }_{k+1}^z-{\sigma }_k^z)$$

$$e_k(\eta )=2J\left({\sigma }_k^{+}{\sigma }_{k+1}^{-}+{\sigma }_k^{-}{\sigma }_{k+1}^{+}-{\mathrm{e}}^{i\eta }{\hat{n}}_k{\hat{v}}_{k+1}-{\mathrm{e}}^{-i\eta }{\hat{v}}_k{\hat{n}}_{k+1}\right).$$

$$H=G(\eta )+iJ\sin \eta ({\sigma }_N^z-{\sigma }_1^z)=G(-\eta )-iJ\sin \eta ({\sigma }_N^z-{\sigma }_1^z)$$

$$\eta =\phi$$

$$e_{k}(\varphi)\mbox{$|\,{\Phi}\,\rangle$}=0,\quad\mbox{$\langle\,{\Phi}\,|$}e_{k}(-\varphi)=0.$$

$$H\mbox{$|\,{\Phi}\,\rangle$}=iJ\sin{\varphi}(\sigma^{z}_{N}-\sigma^{z}_{1})\mbox{$|\,{\Phi}\,\rangle$},\quad\mbox{$\langle\,{\Phi}\,|$}H=-iJ\sin{\varphi}\mbox{$\langle\,{\Phi}\,|$}(\sigma^{z}_{N}-\sigma^{z}_{1}).$$

Full text

${}^{~{}1}$Helmholtz-Institut für Strahlen-und Kernphysik, Universität Bonn, Nussallee 14-16, 53119 Bonn, Germany

Email: popkov@uni-bonn.de

${}^{~{}2}$Institute of Complex Systems II, Forschungszentrum Jülich, 52425 Jülich, Germany

Email: g.schuetz@fz-juelich.de

1 Introduction

A question of considerable interest in the context of one-dimensional transport phenomena is the magnitude of stationary currents in boundary-driven quantum spin systems as a function of system size $N$. In the case of normal (diffusive) transport a current $j$ is asymptotically proportional to $1/N$, while for ballistic transport the current approaches a non-zero constant even in the thermodynamic limit $N\to\infty$. In one dimension this behavior is a hallmark of integrable systems and manifests itself in a finite Drude weight [1, 2]. A way to measure this quantity experimentally in such systems has been proposed recently [3].

We address the relationship between the nature of the boundary driving, integrability and transport properties by studying boundary-driven quantum spin chains in the by now theoretically well-established and experimentally accessible framework of non-equilibrium Lindblad dynamics. This approach models a dissipative coupling of a quantum system to its environment and thus allows for the description of stationary current-carrying quantum states. We explore conditions on the boundary driving under which ballistic transport may occur in a quantum spin system. It turns out that such behavior arises in stationary states in which the ballistic current is associated with a spin rotation along the direction of driving. We shall call such superdiffusive nonequilibrium stationary states “spin helix states” (SHS), in analogy to phenomena in spin-orbit-coupled two-dimensional electron systems [4, 5, 6]. We focus on one-dimensional spin chains, which are of great current interest. However, it will transpire that analogous SHS will appear also in higher dimensions with an appropriate choice of Lindblad boundary driving.

The 1-d SHS generalizes the asymptotic state in the isotropic Heisenberg chain ($XXX$-chain) in the thermodynamic limit $N\to\infty$ that was found recently [7, 8] which is, in turn, reminiscent of the helical ground state of the classical isotropic Heisenberg spin chain with boundary fields and its formal analog of ferromagnetic quantum domains in the Heisenberg quantum chain [9, 10]. The novelty of the SHS is the occurrence of a non-zero winding number in the helical state that turns out to be responsible for the ballistic transport.

Mainly we are interested in exact SHS’s in the experimentally relevant chains of finite length. However, we shall also present numerical results away from the exactly solvable points that highlight the specific features of the exact SHS. Interestingly, these SHS are pure states, which is unusual for solutions of a many-body Lindblad equation. These states arise in the regime $|\Delta|<1$ for the anisotropy parameter of the spin-$s$ chain. For the ground state of the spin-1/2 XXZ Heisenberg chain this is the quantum critical regime, unlike the ferromagnetic regime $\Delta\geq 1$ studied in [10], which exhibits a mathematically somewhat analogous but physically very different behavior. Notice that the nonequilibrium stationary state of a dissipatively boundary driven $XXZ$-chain was argued to converge to the SHS in the Zeno limit of infinitely large boundary dissipation [11, 12]. Here we show how the SHS is produced at arbitrary finite dissipative strength.

The paper is organized as follows. To be concrete, we first consider in Sec. 2 the anisotropic spin-1/2 Heisenberg chain. We define the SHS and derive the conditions under which exact SHS’s arise with judiciously chosen Lindblad dissipators. In Sec. 3 we discuss in some detail transport properties of the spin-1/2 SHS and compare with transport in non-SHS states. Then we go on to generalize the approach to higher-spin chains (Sec. 4) and discuss some classical analogies. In Sec. 5 we draw some conclusions.

2 Spin helix states in the spin-1/2 $XXZ$-chain

The spin-1/2 $XXZ$-chain is defined by the Hamiltonian [13]

[TABLE]

with local interaction matrices $h_{k}$ given in terms of Pauli spin-1/2 matrices by

[TABLE]

Here $\Delta=\cos{\eta}$ is the anisotropy parameter, and in the second representation we have used the local projectors

[TABLE]

and the spin raising and lowering operators $\sigma_{k}^{\pm}=(\sigma_{k}^{x}\pm i\sigma_{k}^{y})/2$. We recall that the Pauli matrices satisfy the $SU(2)$ commutation relations \mbox{[{\sigma^{\alpha}{k}},,{\sigma^{\beta}{l}}]}=2i\delta_{k,l}\sum_{\gamma=1}^{3}\epsilon_{\alpha\beta\gamma}\sigma^{\gamma}_{k} where $\epsilon_{\alpha\beta\gamma}$ is the totally antisymmetric Levi-Civita symbol with $\epsilon_{123}=1$.

The object of interest is the density matrix $\rho$ in a boundary-driven non-equilibrium situation where stationary currents arise from the coupling of the left and right boundary sites $1$ and $N$ to an environment which projects the boundary spins in different directions. The density matrix $\rho$ of the non-equilibrium steady state (NESS) is determined by the stationary Lindblad equation [14, 15]

[TABLE]

with boundary dissipators $\mathcal{D}_{j}$, $j\in\{L,R\}$ acting on the density matrix as

[TABLE]

The Lindblad operators $D_{j}$ which encode the nature of the boundary driving will be specified below. Stationary expectations $\langle\,{O}\,\rangle$ of physical observables $O$ are then given by the trace \mbox{\langle,{O},\rangle}=\mathop{\mathrm{Tr}}\nolimits(O\rho). Our main interest will be in the magnetic moments $\vec{m}_{k}$ at site $k$ of the chain. For convenience we ignore material-dependent factors and choose units such that \vec{m}_{k}=\mbox{\langle,{\vec{\sigma}_{k}},\rangle}.

In the absence of the unitary part given by the spin chain Hamiltonian $H$, the non-unitary dissipative part given by the dissipators $\mathcal{D}_{j}$ forces the system locally at the respective left (L) or right (R) boundary site into some target state. Thus, if the two target states are different, stationary currents associated with local bulk-conserved degrees of freedom are generally expected to flow due to the action of the unitary bulk part of the Lindblad equation.

2.1 The spin-1/2 helix state

For many problems of interest the quantum master equation (5) admits an exact solution in which the stationary density matrix is expressed in matrix product form [16, 17]. Here we take a different approach and make a pure-state ansatz

[TABLE]

with the product state

[TABLE]

This means that we can write

[TABLE]

We take the basis where the $z$-components $\sigma_{k}^{z}$ of the local spin operator are all diagonal and choose

[TABLE]

with the local phase angle

[TABLE]

where $0\leq\varphi<2\pi$.

With the parametrization $a=\mathrm{e}^{i\varphi_{B}/2}$, $b=r\mathrm{e}^{-i\varphi_{B}/2}$ the magnetization profiles m^{\alpha}_{k}:=\mbox{\langle,{\sigma^{\alpha}_{k}},\rangle}/2, i.e., the $\alpha$-components of the dimensionless magnetic moments, are given by

[TABLE]

One recognizes in $\varphi$ the twist angle between neighbouring spins in the $xy$-plane. Therefore we refer to the pure density matrix (9) specified by the properties (10) and (11) as spin helix state (SHS).

The quantity $\varphi(N-1)$ yields the twist angle between boundary target polarizations in the $xy$-plane. Hence any $\varphi\in[0,2\pi[$ of the form

[TABLE]

with $0\leq\Phi<2\pi$ and $0\leq K<N-1$ gives rise to the same spin rotation between the boundary spins by the angle $\Phi$ in the $xy$-plane. We shall refer to $\Phi$ as the boundary twist and to $K$ as the (clockwise) winding number of the spin helix [18]. Without loss of generality we fix the phase $\varphi_{B}=\varphi$ which corresponds to a choice of the coordinate system such that the planar spin component at site 1 points into the $x$-direction. The left target state at site 1 is then the local density matrix $\rho_{L}=(\hat{v}+r^{2}\hat{n}+r\sigma^{x})/(1+r^{2})$. and the right target state is given by $\rho_{R}=(\hat{v}+r^{2}\hat{n}+r\cos{(\Phi)}\sigma^{x}+r\sin{(\Phi)}\sigma^{y})/(1+r^{2})$. For $r=1$ the SHS is fully polarized in the $xy$-plane with perpendicular magnetization $m^{z}_{k}=0$ along the chain. Due to the factorized structure of the SHS there are no spin-correlations between different sites.

Thermal-like properties of this NESS can be characterized by the bond energy density \varepsilon_{k}:=\mbox{\langle,{h_{k}},\rangle}. From the factorization property (9) and the explicit form of the local magnetizations (12) one finds that the bond energy density is spatially constant and given by

[TABLE]

Due do the factorized structure of the SHS there are no energy correlations between non-neighbouring bonds.

The complete absence of correlations in the SHS is reminiscent of very high temperatures. We caution, however, not to interpret this lack of correlations and the flat energy profile along the chain as indicating proximity to some equilibrium state $\rho\propto\exp{(-\beta_{eff}H)}$ with an effective temperature given by (14), not even if $\varphi=0$ when also the magnetization profile is flat. For $\varphi=0$ one can write $\rho\propto\exp{(-\beta_{eff}H_{eff})}$ with an effective Hamiltonian of the form $H_{eff}=\sum_{k}(\sigma_{k}^{z}+u\sigma_{k}^{x}+w)$. Such a non-interacting Hamiltonian corresponds to a subspace of $H$ for $\Delta=0$ [19], but does not in general capture any significant physical property of the thermal density matrix $\rho\propto\exp{(-\beta H)}$ for any finite temperature at any value of $\Delta$.

2.2 Construction of the boundary dissipators

Now we aim at deriving boundary dissipators which allow for maintaining the SHS stationary in the finite $XXZ$-chain. To this end we first make a remark on pure-state solutions of a general stationary Lindblad equation

[TABLE]

where here $j$ belongs to some index set (not necessarily just $L$ and $R$). Let a pure state \rho=\mbox{|,{\Psi},\rangle}\mbox{\langle,{\Psi},|} be the solution of (15). Then $|\,{\Psi}\,\rangle$ is an eigenvector of all the Lindblad operators $D_{j}$ and the Lindblad equation turns into the set of eigenvalue problems

[TABLE]

with (in general complex) eigenvalues $\lambda_{j}$ and (real) eigenvalue $\mu$ of the shifted Hamiltonian

[TABLE]

This can be seen as follows [20, 21]. Sandwich the Lindblad equation (15) with $|\,{\Psi}\,\rangle$. Then the unitary part involving the commutator with $H$ vanishes identically and one gets

[TABLE]

for the dissipative part. By the Schwarz inequality (which generally gives $\geq 0$ for the l.h.s.) the equality is realized if and only if the eigenvalue property

[TABLE]

holds for each dissipative term. Then the Lindblad dissipator can be written as a commutator

[TABLE]

and the Lindblad equation becomes

[TABLE]

Consider now the commutator \mbox{[{A},,{\sigma}]}=0 with a general tensor matrix \sigma=\mbox{|,{\Psi},\rangle}\mbox{\langle,{\Psi^{\prime}},|} such that \mbox{\langle,{k},|,{\Psi},\rangle}\neq 0 and \mbox{\langle,{\Psi^{\prime}},|,{l},\rangle}\neq 0 for all orthonormal basis vectors $|\,{k}\,\rangle$, $|\,{l}\,\rangle$ of the separable Hilbert space to which $|\,{\Psi}\,\rangle$ and $|\,{\Psi^{\prime}}\,\rangle$ belong. Sandwiching with $\langle\,{k}\,|$ and $|\,{l}\,\rangle$ yields

[TABLE]

or, equivalently,

[TABLE]

Hence

[TABLE]

with the same constant $\mu$. This implies

[TABLE]

This proves (16) for any pure state. Conversely, if (16) holds for some vector $|\,{\Psi}\,\rangle$ then the pure state \rho=\mbox{|,{\Psi},\rangle}\mbox{\langle,{\Psi},|} is a solution of the original Lindblad equation (15).

Now we apply this property to the SHS defined by (9) with (10), (11) which we require to satisfy the stationarity condition (5) with boundary Lindblad operators $D_{L,R}$. Notice that one can write the interaction terms $h_{k}$ of the $XXZ$-Hamiltonian (1) as

[TABLE]

with

[TABLE]

This fact allows us to write

[TABLE]

with $G(\eta)=\sum_{k=1}^{N-1}e_{k}(\eta)$.

Remarkably, for the relation

[TABLE]

between the twist angle $\varphi$ of the SHS and the anisotropy $\eta$ of the $XXZ$-chain one has

[TABLE]

This implies G(\varphi)\mbox{|,{\Phi},\rangle}=0 and \mbox{\langle,{\Phi},|}G(-\varphi)=0 and therefore

[TABLE]

To proceed and construct suitable Lindblad operators $D_{L,R}$ it is convenient to define for subscript $j\in\{L,R\}$ the shifted Lindblad operators

[TABLE]

We also note that we can write the shifted Hamiltonian (17) as

[TABLE]

The constants $\lambda_{j}$ are to be determined. According to (16) this implies that one has to solve

[TABLE]

and

[TABLE]

with $\mu\in{\mathbb{R}}$. Here we used that (34) is equivalent to \mbox{\langle,{\Phi},|}\tilde{D}^{\dagger}_{i}=0. This allows us to split these four equations into two pairs of equations for each boundary

[TABLE]

[TABLE]

with $\mu_{L}=(\mu+i\nu)/2$ arbitrary and $\mu_{R}=\bar{\mu}_{L}$ so that $\mu_{L}+\mu_{R}=\mu\in{\mathbb{R}}$ as required by (16). The real-valued constants $\mu,\nu$ can be computed by multiplying from the right by $|\,{\Phi}\,\rangle$. Using (12) yields

[TABLE]

and therefore $\mu=0$, $\nu=j^{z}$. For full planar polarization this reduces to $\mu_{L}=\mu_{R}=0$.

Requiring the left dissipator $D_{L}$ to act non-trivially on the left boundary site 1 one finds from the first eigenvalue equation in (36) that

[TABLE]

with arbitrary constants $\alpha_{L},\beta_{L}$. Then the second equation in (36) is solved by

[TABLE]

For the right boundary the eigenvalue equation \tilde{D}_{R}\mbox{|,{\Phi},\rangle}=0 in (37) gives

[TABLE]

with arbitrary constants $\alpha_{R},\beta_{R}$. From the second equation in (37) one then obtains

[TABLE]

Thus the SHS is stationary under the action of a two-parameter family of boundary dissipators with Lindblad operators $D_{j}=\tilde{D}_{j}+\lambda_{j}$.

3 Transport properties of the SHS

We treat both spin and energy transport, the emphasis being on spin transport.

3.1 Spin transport in the SHS

The $z$-component of the total magnetization is conserved under the unitary part of the time evolution. The associated conserved spin current is defined by the continuity equation through the time derivative of the magnetization profile $\dot{m}^{z}_{k}=j^{z}_{k-1}-j^{z}_{k}$. Since \dot{m}^{z}_{k}=i\mbox{\langle,{\mbox{$[{H},\,{\sigma^{z}_{k}}]$}},\rangle}/2 one gets from the commutation relations of the Pauli matrices the current operator

[TABLE]

In the stationary state the current j^{z}:=\mbox{\langle,{\hat{\jmath}^{z}_{k}},\rangle} does not depend on $k$ and it is of interest to investigate its properties in the SHS. Strictly speaking, the SHS as defined above arises as stationary solution of the Lindblad equation for a finite chain only in the regime $|\Delta|<1$ of the $XXZ$-chain. However, as shown below, it appears asymptotically also in the isotropic Heisenberg chain with $\Delta=1$ and it has a (non-helical) analog in the ferromagnetic regime $\Delta>1$. We discuss these cases separately.

3.1.1 Helical regime $|\Delta|<1$

The factorized form of the SHS defined by (9) - (11) yields

[TABLE]

which even in a large system is of order 1 for macroscopic winding numbers of order $N$. Interestingly, in contrast to the classical relation between a locally conserved current and boundary gradients of the associated conserved quantity, for any winding number there is a current even though there is no gradient $\Delta m^{z}:=m^{z}_{1}-m^{z}_{N}=0$ between the $z$-magnetizations of the boundaries. Moreover, the behaviour of the SHS is also in contrast to the situation where the $XXZ$-chain is driven by two Lindblad operators at each boundary into a state close to an infinite-temperature thermal state [22]. In this case, the effective diffusion coefficient $D^{z}_{eff}\propto Lj^{z}/\Delta m^{z}$ was found numerically for chains up to more than 200 sites to be proportional to $L$ (corresponding to ballistic transport) with a coefficient of proportionality that depends on the anisotropy $\Delta$. Theoretically, a ballistic spin current in this regime was proved by calculating the lower bound for a respective Drude weight, see [2].

The spin transport of the SHS is, in fact, reminiscent of the persistent current $j$ in a mesoscopic ring threaded by a magnetic flux $\Phi$ [23, 24]. At zero temperature one has

[TABLE]

and the Drude weight is given by the spin stiffness [28]

[TABLE]

where $E_{0}$ is the ground state energy and $\Phi_{m}$ is the value of $\Phi$ that minimizes $E_{0}(\Phi)$. Substituting the ground state energy $E_{0}$ of the ring by the energy density (14) times the chain length $L=N-1$ (in lattice units) of the SHS, i.e., $E_{0}\to(N-1)\varepsilon$, identifying the flux $\Phi$ with the magnitude of the boundary twist, and keeping $\Delta$ fixed when taking the derivative w.r.t. $\Phi$ one finds from (48) that $j=j^{z}$ as given by (47) and then (49) gives $D_{SHS}=|J|>0$, indicating infinite DC conductivity.

Expressions for finite temperature analogous to (48) and (49) are derived in [25] and it was conjectured that a finite Drude weight at non-zero temperature is a generic property of integrable systems. Thus the non-thermal (but certainly not zero-temperature) SHS of the integrable XXZ-chain appears to fit into the picture relating the Drude weight obtained via (49), infinite DC conductivity and integrability [1, 2, 26, 27]. The Drude weight $D_{SHS}$, however, does not depend on the anisotropy $\Delta$ unlike the thermal Drude weight [28, 29, 30]. More significantly, however it will be shown below that the ballistic transport in the SHS is, in fact, unrelated to integrability.

3.1.2 Isotropic point $|\Delta|=1$

At the isotropic point $\Delta=1$ where $\eta=0$ and the matching condition (29) yields a trivial constant SHS with twist angle $\Phi=0$ and winding number $K=0$. However, it is interesting to look at the magnetization profiles (12) and the spin current (47) with the boundary driven isotropic $XXX$-chain, corresponding to non-zero boundary twist $\theta\neq 0$ in the $xy$-plane. It was shown in [7, 8] that the boundary target states and the magnetization profiles for large $N$ are of the form (12) with $\varphi=\theta/(N-1)$ and $r=1$. Thus this non-equilibrium steady state of the $XXX$-chain is a SHS in the thermodynamic limit with winding number $K=0$ and boundary twist $\Phi=\theta$.

The $z$-component of the spin current in the $XXX$-chain is asymptotically given by $j^{z}\sim J\theta/N$ [8], which agrees with (47) for $\varphi=\theta/(N-1)$ and large $N$ [31]. Moreover, one can show that in the $XXX$-case one has $\Delta m^{z}:=m^{z}_{1}-m^{z}_{N}=O(1/N)$, indicating ballistic transport of the $z$-component of the spin in the $XXX$-chain since the effective diffusion coefficient $D^{z}_{eff}=Nj^{z}/(\Delta m^{z})$ is proportional to system size $N$. This is consistent with the observation of infinite conductivity in the SHS of the XXZ-chain obtained above from the Drude weight (49) which is finite also for $\Delta=1$ [32].

However, the ballistic transport in the SHS of the $XXX$ chain is in contrast to the transport properties both of the canonical ensemble for which it has been shown that the spin stiffness of the periodic $XXX$-chain at zero $z$-magnetization vanishes at any positive temperature [33] and of the “infinite-temperature” $XXX$-chain with two Lindblad operators at each boundary, reported in [34]. According to exact numerical calculations for short chains up to approx. 10 sites the diffusion coefficient seems to diverge superdiffusively with system size as $D^{z}_{eff}=\propto N^{1/2}$ in this rather different setting. This is remarkable as it implies that the microscopic details of the Lindblad boundary dissipators may determine fundamentally qualitative properties of the bulk.

3.1.3 Ferromagnetic coupling $\Delta>1$

The Heisenberg Hamiltonian with $J<0$ and $\Delta>1$ (corresponding to a purely imaginary anisotropy parameter $i\eta$) has a degenerate ferromagnetic ground state with all spins aligned in positive or negative $z$-direction, corresponding to the SHS with $r=0$ or $r=\infty$ respectively. We note, however, that the SHS with $r$ finite can be defined also for purely imaginary $\varphi$ and therefore the matching condition (29) can be met for $\Delta>1$. However, this state is not a helix state. Substituting $\varphi\to i\eta$ and parametrizing $r=\exp{(u^{\ast}N\eta+i\phi_{0})}$ one obtains for the Heisenberg chain (1) with $\Delta=\cosh{\eta}$ a fully polarized state with vanishing spin current $j^{z}$ and the magnetization profiles given by

[TABLE]

where $\tilde{k}=k-u_{0}N$.

This is the domain wall state of the XXZ-chain with opposite boundary fields in $z$-direction [10] with a left domain of negatively aligned spins and a right domain with positively aligned spins. For $N\gg 1/\eta^{2}$ the domain wall between positive and negative aligned spins is located at $u_{0}N$, provided that $0<u_{0}<N$. Otherwise one has a boundary layer with a width of order $1/\eta$. Only in a region of size $O(1/\eta^{2})$ near the domain wall one has for large $N$ a non-negligible transverse magnetization $m_{k}^{x,y}$. This domain wall state has a direct classical analog as stationary traffic jam state of the asymmetric simple exclusion process with reflecting boundary conditions [35, 36] since for $\Delta>1$ the XXZ-Hamiltonian coincides with the generator of this stochastic interacting particle system [37]. Note that also the state (50) can be dissipatively obtained for infinite dissipation strength in a $XXZ$ chain with fine-tuned anisotropy $\Delta=\cosh\eta$ [11].

3.2 Energy transport in the SHS

The operator for the locally conserved energy current $\hat{\jmath}^{E}_{k}$ associated with bond $(k,k+1)$ is defined by the continuity equation \dot{h}_{k}=i\mbox{[{H},,{h_{k}}]}=\hat{\jmath}^{E}_{k}-\hat{\jmath}^{E}_{k+1} which yields \hat{\jmath}^{E}_{k}=i\mbox{[{h_{k-1}},,{h_{k}}]} [38, 39]. Using the commutation relations of the Pauli matrices one finds

[TABLE]

The energy current j^{E}=\mbox{\langle,{\hat{\jmath}^{E}_{k}},\rangle} then follows from the factorized structure (8) of the SHS and the magnetization profiles (12).

Somewhat surprisingly

[TABLE]

since $\Delta=\cos{\varphi}$ in the SHS. This is consistent with the constant bond energy along the chain (implying the absence of a energy gradient between the boundaries), but nevertheless not completely obvious since (a) from a microscopic perspective it is not a priori clear that the dissipators would not generate an energy current and (b) the total energy current $\sum_{k}\hat{\jmath}^{E}_{k}$ in a periodic chain is a conserved charge of the integrable periodic $XXZ$-chain [38, 39] and hence ballistic transport of energy is generic.

3.3 Numerical results

Now we explore numerically on a concrete example the predicted special properties of the spin helix state as opposed to a generic non-equilibrium state that arises as a solution of the Lindblad equation (5) with Lindblad operators whose parameters do not satisfy the matching condition (29) and conditions (39) - (45) for the Lindblad operators. We focus on the fully polarized SHS with $r=1$ and fix the Heisenberg exchange coupling $J=1$.

For the numerically exact solution of the Lindblad equation we consider an $XXZ$-chain of four sites. For the Lindblad operators we take $\alpha_{L}=\beta_{L}=\alpha_{R}=\beta_{R}=\sqrt{\Gamma}>0$ so that

[TABLE]

For $N=4$ we take $\varphi=2\pi/3$ corresponding to winding number $K=2$ and a zero boundary twist angle $\Phi=0$ in the $xy$-plane. By fixing $\epsilon_{R}=-\epsilon_{L}=0.05$ the variable $\Gamma$ becomes a measure for the dissipative strength. The pure SHS (9) - (11) is then a stationary solution of the Lindblad equation (5) for

[TABLE]

For the purpose of the numerical investigation we do not require these equations to be satisfied and study the purity of the solution of (5) and the corresponding stationary current $j^{z}$ as a function of the anisotropy $\Delta=\cos{\eta}$ and the dissipative strength $\Gamma$.

As a measure for the purity of the nonequilibrium steady state (NESS) $\rho$, we choose the von Neumann entropy $S=-\mathop{\mathrm{Tr}}\nolimits(\rho\log_{2}\rho)$. Notice that $S=0$ if and only if the NESS is a pure state. From the exact numerical solution of (5) with $\eta=\varphi$ one sees that indeed for the value of $\Gamma$ predicted by (54) the NESS becomes pure (Fig. 1). The spin current is maximal in amplitude near this point, but remains approximately equally strong for all $\Gamma\gtrsim 4$.

It is also instructive to look at the NESS as a function of the anisotropy $\Delta=\cos{\theta}$, i.e., now we assume the dissipative strength to satisfy (54), but not $\eta$. In this way, we see a resonance-like behaviour of various system observables around the critical value of the anisotropy $\Delta=\cos{\varphi}$. Even for a small chain of only 4 sites the spin current $j^{z}$ increases by an order of magnitude and changes its sign near the critical anisotropy, see Fig. 2. The von-Neumann entropy vanishes at $\Delta=\cos{\varphi}$, as expected. At the $XXX$-point $\Delta=1$ the von-Neumann entropy is small, but non-zero, in agreement with the notion that the SHS is attained only asymptotically. Also the current at this point as expected from the exact result [8]. For non-zero boundary twist $\Phi$ one obtains qualitatively similar behavior (data not shown).

In order to get some insight in the resonance-like behaviour we note the following. For large amplitude $\Gamma$, the dissipative part of the dynamics, which is quadratic in amplitudes, becomes much larger than the unitary Hamiltonian part of the dynamics, and as a result the boundary spins $1,N$ “freeze” for any $\Delta$. By this we mean that the states to which the dissipation projects the boundary spins, which are mixed states, become very close to completely polarized pure states. At the left boundary, the spin $1$ fixates approximately along the vector $(1,0,0)$ and at the right boundary approximately in the direction $(\cos\varphi(N-1),\sin\varphi(N-1),0)=(\cos{\Phi},\sin{\Phi},0)$. Indeed, analyzing the kernel of the left dissipator, we find that the distance from the actually targeted state and the pure fully polarized state at the left boundary, characterized via $\epsilon:=1-Tr(\rho_{1})^{2}$ with the reduced density matrix $\rho_{1}=Tr_{2,3,\ldots N}\rho$ is proportional to $\varepsilon\sim\Gamma^{-4}$ for large $\Gamma$. The same is true for the right boundary. Now, if the polarization of the leftmost and rightmost spins in the chain differ only slightly (in our example this boundary twist angle is actually zero $\Phi=0$), then one expects almost no current in the system for any $\Delta$ since it will generically favor a homogeneous spin configuration, the neighbouring spins at sites $k,k+1$ being almost collinear. This picture is well borne out by Fig. 2, except close to the critical value $\Delta=\cos{\varphi}$. At this point the spins arrange in the helix structure with a non-zero winding number (2 in our case) which gives rise to the resonance. For the exact helix spin state the spin current takes the value $j_{z}=\sin\varphi\approx-0.866$, close to the maximal possible spin current $|j^{z}_{max}|=1$.

4 Higher-spin chains

The above results can be generalized to the case of spin $s$ with maximal $z$-component $s^{z}=s=(n-1)/2$. We focus on spin chains with conserved $z$-component of the total spin.

4.1 Spin-$s$ chains with conserved $S^{z}$ component

In order to define the Hamiltonian $H$ we introduce the $n$-dimensional matrices $E^{pq}$ with matrix elements $(E^{pq})_{mn}=\delta_{p,m}\delta_{q,n}$. They satisfy the quadratic algbra

[TABLE]

From these we build the local operator

[TABLE]

for the $z$-component of the local spin as well as the total $z$-component

[TABLE]

We assume a local nearest neighbour interactions between spins, i.e.,

[TABLE]

This notation means that the nearest neighbour interaction matrix

[TABLE]

of dimension $n^{2}$ has matrix elements $h_{pn+q+1,p^{\prime}n+q^{\prime}+1}=c^{pq}_{p^{\prime}q^{\prime}}$. The coupling constants satisfy $c^{pq}_{p^{\prime}q^{\prime}}=\bar{c}^{p^{\prime}q^{\prime}}_{pq}$ since $H$ is hermitian. Moreover, we impose the ice rule [13]

[TABLE]

and the symmetry relation

[TABLE]

The ice rule (61) ensures conservation \mbox{[{H},,{\hat{S}^{z}}]}=0 of the $z$-component of the total magnetization and (62) corresponds to lattice reflection symmetry $k\leftrightarrow N+1-k$. We shall also investigate the special case of spin-flip symmetry

[TABLE]

which is the invariance under $S^{z}\leftrightarrow-S^{z}$. Requiring in addition time-reversal symmetry gives the constraints

[TABLE]

on the phases of the coupling coefficients.

4.2 Spin-$s$ helix state

We target a NESS in the form of a pure SHS |\,{\Psi}\,\rangle$$\langle\,{\Psi}\,| with \mbox{|,{\Psi},\rangle}=\mbox{|,{\Psi_{1}},\rangle}\otimes\dots\otimes\mbox{|,{\Psi_{N}},\rangle} and

[TABLE]

with non-zero constants $r_{i}$ that can be complex. In order to achieve this state in a similar fashion as discussed above for $s=1/2$, it is sufficient to require the generalization

[TABLE]

of the telescopic property (31) with diagonal matrices $F_{k}=\sum_{p=0}^{2s}f_{p}E^{pp}_{k}$.

This condition will be satisfied if

[TABLE]

is satisfied for all $k$. In order to see what this implies for the coupling constants $c^{pq}_{p^{\prime}q^{\prime}}$ we define the gauge transformation

[TABLE]

and rewrite the SHS in the form

[TABLE]

where $|\,{\Psi_{0}}\,\rangle$ represents the constant wave function. Consequently, multiplying (67) by $V_{\varphi}$ from the left and noting that $V_{\varphi}$ and $F$ are diagonal matrices, we obtain

[TABLE]

for all $k$. From the definition one finds $V_{\varphi}E_{k}^{pp^{\prime}}V_{\varphi}^{-1}=\mathrm{e}^{ik(p^{\prime}-p)}$ and therefore, using the ice rule,

[TABLE]

Moreover, one has

[TABLE]

Therefore

[TABLE]

On the other hand,

[TABLE]

Thus

[TABLE]

determines the coupling constants of the spin-$s$ chain (58).

This linear system of equations for the coupling constants of the Hamiltonian can be easily solved which we demonstrate for the first non-trivial case $s=1$. Notice that the case $s=1/2$ reproduces the $XXZ$-Hamiltonian discussed earlier.

4.3 Spin-1 chain

The ice rule (61) allows for 19 non-vanishing coupling constants. Hermiticity and reflection symmetry (62) leave as free parameters the real-valued diagonal elements $a_{p}:=c^{pp}_{pp}$, $b_{1}:=c^{01}_{01}=c^{10}_{10}$, $b_{2}:=c^{02}_{02}=c^{20}_{20}$, $b_{3}:=c^{21}_{21}=c^{12}_{12}$ and the spin-flip coefficients $c_{1}:=c^{01}_{10}=c^{10}_{01}\in{\mathbb{R}}$, $c_{2}:=c^{02}_{20}=c^{20}_{02}\in{\mathbb{R}}$, $c_{3}:=c^{12}_{21}=c^{21}_{12}\in{\mathbb{R}}$, $d:=c^{11}_{02}=c^{11}_{20}$, $\bar{d}:=c^{02}_{11}=c^{20}_{11}$. Requiring also spin-flip symmetry (63) leads to the further relations $a_{3}=a_{1}$, $b_{3}=b_{1}$, $c_{3}=c_{1}$. Time-reversal symmetry then implies $\bar{d}=d$.

4.3.1 Computation of $h$ for helix states

We define

[TABLE]

The parameters $\varphi,\zeta$, or equivalently $\delta,\zeta$, characterize the spin-1 helix state. In particular, one has \mbox{\langle,{S_{k}^{x}},\rangle}=2\sqrt{2\zeta}/(1+2\zeta)\cos{(\varphi(k-1))}, \mbox{\langle,{S_{k}^{y}},\rangle}=2\sqrt{2\zeta}/(1+2\zeta)\sin{(\varphi(k-1))}, \mbox{\langle,{S_{k}^{z}},\rangle}=0, and the amplitude attains its maximum of full polarization at $\zeta=1/2$. We exclude from the discussion the non-helical zero-current states $\varphi=0,\pi$ corresponding to $|\delta|=1$ and the non-helical states $\zeta=0,\infty$ with vanishing spin polarization \mbox{\langle,{\vec{S}_{k}},\rangle}=\vec{0}.

The full set of equations (75) for the spin-1 SHS reads

[TABLE]

Therefore

[TABLE]

and $a_{1}=-2d\zeta\delta$, $b_{2}=-c_{2}\cos{(2\varphi)}-\bar{d}\zeta^{-1}\delta$.

Since $b_{2}$ and $c_{2}$ are both real we conclude that also $d\zeta$ and $\bar{d}\zeta^{-1}$ must be real which implies that $d$ has the negative phase of $\zeta$ plus a multiple of $\pi$. For the coefficients $f_{i}$ one finds

[TABLE]

In addition we have

[TABLE]

which yields the consistency condition $c_{2}\sin{(2\varphi)}=(c_{1}+c_{3}-\bar{d}\zeta^{-1})\sin{\varphi}$ which is automatically satisfied for the irrelevant cases $\varphi=0,\pi$ and which otherwise yields

[TABLE]

Thus all parameters are expressed in terms of $\zeta,\varphi$ characterizing the helix state and the three real-valued parameters $c_{i}$ that can be chosen freely.

With the shorthand $h_{k}\equiv h_{k}(c_{1},c_{2},c_{3};\zeta,\varphi)$ we arrive at

[TABLE]

We also note that

[TABLE]

The constant $f_{1}$ is arbitrary since only the difference $F_{k+1}-F_{k}$ and the telescopic sum $\sum_{k=1}^{N-1}(F_{k+1}-F_{k})=F_{N}-F_{1}$ appear in calculations. Hence we can set $f_{1}=0$.

For spin-flip symmetry and time-reversal symmetry where $c_{3}=c_{1}$ and $\bar{\zeta}=\zeta$ the local interaction reduces to

[TABLE]

where $h^{\ast}_{k}(c_{1},c_{2};\zeta,\varphi):=h_{k}(c_{1},c_{2},c_{1};\zeta,\varphi)$. The corresponding divergence term is given by

[TABLE]

4.3.2 Integrable spin-1 chains with helix states

The local Hamiltonian (93) is a special case of the family of spin-1 chains surveyed in [40]. For general parameter values the Hamiltonian built from the local Hamiltonians (93) is not integrable which proves that the phenomenon of ballistic transport in the helix state is not related to integrability. However, on a submanifold in parameter space one can identify two integrable families which are special cases of the $U_{q}[\mathfrak{sl}(2)]$-symmetric Hamiltonian [41]

[TABLE]

where

[TABLE]

with the spin-1 representation of $SU(2)$ and deformation parameter $q=\mathrm{e}^{i\lambda}$.

Comparing coefficients one finds

[TABLE]

which is the integrable Zamolodchikov-Fateev Hamiltonian [42]. Moreover, one has

[TABLE]

which is the bi-quadratic Hamiltonian of [43, 44]. It is remarkable that there is no significant difference in the properties of the helix states for the integrable and the non-integrable cases. The integrable models, however, are of particular interest as they allow for a more detailed study, including transport properties in the pure quantum case and possibly the construction of non-local conserved quantities that are relevant for the derivation of transport properties of these models [26].

5 Concluding remarks

We have defined a family of spin helix states (SHS) with twist angle $\varphi$ in the $xy$-plane between neighboring spins and shown that these states arise as the exact stationary solution of open spin-1 quantum chains with bulk conservation of the $z$-component of the magnetization, but boundary dissipation given by a suitably chosen two-parameter families of Lindblad operators. These helix states are not in any sense close to the quantum ground states of these spin chains. Nevertheless, they are stationary under the Lindblad boundary driving that targets the boundary spins in different directions, with a boundary twist angle $\Phi=(N-1)\varphi\mod{}2\pi$. A non-zero winding number $K$ determined by $\varphi=(\Phi+2\pi K)/(N-1)$ allows for a stationary spin-current $j^{z}$ of order 1.

Specifically, for the spin-1/2 Heisenberg chain with anisotropy parameter $\Delta=\cos(\eta)$ the SHS occurs when $\eta=\varphi$. As a function of $\eta$ the stationary current $j^{z}$ for fixed $\varphi$ shows a resonance-like peak at the SHS value $\eta=\varphi$. If this matching condition is satisfied then for any fixed anisotropy parameter $\Delta=\cos(\eta)$ the SHS carries a spin current $j^{z}=J\sin{(\eta)}$. This corresponds to ballistic transport, i.e., the current does not depend on system size, since for any $N$ one can find a boundary twist angle $\Phi\in[0,2\pi[$ that supports this current. In fact, even when the boundary twist $\Phi$ is zero the SHS carries a current of order 1 at anisotropies of the form $\Delta=\cos{2\pi K/(N-1)}$. This is reminiscent of a result for the XXZ-chain with different Lindblad operator where the Drude weight has peaks at anisotropies $\Delta=\cos{2\pi m/n}$ ($m,n$ being integers), leading to an overall fractal behaviour of the Drude weight as a function of $\Delta$ in the thermodynamic limit $N\to\infty$ [26]. Whether this Drude weight is related to an SHS is an open question.

We generalized the construction to higher spins. For spin 1 we have derived Hamiltonians which allow for the existence of stationary spin-1 SHS under suitable dissipative dynamics at the boundaries. There Hamiltonians include the integrable Zamolodchikov-Fateev chain [42] and also the bi-quadratic Hamiltonian of [43, 44]. We stress, however, that the existence of SHS is not in any way related to integrability. Our solution includes non-integrable spin chains. Moreover, since the construction relies on a local divergence condition when applying the local Hamiltonian on the SHS, it can be generalized to any lattice that allows for the cancellation of all these terms in the sum of the local Hamiltonians over the lattice. So, in particular, one can construct SHS for two- and three-dimensional cubic lattices. By the same token, we expect that one can generalize the approach to Hamiltonians with next-nearest neighour interactions and to Hamiltonians with valence-bond eigenstates.

Generally, the properties of the SHS show, by comparing with known results for other boundary driving mechanisms, that the transport properties of spin chains depend qualitatively on the choice of Lindblad operators. This is somewhat puzzling as the ballistic or other superdiffusive transport is expected to be a bulk property of the chain, not a boundary property. This is reminiscent of boundary-induced phase transitions in classical stochastic particle systems [45, 46]. Whether there is a deeper link is a further open question.

Acknowledgements

Financial support by DFG is gratefully acknowledged. GMS thanks F.C. Alcaraz and D. Karevski for stimulating discussions and the University of São Paulo and the University of Lorraine for kind hospitality.