Transport and Energetic Properties of a Ring of Interacting Spins Coupled to Heat Baths

We study the heat and spin transport properties in a ring of interacting spins coupled to heat baths at different temperatures. We show that interactions, by inducing avoided crossings, can be a means to tune both the total heat current flowing between the ring and the baths, and the way it flows through the system. In particular, we recognize three regimes in which the heat current flows clockwise, counterclockwise, and in parallel. The temperature bias between the baths also induces a spin current within the ring, whose direction and magnitude can be tuned by the interaction. Lastly, we show how the ergotropy of the nonequilibrium steady state can increase significantly near the avoided crossings.


Introduction
Understanding transport properties in quantum systems can lead to various interesting applications such as quantum rectifiers [1][2][3][4][5][6][7][8][9][10][11][12][13], transistors [14], engines [15][16][17][18][19][20][21][22][23][24][25][26][27][28][29][30][31], refrigerators [32][33][34][35][36][37], and batteries [38][39][40]. In thermoelectric systems, one main goal is to achieve an electric current due to a temperature bias. It is thus important to study fundamental aspects of the conversion of heat currents into particle/spin currents. In the following, we will consider a circuit made of a ring of spins which we couple to baths at different temperatures to study whether a spin current can be induced by them and how this current can qualitatively change depending on the system parameters. We will evaluate the maximum amount of work that the steady state can produce by computing its ergotropy [41]. In particular, we focus on the interplay between interactions in the system and the external driving due to the two different thermal baths. To find more generic results, we focus on a minimal model as the study of minimal models allows one to extract the key ingredients necessary to obtain a given physical phenomenon.
The minimal model we consider is a ring of four spin 1/2 with XXZ interaction (this can be mapped to spinless fermions with nearest neighbor interaction via Jordan-Wigner transformation). The ring is coupled at two opposite sites to two baths at different temperatures, as shown in Figure 1a.
In order for the current to have a preferred direction, we need to break the symmetries of the model, which can be done by applying different local magnetic fields. The ring only exchanges energy in the form of heat with the baths and there is no transfer of spins between them. In this way, we can clearly analyze whether temperature biases can generate a spin current.
There can be two types of currents in the system: the heat current and the spin current. As pictorially represented in Figure 1b, we will show that heat currents can flow in three different ways: the heat can flow in a parallel manner (yellow arrows) in the upper and lower part of the ring, or it can flow in an anti-parallel fashion. In the latter case, the flow can be either clockwise (green arrows) or counterclockwise (red arrows). As for the spin current, since no spin is exchanged with the baths, it can only flow clockwise or counterclockwise. We will show that, in our setup, the spin current inversion occurs together with a significant change in the heat exchanged with the baths. They both occur in the proximity of an avoided crossing which is induced by the interaction between the spins. Since the steady state is not in equilibrium, we also investigate the maximum energy that can be extracted from it via unitary processes, i.e., the ergotropy [41]. We show what is the main contributor to the ergotropy of the steady state and also find that it is significantly enhanced near interaction induced avoided crossings.
The article is organized as follows: in Section 2, we briefly describe the interacting spin ring model coupled to heat baths. In Section 3, we discuss the Markovian Redfield master equation which we use to derive our results. In Section 4, we study how the interaction can result in different scenarios of local heat and spin currents and, in Section 5, we analyze its effects on the ergotropy. In Section 6, we draw our conclusions. Figure 1. (a) schematic representation of a ring of interacting spins with local magnetic fields coupled to two heat baths; (b) three possible heat current modes in the system: clockwise, parallel, and counterclockwise. Since there is no spin exchange with the baths, the spin current can only flow in clockwise or counterclockwise directions.

Model
We aim to study the effect of interactions on both the heat and spin currents. We thus consider a prototypical model, i.e., an XXZ chain with periodic boundary conditions. The Hamiltonian is given byĤ whereσ u l = (σ z l + 1)/2 and where site l = 5 coincides with site one (We have written the interaction asσ u lσ u l so that, when converting the spins to fermionic particles, this interaction term corresponds to only nearest neighbor interaction without a ∆-dependent local chemical potential). Theσ a l with a = x, y, z are the operators corresponding to the Pauli matrices.
In order to generate internal currents in the system, it is necessary to break the reflection and rotational symmetries. This can be done by applying an external perturbation, or disorder, in terms of local magnetic field. The minimum type of perturbation needed to break the symmetries is to apply a different local magnetic field to two consecutive spins. Here, local magnetic fields h z l are applied to the second and third spins as depicted in Figure 1a. We seth = k B = J = 1 throughout.
The two spins at sites l = 1, 3 are coupled to two heat baths as shown in Figure 1a. The left and right heat baths are infinite collections of harmonic oscillators with their respective HamiltoniansĤ L andĤ R given byĤ We consider baths characterized by the spectral density given by J(ω) = γω exp (−ω/ω c ), where ω c is the cutoff frequency and γ is the system-bath coupling strength. We use a cutoff of ω c = 10 for both baths. The system-bath coupling strength γ is chosen as 0.01, which is sufficiently weak as compared to the level spacings, even in the presence of avoided crossings.
The ring and the baths interact via system operatorsσ z 1 andσ z 3 aŝ Since the operatorσ z l conserves the total number of spins in the system, there is no spin current between the baths and the system. However, as we show later, the heat current can induce a spin current within the system in a direction that depends on the strength of the interaction.

Minimality of the Model
We have considered a circuit of four spins, and this is the minimal number required to observe the effects of interaction on the ring current. With two spins, it is not possible to form a close circuit. With three spins, it is possible to form the smallest circuit; however, in such a system, the interaction term would be equivalent to a global energy shift. This can be shown quite simply: since the total magnetization is conserved, for three spins, it is only possible to have a current only with one spin up and two spins down or vice versa. However, in either case, the number of possible configurations is three and in all configurations, the number of parallel and anti-parallel neighboring spins is invariant. Hence, the interaction is effectively a homogeneous local potential. With four spins instead, in the magnetization sector with two spins up and two spins down, it is possible to have configurations in which the number of parallel or anti-parallel neighboring spins is different, e.g., alternating spins up and down or two spins up followed by two spins down. Hence, four spins is the smallest size of a circuit which can be used to explore the effects of interactions.

Methods
To study the nonequilibrium properties of the steady state, we use the Redfield master equation [42] which can be derived from Equations (1)-(3). The Redfield master equation is a second-order perturbative master equation that produces accurate results for the occupation of the energy levels at zeroth order in system-bath coupling. The off-diagonal elements in the energy eigenbasis are, on the other hand, accurate up to the second order [43][44][45]. One important advantage of using this master equation is that it does not require the secular approximation which can result in the vanishing local currents within the system [46][47][48][49]. The Redfield master equation is often criticized for producing negative probabilities [50]. We emphasize here that a strong Markovian environment, i.e., a fast decaying correlation function, usually helps to avoid such a problem. Moreover, the emergence of negative probabilities is a clear warning that the master equation is being used beyond its regime of validity (a check which is not present in master equations in Gorini-Kossakowski-Sudarshan-Lindblad (GKSL) form [51,52]). In the simulations presented here, we have not observed negative probabilities.
To compute the nonequilibrium steady-state (NESS) density operatorρ NESS =ρ(t = ∞), we use the Redfield master equation given by [42,53], where the dissipative part R l are and where, for our setup,Ŝ L =σ z 1 andŜ R =σ z 3 , and we associate the subindex L, of the left bath, with site l = 1 and subindex R, of the right bath, with site l = 3. Since the upper limit of the integral in Equation (6) is τ = ∞, Equation (4) is also referred to as the Markovian Redfield master equation. Note that the action of the baths on the system is considered to be additive, which is in general a good approximation for sufficiently weak coupling and Markovian baths [54]. It should also be remarked that this approach considers the full Hamiltonian of the systemĤ S in derivingŜ l , and not a local approximation of it. This is important because the local Hamiltonian approximation can result in a failure to capture the dependence on many-body interactions [55] and can be thermodynamically inconsistent [56]. In Equations (5) and (6), R l [·] is the dissipator that contains all the bath information through the transition operatorŜ l . The bath correlator C l (τ) is explicitly given by where J(ω) is the spectral density specified in Section 2 [53] and T l is the temperature in either the left (T L ) or right (T R ) bath (see also Figure 1). The Markovian Redfield master equation in the energy eigenbasis results in the following equation for the density matrix ρ α,β , where ∆ α,β = E α − E β is the energy difference between the energy levels α and β, while R α ,β l,α,β is a tensor acting on the density matrix ρ α,β given by with δ α,β the Kronecker delta. Here, S † l,β,β and S † l,β,β are derived, respectively, from the coefficients of the operatorsŜ l andŜ l in the energy eigenbasis.
From Equation (6), we get a transition matrix in the energy basis given by which can be evaluated using Plemelj formula ∞ 0 dte ±i t = πδ ( ) ± i P . The imaginary part, usually known as the Lamb shift, is a small perturbation to the system Hamiltonian and it is therefore neglected.
In case the secular approximation had been taken, as for master equations in GKSL form, then all the terms in R α ,β l,α,β in Equation (9), with α = α and β = β will be omitted. As a consequence, the resultant master equation would give vanishing off-diagonal elements for the reduced density matrix [48,57].

Local Current Modes
We are interested in the local heat and spin currents within the system. To obtain the expressions for the spin and heat currents, we use the conservation laws of, respectively, the local spinσ u l and the local energy on a bondĥ l,l+1 wherê Hence, the local heat current operator on the l-th site, except for the sites in contact with the baths, is given byĵ Similarly, for the local spin current operator, we havê which is uniform through the system because the total magnetization is conserved, and it can only be in the clockwise or counterclockwise direction. In our setup, the heat current in the upper part of the ring is given by j H 2 = ĵ H 2 , while the one in the lower part of the ring is given by j H 4 = ĵ H 4 . Here and henceforth, we use the notation · = Tr ( ·ρ NESS ) to indicate the trace over the steady-state density operator of the systemρ NESS . The total heat current is given by the difference of the current in the upper and lower part of the ring (Current operators are defined as positive in the clockwise direction), and it is also given by the energy exchanged with each bath, which is computed by j H L = Tr Ĥ S R L [ρ NESS ] , i.e., the energy drawn from the left bath, or equivalently from j H R = Tr Ĥ S R R [ρ NESS ] , i.e., the energy drawn from the right bath. Indeed, we have that j H L = −j H R = j H 2 − j H 4 . We use j S = ĵ S l to denote the local spin current. Due to the conservation of total magnetization (number of spins up), the system has multiple invariant subspaces. When there is zero or one spin up, the interaction ∆ does not play any role in the system. When there are three or four spin ups, the interaction ∆ in Equation (1) acts effectively as local fields. For two spins up, ∆ differentiates cases in which the two spins are next to each other or not, hence acting as a nearest neighbor interaction. We thus study the system in the symmetry sector with two spins up and contrast it with results from the sector with three spins up.
In Figure 2a, we consider the case with three spins up. The local heat currents j H 2 and j H 4 are depicted by the orange dashed and green dot-dashed lines, respectively. While the total heat current (red solid line) remains unchanged when ∆ increases, the local heat currents change linearly with ∆ because of the linear shift of the eigenenergies. Because of this linear change, the local heat currents geometry changes from counterclockwise (pink shaded region in Figure 2) to parallel (yellow shaded region) and finally to clockwise (green shaded region). The local spin currents, on the other hand, do not vary with respect to ∆ (see Figure 2c). This is due to the fact that, in this sector with three spins up, ∆ α,β in Equations (8) and (10) is independent from ∆ and henceρ NESS is also invariant with it. This explains why the total heat current and the local spin current do not change. However, the local heat current operator depends on ∆ via Equations (11) and (12), and hence j H changes with ∆.
In Figure 2b,d we show the effect of the interaction ∆ on the currents in the sector with two spins up. In this case, the local heat currents also demonstrate a transition between the three different geometries, from parallel (yellow shaded regime), to clockwise direction (green shaded regime), and to counterclockwise direction (red shaded region). However, the total heat currents vary in a non-monotonous way, with a minimum at ∆ ≈ 2.19. This minimum corresponds to a sharp change in the geometry of the local heat currents and the local spin current which experience an interaction induced inversion of the direction, from clockwise to counterclockwise. This current inversion occurs because of the presence of a small avoided crossing near ∆ ≈ 2.19 which connects two different energy eigenstates at sufficiently low energy such that they have a large enough weight to play a significant role in the properties of the system. We will discuss in more detail the role of avoided crossings in the next section.

Ergotropy
We now show how interactions can significantly affect the possibility of extracting energy from the system via unitary processes U. In this section, we focus on the sector with two spins up. The maximum energy that can be extracted via a unitary process is quantified by the ergotropy E [41], which is calculated for the steady state using whereρ passive is the corresponding passive state [58] which is built from the eigenvalues p k ofρ NESS , and from the system Hamiltonian in its eigenbasisĤ More precisely,ρ which means that higher energy levels are less populated. Given a state, the passive state is the lowest energy state that can be reached via unitary transformations. Examples of passive states are thermal states, from which it is not possible to extract energy (work) with unitary operations.
In Figure 3a, we study the ergotropy extraction within the system due to the temperature bias and the interaction strength. Increasing temperature bias, ∆T = T L − T R , drives the system more out of equilibrium while interaction can be used to tune the energy level differences. For illustration purposes, we consider a fixed average temperature bias (T L + T R )/2 = 1.5. Various hot spots for ergotropy are identified far from equilibrium as shown in Figure 3a. Note that the ergotropy is not symmetric with respect to the interaction ∆ or temperature bias ∆T, indicating a rectification of ergotropy in the system. In general, ergotropy can be non-zero for two main reasons: an occupation of higher energy levels due to a strong temperature bias, or the presence of coherence in the nonequilibrium steady state. From Figure 3b,c we observe that in our setup the main cause is the change in the occupation of the energy levels, which in the following we refer to as "population inversion". We note that this is expected within the weak system-bath coupling limit and our perturbative approach which can accurately study systems with small off-diagonal terms in the energy eigenbasis [43][44][45]. Figure 3b,c shows the occupations of different energy levels ρ k,k (where k increases for increasing energy) versus ∆. When comparing Figure 3b with the line ∆T = 1 in Figure 3a, or Figure 3c with the line ∆T = −1 in Figure 3a, we observe that the ergotropy maxima correspond to the regions with population inversion, highlighted by the light blue shadings in Figure 3b,c. In particular, in Figure 3b, there is a cross over between ρ 2,2 (orange line) and ρ 3,3 (green line) for larger ∆ and between ρ 4,4 (red line) and ρ 5,5 (purple line) for lower ∆. In Figure 3c, the population inversion is between ρ 3,3 (green line) and ρ 4,4 (red line). This confirms that the leading contribution to ergotropy is from population inversion.
In the following, we show that the population inversion occurs close to avoided crossings due to the presence of smaller energy gaps (We emphasize that the avoided crossings and the ergotropy maxima do not exactly coincide due to small corrections from other terms in the master Equation (8)).
This gives a mechanism to use interactions to tune the amount of ergotropy in the system. In Figure 4a,b, we show the ergotropy E and the total heat current j H versus ∆ for different values of the local magnetic fields h z l . For both E and j H , we observe peaks which shift due to the interplay between the interaction and the local magnetic fields. With Figure 4c, we can associate these peaks with the avoided crossings. In Figure 4c, we show the difference in energy between levels E 3 and E 2 (orange lines), levels E 4 and E 3 (purple lines) and levels E 5 and E 4 (blue lines) versus ∆ and for the same local magnetic fields used in Figure 4a,b. It is clear that the maxima of ergotropy and the minima of heat current occur together with the small avoided crossings. Moreover, Figure 3b shows that avoided crossings in low energy states (ρ 3,3 and ρ 4,4 ) result in larger ergotropy as these levels are more occupied than those with higher energy.  (c) energy difference as a function of interaction ∆ for low energy states E 3 − E 2 (orange lines), intermediate energy states E 4 − E 3 (purple lines), and high energy states E 5 − E 4 (blue lines). In (a-c), color gradient indicates the strength of the local magnetic field. h z 2 /h z 3 = 2 where h z 2 = 2, 1.8, 1.6, 1.4, 1.2 for increasingly dark colors. T L = 2, and T R = 1. The system-bath coupling γ = 0.01 for both baths, and we are considering the sector with two spins up.

Conclusions
We have considered a minimal system to study the interplay between heat and spin currents. In particular, our setup is composed of four spins in a ring which is connected to two baths at different temperatures. Each spin is interacting with its nearest neighbor and a position-dependent magnetic field breaks the reflection and rotation symmetries of the system. This setup could be implemented with ultracold ions realizations of spin-1/2 systems [59][60][61], or with quantum dots circuits [62][63][64]. Even without interactions, it is possible for the heat current in the system to flow in three different ways: clockwise, counterclockwise, and in parallel flows. However, we show here that interactions give us a means to cause sizeable changes in the total heat current and induce an inversion in the spin current. When coupling the system with an ancilla, such changes in internal currents could be used for sensing or activation. The internal currents can also be used as indicators of quantum criticality [65,66]. Coupling such a system to two baths brings it to a nonequilibrium steady state from which it is possible to extract work with unitary processes, i.e., the system can have non-zero ergotropy. We show that the ergotropy is, in general, non-symmetric with respect to the temperature bias ∆T, and it can be significantly enhanced close to avoided crossings due to population inversions.
For larger rings, the number of relevant avoided crossings can increase; however, we expect the same qualitative behavior. Systems with more complex topologies, e.g., two or more rings, can present a larger variety of behaviors which is worth investigating further. In future works, similarly to [20], we may consider the effect of time-independent and time-dependent gauge fields and consider the thermodynamic properties of such systems. In our setup, the weak coupling between the system and the baths limits the contribution of coherence to the ergotropy. It would thus be interesting to study the effects of stronger system-bath coupling on heat currents and ergotropy.