# Heat current magnification in Classical and Quantum spin networks

Vipul Upadhyay<EMAIL_ADDRESS>Department of Physics, Indian
Institute of Technology Delhi, Hauz Khas 110 016, INDIA  Poshika Gandhi
Institute of Physics, University of Freiburg, Hermann-Herder-Str. 3, 79104
Freiburg, Germany  Rohit Juneja  Department of Physics, Indian Institute of
Technology Delhi, Hauz Khas 110 016, INDIA  Rahul Marathe
<EMAIL_ADDRESS>Department of Physics, Indian Institute of
Technology Delhi, Hauz Khas 110 016, INDIA

###### Abstract

We investigate heat current magnification due to asymmetry in the number of
spins in two-branched classical and quantum spin systems. We begin by studying
the classical Ising like spin models using Q2R and CCA dynamics and show that
just the difference in the number of spins is not enough and some other source
of asymmetry is required to observe heat current magnification. Unequal
spin–spin interaction strength in the upper and lower branch is employed as a
source of this asymmetry and it proves adequate for generating current
magnification in both the models. Suitable physical motivation is then
provided for current magnification in these systems, along with ways to
control and manipulate magnification through various system parameters. We
also study a five spin Quantum system with modified Heisenberg XXZ interaction
and preserved magnetisation using the Redfield master equation. We show that
it is possible to generate current magnification in this model by the
asymmetry in the number of spins only. Our results indicate that the onset of
current magnification is accompanied by a dip in the total current flowing
through the system. On analysis it is revealed that this dip might occur
because of the intersection of two non-degenerate energy levels for certain
values of the asymmetry parameter in the modified XXZ model. We deduce that
the additional degeneracy and the ergodic constraint due to fixed
magnetisation in the system are the main reasons for current magnification and
other atypical behaviors observed. We then use the concept of ‘ergotropy’ to
support these findings. Finally, for both the classical and quantum models, we
see that current magnification is only observed when temperature gradient and
intra-system interaction strength have similar order of energy.

##  I Introduction

Heat currents can influence systems in a significant manner at the nano-scale
and as a result, their control and manipulation at small scale is currently a
subject of intense research in both the classical as well as quantum domain.
These studies may be broadly divided into two categories, one dealing with
finding ways to control heat currents by modelling and manufacturing of small
scale thermal devices like thermal diodes Kaushik  _et al._ (2018); Upadhyay
_et al._ (2021); Landi  _et al._ (2014); Jiang  _et al._ (2010); Ordonez-
Miranda  _et al._ (2017); Segal and Nitzan (2005); Karg ı  _et al._ (2019);
Segal (2008); Zhang  _et al._ (2013); Terraneo  _et al._ (2002), thermal
transistors Saira  _et al._ (2007); Chung Lo  _et al._ (2008); Ruan and Wang
(2020); Naseem  _et al._ (2020); Jiang  _et al._ (2015); Zhang  _et al._
(2018); Sánchez  _et al._ (2017); Messina  _et al._ (2012); Hwang  _et al._
(2009) etc., and the other more direct way of explicitly studying heat current
transport and distribution in various systems, and developing transportation
theories Dhar (2008); Poulsen  _et al._ (2022a, b); Chioquetta  _et al._
(2021); Li  _et al._ (2004); Cho and McKenzie (2005a); Lepri  _et al._ (2003);
Li  _et al._ (2012); Pereira (2017); Kosloff (2013); Kirchberg and Nitzan
(2022). Interestingly, for a multi branched system kept between two heat baths
at different temperatures, it is observed that under certain conditions heat
current in one of the branches may become larger than the total current
flowing between the baths, thus leading to the phenomenon called the current
magnification or circular current Xu  _et al._ (2019); Dugar and Chien (2022);
Roy (2007); Rai  _et al._ (2010); Marathe  _et al._ (2010); Bandopadhyay  _et
al._ (2004); Jayannavar  _et al._ (1995); Cho and McKenzie (2005b); Xiao  _et
al._ (2013); Nakanishi and Tsukada (2001); Eckmann and Zabey (2004); Roy
(2007). Current magnification has been shown to exist in a wide variety of
physical systems like spins Xu  _et al._ (2019), molecules Rai  _et al._
(2010); Nakanishi and Tsukada (2001), classical harmonic chains Marathe  _et
al._ (2010); Eckmann and Zabey (2004), and metallic rings Jayannavar and
Singha Deo (1995); Roy (2007). In general, it is observed that to get current
magnification, the system must possess either rotational or reflection
asymmetry Xu  _et al._ (2019). For the simplest case of a system containing
two branches, the heat current can flow in three possible ways, namely,
parallel currents in the branches, clockwise circulating current, and
anti–clockwise circulating current (see Fig. 1).

Figure 1: Possible ways for current to flow in a two-branch system (motivated
from Xu  _et al._ (2019))

In a recent study Xu  _et al._ (2019), it was observed that it is possible to
get circulating heat currents in Quantum spin systems with modified Heisenberg
exchange interaction if the onsite magnetic field is inhomogeneous and the
total magnetisation of the system is conserved. However, it is difficult to
realise this system experimentally as applying different magnetic fields at
different sites for such small systems is in general a difficult task.
Motivated by this and earlier studies on heat current magnification, here we
study some simple classical and quantum models for heat current magnification
in two-branched spin systems. In particular, we study the heat current flow in
the classical Creutz Cellular automaton (CCA) and Q2R models described in
Vichniac (1984); Saito  _et al._ (1999); Creutz (1986); Schulte  _et al._
(1987); Stauffer (2000). We consider systems with unequal branch spin numbers
and with equal as well as unequal branch spin–spin interaction strengths. We
find the required conditions for getting current magnification, and the range
of parameters that optimise it. We also study a five spin Quantum system
similar to Xu  _et al._ (2019) but with spin number asymmetry in two branches.
We use the Redfield master equation, and explore the heat flow and energetics
in it. Similar to the classical models we find the parameter range suitable
for getting current magnification in this system. The models that are studied
in this manuscript can have possible experimental realizations for example in
systems like quasi-1D Ising chains Coldea  _et al._ (2010); Wolf (2000) or
quantum simulators like NMR Micadei  _et al._ (2019); Zhang  _et al._ (2009),
quantum dots Seo  _et al._ (2013); Rogge and Haug (2008), trapped ions
Bermudez  _et al._ (2013); Porras and Cirac (2004).

The paper is organized as follows. In Sec. II we consider the classical spin
models and provide analytical and numerical results. In Sec. III we discuss
the numerical Pollard and Friesner (1994) implementation of branch spin number
asymmetry in the Quantum case by analysing our system by the Redfield master
equation and highlight the relative ease for achieving current magnification
in our model. We summarise our findings and conclude in Sec. IV.

##  II Model and Classical Analysis

Figure 2: Schematic representation of the model system

We consider an Ising-like two-state spin chain with periodic boundary
conditions as shown in Fig. 2. This system is kept between two heat baths at
temperatures $T_{L}$ and $T_{R}$ to its left and right respectively. The
system interacts with the heat baths at two individual nodes such that only
the nodal spins are in direct contact with them while the bulk spins are not.
This results in the formation of two spin branches between the baths as shown
in Fig. 2. We use $N_{U}(N_{D})$ for the number of spins in the upper (lower)
branch which are not in contact with the heat baths. So the total number of
spins is given by $N_{S}=N_{U}+N_{D}+2$. The spins are numbered in clockwise
direction, i.e., starting from the spin in contact with the left bath followed
by the upper branch spins and then the lower branch spins implying that the
spins numbered $`1^{\prime}$ and $`N_{U}+2^{\prime}$ interact with the left
and right bath respectively. The state or configuration of our system is
described as:

$\eta=\\{\sigma_{1},\sigma_{2},\ldots,\sigma_{N_{u}+2},\ldots\sigma_{N_{S}}\\},$
(1)

where $\sigma_{i}$ denotes the state of $i^{th}$ spin and
$\sigma_{i}\in\\{-1,1\\}$. The energy corresponding to a particular
configuration is decided by the Hamiltonian of the system. For the cases that
we study in this paper, only the nearest neighbor spin–spin interaction are
considered. Thus the Hamiltonian is given by

$\displaystyle H^{I}_{S}=-\sum_{i}J_{i}\sigma_{i}\sigma_{i+1},$ (2)

where $H^{I}_{S}$ is the spin–spin interaction term in the total Hamiltonian,
$J_{i}$ is the interaction strength between the $i^{th}$ and $(i+1)^{th}$ spin
and $\sigma_{N_{S}+1}=\sigma_{1}$. This means that the energy cost for
flipping a spin $\sigma_{i}$ is,

$\displaystyle\Delta
E_{i}=2J_{i}\sigma_{i}\sigma_{i+1}+2J_{i-1}\sigma_{i}\sigma_{i-1},$ (3)

subject to the periodic boundary conditions. Due to the interaction with the
baths as well as the internal interaction between the spins, the configuration
of the system may change resulting in the flow of heat currents in our system.

###  II.1 Time evolution dynamics

To carry our analysis forward, we first need to define the appropriate time
evolution dynamics for our system. The master equation for our system is
characterised by configuration pairs $(\eta,\eta^{i})$ where $\eta^{i}$ is
same as $\eta$ except for the flip of $i^{th}$ spin and is given as:

$\displaystyle\frac{dP(\eta)}{dt}=\sum_{i=1}^{Ns}r_{\eta^{i}\to\eta}P(\eta^{i})-\sum_{i=1}^{Ns}r_{\eta\to\eta^{i}}P(\eta),$
(4)

where $P(\eta)$ is the probability of finding the system in the configuration
$\eta$ at time $t$ and $r_{\eta^{i}\to\eta}$ is the transition rate for
flipping the spin $i$ and taking the system from configuration $\eta^{i}$ to
$\eta$. As already discussed above, for our model the nodal spins are in
direct contact with the baths while the bulk spins are not. This suggests that
to study our system we need hybrid dynamics with separate transition rates for
nodal and bulk spins.

Firstly, because of their interaction with the baths, the dynamics of the
nodal spins can be studied using the Metropolis algorithm, for which a spin
flips according to the following transition rate Kaushik  _et al._ (2018)

$\displaystyle
r_{\eta\to\eta^{i}}^{\alpha}=\text{min}(1,~{}e^{-{\beta_{\alpha}\Delta
E_{i}}}),$ (5)

where $\beta_{\alpha}=1/{k_{B}T_{\alpha}}$, $\alpha\in\\{L,R\\}$,
$i\in\\{1,N_{U}+2\\}$ and we work in units where $k_{B}=1$. Since the bulk
spins are not in contact with any baths, a dynamics involving random numbers
is not suitable for them and a deterministic energy-conserving dynamics is
needed for understanding their time evolution. Such kind of dynamics generally
fall under Cellular automaton Wolfram (1983); Vichniac (1984). We use two
different types of cellular automaton for our problem, first the Q2R dynamics
in the subsection II.4, and then the CCA dynamics in II.7. We now discuss the
current difinitions required for the classical analysis.

###  II.2 Heat Current Definitions

To define the heat current flowing out of a bath, we look at the nodal spin
associated with it and identify all the configuration pairs $(\eta,\eta^{i})$
which differ by its flip. The corresponding heat energy flow is then given as:

$\displaystyle I_{\alpha}=\sum_{(\eta,\eta^{i})}\Delta
E_{i}~{}(r_{\eta\to\eta^{i}}^{\alpha}P(\eta)-r_{\eta^{i}\to\eta}^{\alpha}P(\eta^{i})),$
(6)

where the sum runs over all the configuration pairs that differ by a flip of
the $i^{th}$ spin with $\alpha\in\\{L,R\\}$ for $i\in\\{1,N_{U}+2\\}$. Since
the bulk spins are not in contact with the baths, their flipping should cause
no change in total energy of the system for any of the dynamics used. As a
result, the Eq. (6) gives us zero when applied for flipping of bulk spins.
However, the flip of bulk spins still leads to a transfer of energy from one
spin to other over a bond Saito  _et al._ (1999). To better understand this,
let us consider the interaction energy associated with a spin $i$,

$\displaystyle
E_{i}=-J_{i}\sigma_{i}\sigma_{i+1}-J_{i-1}\sigma_{i}\sigma_{i-1}.$ (7)

The above expression contains the two bond energy terms associated with the
spin $i$. Since this spin does not interact with the baths, it can flip only
when the associated energy cost given in Eq. (3) is either zero or can be
compensated by some other non-interacting internal energy of the system. In
both these cases, the individual bond energies change, signifying a transfer
of energy through the bonds and hence a heat current. Focusing on the
$`i\text{--}i+1^{\prime}$ bond, the energy transferred through it on the
flipping of $i^{th}$ spin is:

$\displaystyle\Delta E_{i,i+1}=2J_{i}\sigma_{i}\sigma_{i+1}.$ (8)

To find the associated heat current, we identify all the configuration pair
transitions related by the flip of $i^{th}$ spin but focus only on the
corresponding energy change on the $`i\text{--}i+1^{\prime}$ bond. This gives
the following expression for heat current passing through bulk spin $i$

$\displaystyle I_{U(D)}=\sum_{(\eta,\eta^{i})}\Delta
E_{i,i+1}(r_{\eta\to\eta^{i}}P(\eta)-r_{\eta^{i}\to\eta}P(\eta^{i})),$ (9)

where $I_{U(D)}$ is the heat current in the upper (lower) branch depending on
the position of the $i^{th}$ spin. Since, the flipping of the bulk spins is
deterministic, it will either definitely happen or is completely forbidden and
as a result $r_{\eta\to\eta^{i}}\in\\{0,1\\}$, unlike the stochastic flipping
of the nodal spins. Finally, current magnification occurs when the heat
current in one of the branches is larger than the current flowing between the
system and the bath. According to the convention used in the current
definitions above, we get current magnification if both the branch currents
have the same sign.

###  II.3 Information about Methodology

We now discuss the methodology Sandvik (2010) used to find results for the
classical case. All the results shown in this manuscript are for the steady
state. The simulation for the classical case are performed by randomly
selecting a spin and flipping it depending on a suitable dynamics. Here, we
use a hybrid dynamics involving the metropolis algorithm for the nodal spins
and the Q2R or CCA dynamics for the bulk spins. Whichever the case, only one
spin is allowed to flip per time step. The system is allowed to relax to the
steady state and the required quantities are then calculated. The
corresponding numerical and analytical calculations are done by writing the
Master equation and solving the transition matrix like the one given in Eq.
(11). Since, we are only interested in the steady state calculations, we
require the eigenvector corresponding to the eigenvalue $`0^{\prime}$ of the
transition matrix. On applying the probability normalisation condition on this
eigenvector, steady state probabilities can be obtained. These can then be
used to calculate the currents analytically.

Note on experimental parameter range: In the units that we work, for both the
classical and the quantum systems, the effect on system due to the interaction
with the baths is quantified by an expression of the type $e^{J/T}$. This
means that for any meaningful dynamics to take place, the temperature of the
baths should correspond to the thermal energy of the same order as $J$. The
value of $J$ is typically decided by the experimental setup being studied. For
example the experimental setup studied by Coldea et al. in Coldea  _et al._
(2010), the value of $J$ is in few mev which roughly corresponds to the
thermal energy associated with temperatures of order $10$K.

Figure 3:  Illustration of the model system consisting of 3 spins with
$N_{U}=1$ and $N_{D}=0$.

###  II.4 Q2R dynamics with symmetric upper and lower branch interaction
strength

We start our analysis by studying the Q2R dynamics Saito  _et al._ (1999);
Creutz (1986) for the flipping of the bulk spins. This is one of the simplest
model which allows us to have a deterministic mechanism for time evolution of
the bulk spins while preserving the total energy of the system. The
Hamiltonian for this model is the same as the typical Ising Hamiltonian and is
obtained by substituting $J_{i}=J$ in equation (2). According to this
dynamics, a bulk spin $i$ can flip only if its neighboring spins have opposite
orientation i.e.,

$\displaystyle\sigma_{i-1}=-\sigma_{i+1}.$ (10)

This ensures that the total energy change of the system in Eq. (3) is zero on
flipping of the bulk spin.

Figure 4: (Color Online) (a) Comparison of total current ($I_{L}$) between
analytical, numerical and simulation results for the 3 spin system with
$N_{U}=1$, $N_{D}=0$, (b) Variation of Heat Currents ($I$) with $T_{L}$ for
$N_{U}=1$, $N_{D}=0$, (c) Variation of $I_{L}$ with $T_{L}$ for a 8 spin
system but with different branch spin distributions, (d) Variation of Heat
Currents ($I$) with $N_{U}$ for $N_{D}=4,T_{L}=2,T_{R}=0.1$. For all the cases
unless otherwise specified, we use $T_{R}=1$, $J=1$.

To better understand the implications of the Q2R dynamics, we perform an
analytical study for the minimalistic model with just three spins and
$N_{U}=1$ and $N_{D}=0$ as shown in Fig. 3. This being the simplest case,
where the required asymmetric spin branching is possible, offers a good
starting point. The master equation (4) for this system written in the matrix
form is given by,

$\displaystyle\frac{d}{dt}\begin{pmatrix}P(\downarrow,\downarrow,\downarrow)\\\
P(\uparrow,\downarrow,\downarrow)\\\ P(\downarrow,\uparrow,\downarrow)\\\
P(\uparrow,\uparrow,\downarrow)\\\ P(\downarrow,\downarrow,\uparrow)\\\
P(\uparrow,\downarrow,\uparrow)\\\ P(\downarrow,\uparrow,\uparrow)\\\
P(\uparrow,\uparrow,\uparrow)\end{pmatrix}=\begin{pmatrix}-e^{-4J\beta_{L}}-e^{-4J\beta_{R}}&1&0&0&1&0&0&0\\\
e^{-4J\beta_{L}}&-3&0&1&0&1&0&0\\\ 0&0&-2&1&0&0&1&0\\\
0&1&1&-3&0&0&0&e^{-4J\beta_{R}}\\\ e^{-4J\beta_{R}}&0&0&0&-3&1&1&0\\\
0&1&0&0&1&-2&0&0\\\ 0&0&1&0&1&0&-3&{e^{-4J\beta_{L}}}\\\
0&0&0&1&0&0&1&{-e^{-4J\beta_{L}}-e^{-4J\beta_{R}}}\end{pmatrix}\begin{pmatrix}P(\downarrow,\downarrow,\downarrow)\\\
P(\uparrow,\downarrow,\downarrow)\\\ P(\downarrow,\uparrow,\downarrow)\\\
P(\uparrow,\uparrow,\downarrow)\\\ P(\downarrow,\downarrow,\uparrow)\\\
P(\uparrow,\downarrow,\uparrow)\\\ P(\downarrow,\uparrow,\uparrow)\\\
P(\uparrow,\uparrow,\uparrow)\end{pmatrix}.$ (11)

Solving the master equation (as discussed in Sec. II.3), we arrive at the
following steady state probabilities:

$\displaystyle
P(\downarrow,\downarrow,\downarrow)=\frac{e^{4J\beta_{L}+4J\beta_{R}}}{\mathcal{D}_{1}},$
$\displaystyle
P(\uparrow,\downarrow,\downarrow)=\frac{3e^{4J\beta_{L}}+5e^{4J\beta_{R}}}{8\mathcal{D}_{1}},$
$\displaystyle
P(\downarrow,\uparrow,\downarrow)=\frac{e^{4J\beta_{L}}+e^{4J\beta_{R}}}{2\mathcal{D}_{1}},$
$\displaystyle
P(\uparrow,\uparrow,\downarrow)=\frac{5e^{4J\beta_{L}}+3e^{4J\beta_{R}}}{8\mathcal{D}_{1}},$
$\displaystyle
P(\downarrow,\downarrow,\uparrow)=\frac{5e^{4J\beta_{L}}+3e^{4J\beta_{R}}}{8\mathcal{D}_{1}},$
$\displaystyle
P(\uparrow,\downarrow,\uparrow)=\frac{e^{4J\beta_{L}}+e^{4J\beta_{R}}}{2\mathcal{D}_{1}},$
$\displaystyle
P(\downarrow,\uparrow,\uparrow)=\frac{3e^{4J\beta_{L}}+5e^{4J\beta_{R}}}{8\mathcal{D}_{1}},$
$\displaystyle
P(\uparrow,\uparrow,\uparrow)=\frac{e^{4J\beta_{L}+4J\beta_{R}}}{\mathcal{D}_{1}}.$
(12)

where,

$\displaystyle\mathcal{D}_{1}$
$\displaystyle=2e^{2J(\beta_{R}+\beta_{L})}\mathcal{D},$
$\displaystyle~{}\text{with}~{}\mathcal{D}$
$\displaystyle=\left(3\cosh(2J(\beta_{R}-\beta_{L}))~{}+~{}e^{2J(\beta_{L}+\beta_{R})}\right),$
(13)

and the notation $P(\uparrow,\downarrow,\downarrow)$ specifies the probability
of getting the left spin up, middle spin down, and right spin down
respectively. Using the above probabilities and the rates defined in Eq. (5),
the steady state heat current flowing out of the left bath in Eq. (6) becomes,

$\displaystyle I_{L}$
$\displaystyle=\frac{8}{3}J(e^{-4J\beta_{L}}P(\downarrow,\downarrow,\downarrow)-P(\uparrow,\downarrow,\downarrow))$
$\displaystyle=J\sinh(2J(\beta_{R}-\beta_{L}))/\mathcal{D}.$ (14)

Similarly, the branch heat currents are found by using the Eq. (9) for the
bulk spins in upper and lower branch and have the following expressions,

$\displaystyle I_{U}$
$\displaystyle=\frac{J}{3}\sinh(2J(\beta_{R}-\beta_{L}))/\mathcal{D},$
$\displaystyle I_{D}$
$\displaystyle=-\frac{2J}{3}\sinh(2J(\beta_{R}-\beta_{L}))/\mathcal{D}.$ (15)

From the above expressions for total and branch currents, we see that all the
currents go to zero for equal temperatures of baths as expected. For different
bath temperatures, the branch currents add up to give us the total current
i.e. $I_{L}=I_{U}-I_{D}$. The minus sign comes because of the convention used
for defining current direction in the system ($+$ve when current direction is
left to right). Also, the magnitude of the total current is always greater
than the magnitude of the branch currents for any values of system parameters,
so heat current magnification is not possible for this system. Interestingly,
the heat current flowing through a branch is inversely proportional to the
number of spin-spin bonds in that branch. This is similar to Ohm’s law in an
electric circuit with the spin-spin bond behaving analogous to a resistor.

Figure 5:  (Color Online) Variation of Heat Currents ($I$) with (a) $T_{L}$,
with (b) $T_{L}$ for $N_{D}=2,N_{U}=3$, with (c) $T_{R}$, with (d) $J_{2}$,
with (e) $N_{U}$ for $N_{D}=4$ and (f) Variation of Heat Currents ($I$) and
$\langle\sigma_{2}\sigma_{3}\rangle$ with $T_{L}$. For all the cases, unless
otherwise specified we use
$J_{1}=2,T_{L}=2,T_{R}=0.1,J_{2}=1.9,N_{D}=3,N_{U}=2$. In all the figures
where multiple regions of current circulation directions are present, we use
green, gray and yellow background colors to indicate clockwise, parallel and
anticlockwise circulating currents respectively. For all other cases, we use a
white background.

The corresponding simulation results for this model are shown in Fig 4. We see
in Fig. 4$(a)$ that for the $3$-spin system, our analytical and numerical
results match perfectly with the simulation results. In Fig. 4$(b)$, we plot
the current as a function of the temperature of the left bath $T_{L}$, keeping
the temperature of right bath $`T_{R}^{\prime}$ fixed . We find that the heat
currents are non-monotonic functions of temperature difference, initially
increasing sharply and then saturating for large temperature difference
between the baths. To understand why this is so, we put the limit of large
temperature gradient $`T_{L}>>T_{R}^{\prime}$ in the equations (II.4), (II.4)
and (II.4). Since we have fixed $T_{R}=1$, this means that
$`T_{L}>>1^{\prime}$ and $\beta_{R}\pm\beta_{L}\sim\beta_{R}$. Using this we
get the following asymptomatic expressions for currents:

$\displaystyle I_{L}$
$\displaystyle\sim\frac{J\sinh(2J\beta_{R})}{3\cosh(2J\beta_{R})+e^{2J\beta_{R}}},$
$\displaystyle I_{U}$
$\displaystyle\sim\frac{\frac{J}{3}\sinh(2J\beta_{R})}{3\cosh(2J\beta_{R})+e^{2J\beta_{R}}},$
$\displaystyle
I_{D}\sim\frac{\frac{-2J}{3}\sinh(2J\beta_{R})}{3\cosh(2J\beta_{R})+e^{2J\beta_{R}}}.$
(16)

We see that the above expressions are independent of $\beta_{L}(T_{L})$, hence
the current saturates for large values of $`T_{L}^{\prime}$. We also note from
Fig. 4$(c)$ that the total current flowing in the system does not only depend
upon the total number of spins but also on how they are distributed in
branches. In Fig. 4$(d)$ we see that the current decreases with an increase in
the number of spins which is consistent with our earlier proposition of the
spin-spin bonds behaving similarly to the resistors. From the above analysis,
it is clear that we will not get current magnification just by branch spin
number asymmetry for the Q2R dynamics and a dynamics with additional source of
asymmetry is required. We study one of such dynamics in the next section

###  II.5 Q2R dynamics with asymmetric upper and lower branch interaction
strength

Since, just the branch spin number asymmetry is not sufficient for generating
current magnification. We now employ different spin-spin interaction strength
in the upper and lower spin branch as a new source of asymmetry in our model.
This effectively means that we have two different thermal wires for upper and
lower branch. The Hamiltonian for this system is given as:

$\displaystyle
H_{S}=-J_{1}\sum_{i=1}^{Nu+1}\sigma_{i}\sigma_{i+1}-J_{2}\sum_{i=Nu+2}^{Ns}\sigma_{i}\sigma_{i+1}.$
(17)

Similar to the previous case, it is possible to analytically solve this model
for a $3$-spin system but the expressions are too complicated to be included
here. However, we can still infer the general characteristic of this model by
looking at the simulation results given in Fig. 5. We see in Fig. 5$(a)$ that
it is possible to get current magnification for asymmetric branch interaction
strengths in the Q2R model if the system parameters are in optimal range. This
deviates from our earlier observations and tells us that for different
interaction strengths in upper and lower branch we can no longer assume the
bond as being analogues to a resistor. In Fig. 5$(b)$ we see that
interchanging the number of spins in the branches still gives us current
circulation in the same direction but the magnitudes of currents change. We
also note (see Fig. 5$(c)$) that we need the temperature of one of the heat
baths to be an order below the energy scale of the bond interaction energy to
achieve significant current magnification. Interestingly, as seen Fig. 5$(d)$,
increasing the difference between interaction strength of the branches does
not necessarily increase the relative current magnification and a possibility
for optimisation of current magnification exists. We see that the direction of
current circulation changes depending on the relative magnitude of interaction
strength in upper and lower branch. This may result in a scenario where no
current flows through one spin branch and total current equals one of the
branch currents. This happens for $J_{2}\sim 2$ in Fig. 5$(d)$. We also see in
Fig. 5$(e)$ that the magnitudes of the currents decrease with an increase in
the number of spins and current magnification is possible even for same number
of spins in the branches. To check how the current magnification depends on
the spin-spin correlation defined as:

$\displaystyle\langle\sigma_{i}\sigma_{j}\rangle=\sum_{\\{\sigma\\}}P(\\{\sigma\\})\sigma_{i}\sigma_{j},$
(18)

where the sum runs over all the system configurations $\\{\sigma\\}$ and
$P(\sigma)$ is the probability of getting a particular configuration $\sigma$.
We plot the currents as well as the correlation between spin $2$ and $3$,
namely $\langle\sigma_{2}\sigma_{3}\rangle$ in the Fig. 5(f) and see that for
the region where we get current magnification, the spin–spin correlation is
high and when it decreases the magnification also decreases. We now discuss a
possible physical mechanism behind this.

###  II.6 Possible Physical Mechanism

Figure 6: Example of a process which result in the transfer of energy
$2(J_{1}-J_{2})$ from left bath to the right bath.

To make sense of the results observed above, we study the example of an energy
transfer process in a 4-spin system with $`N_{U}=1,N_{D}=1^{\prime}$ as shown
in Fig. 6. The configuration change steps corresponding to this process are,

$\displaystyle\\{\uparrow\uparrow\uparrow\downarrow\\}\rightarrow\\{\downarrow\uparrow\uparrow\downarrow\\}\rightarrow\\{\downarrow\downarrow\uparrow\downarrow\\}\rightarrow\\{\downarrow\downarrow\uparrow\uparrow\\}\rightarrow\\{\downarrow\downarrow\downarrow\uparrow\\}$
(19)

where as before the arrows indicate the states of spins numbered from left to
right. The spins numbered $1$ and $3$ interact with the left and right baths
respectively and follow Metropolis dynamics while spins $2$ and $4$ follow Q2R
dynamics. The bond interaction strength is $J_{1}$ and $J_{2}$ for upper and
lower spin branch respectively. The first step in the process involves the
flipping of first spin via the Metropolis algorithm and signifies a transfer
of energy $2(J_{1}-J_{2})$ from the left bath to the system (see Fig. 6(a)).
The next two steps in the process are shown in Fig. 6(b) and involve the
deterministic flipping of the bulk spins via the Q2R dynamics that result in
energy transfer of magnitude $2J_{1}$ in the upper branch and $-2J_{2}$ in the
lower branch successively. The final step involves the flipping of right nodal
spin which indicates an energy transfer of magnitude $2(J_{1}-J_{2})$ from the
upper system branch to the right bath and simultaneously a energy transfer of
$2J_{2}$ from the upper branch to the lower branch (see Fig. 6(d)). Looking at
the cumulative effect of all these steps, we get a process resulting in
transfer of energy $2(J_{1}-J_{2})$ from the left bath to the right bath
accompanied by the flow of $2J_{1}$ energy in the upper branch and $2J_{2}$
energy in the lower branch, hence resulting in current magnification and
clockwise current circulation inside the system. For the circulation of same
energy in the anti-clockwise direction, we will have to transfer the energy
$2(J_{1}-J_{2})$ from the right bath to the left bath which has less
probability because of the lower temperature of the right bath, hence this
process gives us a net current circulation.

Figure 7: (Color Online) Variation of Heat Currents ($I$) with (a) $T_{L}$,
with (b) $J_{2}$ for the system shown in Fig. 6 with $N_{D}=1,N_{U}=1$, For
all the cases, unless otherwise specified we use
$J_{1}=2,T_{L}=2,T_{R}=0.1,J_{2}=1.9.$

Figure 8: (Color Online) Variation of Heat Currents ($I$) for the CCA model
with (a) $T_{L}$ for $J_{1}=J_{2}=1$, with (b) $T_{L}$, with (c) $T_{R}$ and
(d) with $\mu$. For all the cases unless otherwise specified
$J_{1}=2,T_{L}=2,T_{R}=0.1,J_{2}=1,N_{D}=3,N_{U}=2$ and $\mu=4$.

The other way to transfer energy from left to the right bath includes a
process where the energy transfer is of magnitude $2(J_{1}+J_{2})$. This gives
us no current circulation, but since the energy is of higher magnitude this
process has less probability of occurrence. Thus overall the process with
current circulation dominates the energy transfer and hence current
magnification is happens on an average. But if the temperature of right bath
is increased, the probability of the reverse process also increases and
current circulation decreases as was seen in Fig. 5$(b)$ . Increasing the gap
between values of $J_{1}$ and $J_{2}$ increases the energy transfer magnitude
$2(J_{1}-J_{2})$ making the process less and less probabilistic, hence
reducing relative current magnification. Therefore large current magnification
is achieved if the interaction strengths are of similar order of magnitude
(see Fig. 5$(c)$ ). The above mechanism also explains why we will not get any
heat current for symmetric interaction strength in both the branches as the
net current due to this process becomes zero. No such heat transfer process is
possible when we have the two nodal spins directly connected to each other, so
there will be no current magnification in such a case as well. The
corresponding simulation results for the system shown in Fig. 6 are shown in
Fig. 7. We see qualitatively similar results for this case as in Fig. 5.
Though, the value of the circulating current is small pertaining to the
symmetry in the number of spin in upper and lower branch.

###  II.7 CCA Model

In the Q2R dynamics discussed above, the system Hamiltonian consisted of only
the interaction energy terms and as such the heat current in the bulk could
only be transferred through the passing of energy from one bond to the other.
We now work with a model which allows additional ways for the system to
possess energy. To achieve this we use the relatively complex CCA dynamics
Creutz (1986); Saito  _et al._ (1999) for studying the time evolution of bulk
spins. This dynamics allows for a change in interaction energy of the system
by using an additional source of energy at each site, a local “demon” Creutz
(1986). As stated earlier, for studying this dynamics we need to modify the
total Hamiltonian and include an extra energy term to each spin site. This
term corresponds to the presence of kinetic energy at each site and this
results in the following total Hamiltonian:

$\displaystyle
H_{S}=-J_{1}\sum_{i=1}^{Nu+1}\sigma_{i}\sigma_{i+1}-J_{2}\sum_{i=Nu+2}^{Ns}\sigma_{i}\sigma_{i+1}+\mu\sum_{i=1}^{i=Ns}\tilde{\sigma}_{i},$
(20)

where $\mu$ is the scaling parameter for the kinetic energy reserve at each
site and $\tilde{\sigma}_{i}$ represents the state of this new degree of
freedom with $\tilde{\sigma}_{i}\in\\{0,1,2,3\\}$. The presence of additional
energy reserves at each site leads to the possibility of more bulk spins flips
as some of the extra energy changes can be compensated. Since the magnitude of
energy supplied by the kinetic energy reserve is in multiples of $\mu$, the
bulk spin ‘$i$’ can only flip under this dynamics if the following condition
is satisfied:

$\displaystyle\mu\tilde{\sigma}_{i}-\Delta E_{i}\in\\{0,\mu,2\mu,3\mu\\},$
(21)

where, ‘$\Delta E_{i}=2J_{i}\sigma_{i}(\sigma_{i-1}+\sigma_{i+1})$’ and it can
take the values ‘$\Delta E_{U}\in\\{0,\pm 4J_{1}\\}$’ or ‘$\Delta
E_{D}\in\\{0,\pm 4J_{2}\\}$’ for upper and lower branch spins respectively.
The above equation tells us that the bulk spin ‘$i$’ can only flip if its
kinetic energy reserve can compensate for the corresponding energy cost. So,
the dynamics is still deterministic and total energy is conserved for the flip
of a bulk spin. A consequence of the above condition is that for CCA to allow
additional flips than Q2R, the values ‘$\mu$’ and ‘$J_{1},J_{2}$’ must be
related. To illustrate this we give a few particular examples. If we choose,
‘$J_{1}=2,J_{2}=1$’, and ‘$\mu=4$’ then CCA dynamics can be applied for both
the branches. Since, for $\mu=4$, $\mu\tilde{\sigma}\in\\{0,4,8,12\\}$,
$|\Delta E_{U}|\in\\{0,8\\}$ and $|\Delta E_{D}|\in\\{0,4\\}$. However, for
some other value say $\mu=2$, we have $\mu\tilde{\sigma}\in\\{0,2,4,6\\}$,
$|\Delta E_{U}|\in\\{0,8\\}$ and $|\Delta E_{D}|\in\\{0,4\\}$. Hence, for the
upper branch we essentially have the Q2R dynamics because the condition in Eq.
(21) is only satisfied for ‘$\Delta E_{U}=0$’, but for the lower branch CCA
dynamics still applies. We can see that irrespective of the particular values
of $J_{1}$, $J_{2}$ and $\mu$, the above condition is always satisfied for
‘$0$’ energy change and the Q2R dynamics always remains embedded in the CCA
dynamics. In particular the Q2R dynamics studied in the earlier sections is a
special case of the CCA dynamics for $\mu=0$. Since the interaction energy of
the spin bonds still remains the same, flipping of a spin signifies the same
transfer of energy as before but that energy could come either from the
neighboring bond or the kinetic energy reserve. The expressions for heat
currents remain the same and are still given by equations (6) and (9).
Simulating this model gives us the results shown in Fig. 8.

In Fig. 8(a), we see that it is not possible to get current magnification even
in the CCA model just by the branch spin number asymmetry and the total heat
current distribution still has Ohm’s law like characteristics for symmetric
upper and lower branch interaction strength. Similar to the earlier cases, we
see in Fig. 8(b) that we can get current magnification for this model if the
spin–spin interaction strength differs in upper and lower branch. In Fig.
8(c), we see that the current magnification is possible only when the
temperature of one of the baths is an order below the interaction energy of
the bonds. Finally, on studying the variation of current with‘ $\mu$’ in Fig.
8(d), we see that the steady state currents are independent of the values of
$\mu$, even though CCA dynamics allows additional flips when
$\mu\in\\{2,4,8\\}$. This means that we get similar steady state heat currents
for both CCA and Q2R models for similar values of system parameters. The extra
momentum term does not contribute critically to the steady state heat
currents. This in agreement with what was observed by Saito et al. in Saito
_et al._ (1999).

##  III Quantum Analysis

We now want to study the effects of branch spin number asymmetry on the heat
current transport in the quantum domain. To do this, we perform a
complementary study to the one performed by Xu et al. in Xu  _et al._ (2019).
It was shown by them that it is possible to get heat current magnification in
a quantum system with $4$ spins interacting via modified Heisenberg exchange
interaction if the onsite magnetic field is inhomogeneous and the total
magnetisation of the system is conserved. However, it is difficult to realise
this system experimentally as applying different magnetic fields at different
sites for such small systems is in general a difficult task. So in place of
the inhomogeneous magnetic field, we use unequal number of spins in branches
similar to the classical models. This will also help us to distinguish the
current magnification in quantum spin systems from the classical ones.

###  III.1 The model

We study a five spin quantum system containing unequal branch spin numbers as
shown in Fig 9. The spins are interacting with each other via a slightly
modified Heisenberg XXZ interaction Xu  _et al._ (2019) and the system
Hamiltonian is given as:

$\displaystyle\hat{H}_{S}$
$\displaystyle=J\sum_{i}(\hat{\sigma}^{i}_{x}\hat{\sigma}^{i+1}_{x}+\hat{\sigma}^{i}_{y}\hat{\sigma}^{i+1}_{y})+\Delta\hat{\sigma}^{i}_{z^{\prime}}\hat{\sigma}^{i+1}_{z^{\prime}},$
(22)

where, $J$ is the interaction strength between the $x$ and $y$ components of
neighboring spins, $\Delta$ is called the asymmetry parameter of the
Heisenberg XXZ interaction and signifies the interaction strength between $z$
components of the spins, $\hat{\sigma}_{\alpha}^{i}$ denotes the
$\alpha\in\\{x,y,z\\}$ component of the usual Pauli spin-$\frac{1}{2}$ matrix
for the $i^{th}$ spin. The unit operators are denoted by $\hat{I}$ and
$\hat{\sigma}^{i}_{z^{\prime}}=\frac{\hat{\sigma}^{i}_{z}+\hat{I}}{2}$. This
slight modification to the usual XXZ interaction is done so that the
transition to Fermionic system via the Jordan Wigner transformation Coleman
(2015) does not lead to a $\Delta$ dependent local chemical potential Xu  _et
al._ (2019). This system is interacting with two Bosonic baths with the
following Hamiltonian,

Figure 9:  Illustration of the Quantum Spin system consisting of 5 spins with
$N_{U}=2$ and $N_{D}=1$.

$\displaystyle\hat{H}_{B_{i}}$
$\displaystyle=\sum_{n}\omega_{n}\hat{a}_{n}^{i\dagger}\hat{a}^{i}_{n},$ (23)

where $\hat{a}_{n}^{i\dagger}(\hat{a}_{n}^{i})$ is the creation (destruction)
operator of the $n^{th}$ mode of the ‘$i^{th}$’ bath. To achieve the spin
branch configuration given in Fig 9, the spins numbered ‘1’ (‘4’) should
interact with the left (right) baths respectively. The system bath interaction
is chosen to preserve initial magnetisation of the system and has the form:

$\displaystyle\hat{H}_{SB_{L}}$
$\displaystyle=\hat{\sigma}^{1}_{z}\otimes\sqrt{\gamma}\sum_{n}c_{n}^{L}(\hat{a}^{L}_{n}+\hat{a}_{n}^{L\dagger}),$
$\displaystyle\hat{H}_{SB_{R}}$
$\displaystyle=\hat{\sigma}^{4}_{z}\otimes\sqrt{\gamma}\sum_{n}c_{n}^{R}(\hat{a}^{R}_{n}+\hat{a}_{n}^{R\dagger}),$
(24)

where, $\sqrt{\gamma}$ is the coupling strength between system and the baths
and its value is always taken as positive. Though the operator
$\hat{\sigma}_{z}$ can’t change the magnetisation of a spin, it can still lead
to a transfer of energy as the eigenstates of Hamiltonian and
$\hat{\sigma}_{z}$ are not the same.

###  III.2 Master Equation

To study the dynamics of the system, we again impose Markovian time evolution
and hence the density matrix of our system satisfies a master equation of the
type Breuer and Petruccione (2002),

$\frac{d}{dt}\hat{\rho}(t)=\mathcal{L}_{L}\hat{\rho}(t)+\mathcal{L}_{R}\hat{\rho}(t),$
(25)

where $`\hat{\rho}(t)^{\prime}$ is the density matrix of the system of
interest at time ‘t’ and $`\mathcal{L}_{L(R)}^{\prime}$ is the Liouvillian
operator of the left (right) bath. Specifically, for the analysis in this
manuscript, the Redfield master equation Cattaneo  _et al._ (2019); Breuer and
Petruccione (2002) is used. This is because the more simpler Lindblad master
equation doesn’t give the correct coherences in the steady state and hence we
obtain incorrect branch currents in that case Tupkary  _et al._ (2022); Xu
_et al._ (2019). The Redfield master equation is used for the weak coupling
regime between the system and the baths and is correct only upto order
$\gamma$ Cattaneo  _et al._ (2019); Breuer and Petruccione (2002). Its form
for our model in the Schr$\ddot{o}$dinger picture and in the energy basis is
given as (working in the units where $\hbar=1$),

$\displaystyle\dot{\hat{\rho}}_{nm}=$
$\displaystyle-i\Delta_{nm}\hat{\rho}_{nm}+\gamma\sum_{i,j}\big{(}\hat{S}^{L}_{ni}\hat{S}^{L}_{jm}[\tilde{W}^{L}_{jm}+W^{L}_{ni}]$
$\displaystyle-\delta_{mj}\sum_{l}\hat{S}^{L}_{nl}\hat{S}^{L}_{li}W^{L}_{li}(0)-\delta_{ni}\sum_{l}\hat{S}^{L}_{jl}\hat{S}^{L}_{lm}\tilde{W}^{L}_{jl}(0)\big{)}\rho_{ij}$
+Right Bath Terms (26)

The above equation can be divided into two parts, the first part not
multiplied by $\gamma$ is just the typical time evolution of an isolated
quantum system with the Schr$\ddot{o}$dinger equation and the part multiplied
by $\gamma$ corresponds to the time evolution of the system due to interaction
with the baths. Here $\hat{\rho}_{nm}$ is the matrix element of the system
density matrix and $`\Delta_{nm}^{\prime}$ is the energy difference between
the $n^{th}$ and $m^{th}$ energy levels or $\Delta_{nm}=E_{n}-E_{m}$.
$\hat{S}^{L}(\hat{S}^{R})$ is the operator of the system that interacts with
the left (right) bath and $\hat{S}^{L,R}=\hat{\sigma}^{1,4}_{z}$ for our case.
Finally, the $W$ matrices depend on the system energy spectrum and the heat
bath correlations. They are defined as:

$\displaystyle W^{L}_{ij}=\int_{0}^{\infty}dte^{-i\Delta_{ij}t}C^{L}(t),$ (27)
$\displaystyle\tilde{W}^{L}_{ij}=\int_{0}^{\infty}dte^{-i\Delta_{ij}t}C^{L}(-t),$
(28)

with $C^{L}(t)=Tr_{B}\text{ }(\hat{B}_{L}(t)\hat{B}_{L}\hat{\rho}_{B}(0))$,
$\hat{B}_{L}=\sum_{n}c_{n}^{L}(\hat{a}_{n}+\hat{a}_{n}^{\dagger})$ and
$\hat{\rho}_{B}(0)$ is the equilibrium density matrix corresponding to the
bath. Ignoring the Lamb shift terms Cattaneo  _et al._ (2019); Breuer and
Petruccione (2002), $\tilde{W}^{L(R)}_{j,i}=W^{L(R)}_{i,j}$ and we can write,

$\displaystyle
W^{L(R)}_{i,j}=\begin{cases}J(\Delta_{ij})N(\Delta_{ij},T_{L(R)})&\Delta_{ij}>0,\\\
J(|\Delta_{ij}|)\big{(}1+N(|\Delta_{ij}|,T_{L(R)})\big{)}&\Delta_{ij}<0,\end{cases}$
(29)

where, $N(\omega,T_{L(R)})$ is the Bose-Einstein distribution function,

$\displaystyle N(\omega,T_{L(R)})=\frac{1}{e^{\omega/T_{L(R)}}-1}.$ (30)

and we choose,

$J(\omega)=\omega e^{-{\omega}/{\Omega_{C}}},$ (31)

as the spectral density of the baths. We fix $\Omega_{C}=10$ without loss of
generality. Since, the above equation is only correct upto the order $\gamma$,
it means that $`\gamma^{\prime}$ should be smaller than $`1^{\prime}$ for any
meaningful discussion, so we take $\sqrt{\gamma}=0.1$ without loss of
generality. The functional form of the bath spectral density and numerical
values of $\Omega_{C}$ and $\gamma$ are similar to the ones taken in Xu  _et
al._ (2019).

Figure 10: (Color Online) Variation of Heat Currents ($I$) with the asymmetry
parameter $\Delta$ for (a) for setup 5(3$\uparrow$$2\downarrow$) (b) for setup
5($2\uparrow 3\downarrow$), (c) for a 6 spin system with setup
6(4$\uparrow$$2\downarrow$) with $N_{U}=3,N_{D}=1$ and (d) Variation of Heat
Currents ($I$), correlation and concurrence between spins 2 and 3 with $T_{L}$
for $\Delta=2.1$ . For all the cases unless otherwise specified we work with
the setup $5(3\uparrow 2\downarrow)$ and $T_{L}=2$, $T_{R}=1$, $J=1$. In (d)
we scale the correlation and concurrence by a factor $0.2$ for the visual
convenience.

###  III.3 Definitions

We now define the physical quantities required for studying the Quantum model.
Similar to the classical case, all the results are calculated for the steady
state and we denote our steady state density matrix as,

$\displaystyle\hat{\rho}_{S}=\lim_{t\to\infty}\hat{\rho}(t).$ (32)

The expression for the current flowing through a branch is derived by using
the Heisenberg equation of motion Breuer and Petruccione (2002) for the bond
energy operator between spins $l$ and $l+1$. The resulting operator for the
current flowing through the site ‘$l$’ of our system is given as Xu  _et al._
(2019):

$\displaystyle\hat{I}[l]=i[\hat{h}_{l-1,l},\hat{h}_{l,l+1}],$ (33)

with the bond energy operator inferred from the Eq. (22) as,

$\displaystyle\hat{h}_{l,l+1}$
$\displaystyle=J(\hat{\sigma}_{x}^{l}\sigma_{x}^{l+1}+\hat{\sigma}_{y}^{l}\hat{\sigma}_{y}^{l+1})+\Delta\hat{\sigma}_{z^{\prime}}^{l}\hat{\sigma}_{z^{\prime}}^{l+1}.$
(34)

Thus, the average steady state branch heat currents become:

$\displaystyle I_{U(D)}=Tr(\hat{I}[l]\cdot\hat{\rho}_{S}),$ (35)

where $I_{U(D)}$ is the heat current in the upper (lower) branch depending on
the position of the $l^{th}$ spin. Average current flowing out of a bath is
derived from the Liouville operator and is given as Karg ı  _et al._ (2019);
Upadhyay  _et al._ (2021):

$\langle I_{L}\rangle=Tr[\mathcal{L}_{L}\hat{\rho}_{S}\hat{H}_{S}].$ (36)

Magnetisation of the system is defined with respect to the $\hat{\sigma}_{z}$
operator of each spin with the excited eigenstate having magnetisation ‘+1’
and the ground eigenstate having magnetisation ‘-1’. Total magnetisation is
the difference between the spins in excited state and spins in ground state.
To avoid confusion we will use the notation
$N_{s}(N_{E}\uparrow,N_{G}\downarrow)$ where $N_{s}$ is the total number of
spins, $N_{E}(N_{G})$ is the number of spins in excited (ground) eigenstate.
The average magnetisation of a spin ‘i’ is defined as:

$\displaystyle M_{i}=Tr(\hat{\sigma}_{i}^{z}\cdot\hat{\rho}_{S}).$ (37)

The spin-spin correlation along $`z^{\prime}$ direction is given as:

$\displaystyle\langle\hat{\sigma}^{z}_{i}\hat{\sigma}^{z}_{j}\rangle=Tr(\hat{\sigma}_{i}^{z}\hat{\sigma}^{z}_{j}\cdot\hat{\rho}_{S}).$
(38)

We also use the concurrence Wootters (1998) as a measure of entanglement
between two spins, it is given by the following formula,

$\displaystyle C(\hat{\rho})=\text{max}[0,\lambda 1-\lambda 2-\lambda
3-\lambda 4],$ (39)

where $\lambda$’s are the eigenvalues in decreasing order of the matrix
$\hat{R}=\sqrt{\sqrt{\hat{\rho}}\tilde{\rho}\sqrt{\hat{\rho}}}$ and
$\hat{\rho}$ is the reduced density matrix for the two spin subsystem of
interest with
$\tilde{\rho}=(\hat{\sigma}_{y}\otimes\hat{\sigma}_{y})\hat{\rho}^{*}(\hat{\sigma}_{y}\otimes\hat{\sigma}_{y})$
and $\hat{\sigma}_{y}$ being the usual ‘Pauli Y’ matrix. We will discuss in
detail how the current magnification influences the correlations between the
spins and the concurrence later.

Finally, we introduce the concept of ergotropy Allahverdyan  _et al._ (2004),
which is defined as the maximum amount of work extractable from a quantum
state under unitary time evolution. The state which is attainable by the
unitary transformation and through which no work is extractable is called the
passive state (for example an equilibrium state). To find it, we first need to
write the spectral decomposition of the steady state density matrix as well as
the system Hamiltonian in the following manner:

$\displaystyle\hat{\rho}_{S}$
$\displaystyle=\sum_{i}r_{i}|\rho_{i}\rangle\langle\rho_{i}|,$
$\displaystyle\text{with }r_{i}\geq r_{i+1}$ $\displaystyle\hat{H}_{s}$
$\displaystyle=\sum_{i}E_{i}|E_{i}\rangle\langle E_{i}|,$
$\displaystyle\text{with }E_{i}\leq E_{i+1}$ (40)

where ‘$r_{i},E_{i}$’ and ‘$|\rho_{i}\rangle,|E_{i}\rangle$’ are the
eigenvalues and corresponding eigenvectors of the density matrix and the
system Hamiltonian respectively. Given the above decompositions, the passive
state $(\hat{\rho}_{passive})$ is defined as Xu  _et al._ (2019); Allahverdyan
_et al._ (2004),

$\displaystyle\hat{\rho}_{passive}=\sum_{i}r_{i}|E_{i}\rangle\langle E_{i}|.$
(41)

As can be seen from the above expressions, passive state is not uniquely
defined for degeneracies in energy spectrum of the system due to permutations
of degenerate energy states. However, the ergotropy $(\epsilon)$ can still be
uniquely defined and is given by,

$\displaystyle\epsilon$
$\displaystyle=Tr(\hat{H}_{S}\cdot\hat{\rho}_{S})-Tr(\hat{H}_{S}\cdot\hat{\rho}_{passive})$
$\displaystyle=\sum_{i,j}r_{i}E_{j}\big{(}|\langle\rho_{i}|E_{j}\rangle|^{2}-\delta_{i,j}\big{).}$
(42)

From the above discussion, we can infer that ergotropy can be used to quantify
how far a system is from the passive state.

Figure 11: (Color Online) (a) Variation of energy levels with $\Delta$, (b)
Variation of occupation probability ($P$) with $\Delta$, (c) Variation of
individual spin magnetisation ($M_{k}$) with $\Delta$ and (d) Variation of
ergotropy ($\epsilon$) with $\Delta$. For all the cases unless otherwise
specified, we work with the setup $5(3\uparrow 2\downarrow)$ and $T_{L}=2$,
$T_{R}=1$, $J=1$.

###  III.4 Information about numerical methodology

For solving the master equation given in Eq. (III.2), we again write it in the
matrix form, though we use the more effective numerical method for the
Redfield master equation defined in Pollard and Friesner (1994). Numerically
solving this equation by replacing the time derivative terms with
$`0^{\prime}$ and imposing the constraint of fixed initial magnetisation then
gives us the steady state density matrix Sandvik (2010). Once we have the
steady state density matrix, all the quantities of interest can be evaluated
from the definitions given in the subsection III.3. It is important to note
that the Redifield master equation is not guaranteed to be trace preserving
and totally positive Breuer and Petruccione (2002), so we also check if the
obtained solutions have negative probabilities. For the parameter range we
work on, no such regions were found. We now discuss in detail the results for
the Quantum model below.

###  III.5 Results

The results for the quantum model are shown in Fig. 10. We see that large
current magnification can be generated just by having unequal number of spins
in the branches if the total magnetisation is conserved and the asymmetry
parameter ‘$\Delta$’ is in the suitable range. In Fig 10(a), on plotting the
currents with the asymmetry parameter $\Delta$ for the $5$ spin setup
$5$($3\uparrow 2\downarrow$), we see that current magnification occurs for
positive values of asymmetric interaction parameter $`\Delta^{\prime}$. We
also see that the onset of current magnification is marked by a sudden dip in
the total current near $\Delta\sim 0.$ Interestingly on reversing the
magnetisation via a configuration change $5(3\uparrow 2\downarrow)\to
5(2\uparrow 3\downarrow)$, the total current remains the same and current
magnification also occurs for a certain parameter range but the branch
currents drastically change as seen in Fig. 10(b).

Figure 12: (Color Online) Variation of (a) energy levels and (b) currents with
a large range of asymmetry parameter $\Delta$. Note that only the magnitude of
currents are considered for the log plot in (b) but the background colors for
parallel, clockwise or anticlockwise circulating currents are kept. For all
the cases unless otherwise specified, we work with the setup $5(3\uparrow
2\downarrow)$ and $T_{L}=2$, $T_{R}=1$, $J=1$.

In Fig. 10(c), we observe that the dip in the total current is sharper for a
$6$ spin system and is accompanied by current magnification in a narrow range.
Finally on plotting the currents with the temperature of the left bath
‘$T_{L}$’ while keeping the temperature of right bath ‘$T_{R}$’ fixed in Fig.
10(d), we see that the magnitude of currents initially increase non-
monotonically with an increase in the value of ‘$T_{L}$’, but they plateau for
larger temperature differences. On the other hand, the corresponding current
magnification initially increases for a small region, reaches a maximum and
then starts decreasing before finally vanishing at some finite temperature
difference. We also look at the spin-spin correlation along ‘$z$’ direction
and concurrence between spins 2 and 3 in the Fig. 10(d) to check the
dependence of current magnification on correlations. We see that similar to
the current magnification, correlation and concurrence initially increase for
a small region and then start declining with an increase in temperature of the
left bath. Note that $\langle\hat{\sigma}^{z}_{2}\hat{\sigma}^{z}_{3}\rangle$
is negative because of the antiferromagnetic interaction between the spins. We
try to make sense of the above observations in the following discussion.

We discuss only about the $5$ spins setup $5$(3$\uparrow$2$\downarrow$). The
results in Fig. 10 above showed that the onset of current circulation is
accompanied by a dip in total current near $\Delta\sim 0$. To understand why
the total current drops, we study the properties of the occupied energy
levels. It turns out that the non-zero occupancy is seen in only $10$ energy
levels. These levels correspond to $10$ different configurations having the
same magnetisation. Out of these, only $6$ energy levels are non-degenerate,
numbering these energy levels in ascending order of their energies we plot
$E_{1},E_{2},E_{3},E_{4},E_{5},E_{6}$ with $\Delta$ in Fig. 11(a) and observe
that the energy levels $E_{4}$ and $E_{5}$ intersect each other at $\Delta=0$.
The other energy levels are sufficiently far from each other in the plotted
region. This means that the system looses one phonon transfer channel through
which heat current can pass at $\Delta=0$ and as a result the total heat
current decreases. To check what happens to the probability of occupancy
corresponding to these energy levels, we plot the variation of probability of
occupancy with $\Delta$ for the same range in Fig. 11(b). We see that $P_{2}$
suddenly drops and $P_{3}$ increases indicating population inversion between
energy levels $E_{2}$ and $E_{3}$ near $\Delta=0$. The probability of
occupancy for other energy levels behaves in the expected way with higher
energy levels having lesser probability of occupancy and vice versa. This can
only happen if the system is far from the passive state defined in (41) and
seems to suggest that the additional degeneracy in the two energy levels and
the constraint of fixed initial magnetisation is forcing the system to settle
at states that are unfavorable according to a passive state. We also plot the
variation of individual spin magnetisation with $\Delta$ in Fig. 11(c) and see
that though the total magnetisation remains fixed throughout this region,
individual spin magnetisation behave in an atypical manner near $\Delta\sim 0$
with local minima or maxima dipping for each spin. To support our observations
further, we evaluate the ergotropy according to the Eq. (III.3) but only
consider the energy levels with non zero probability of occupancy. Plotting
the variation of ergotropy with $\Delta$ in Fig. 11(d), we see that the
ergotropy suddenly starts rising around $\Delta\sim 0$ indicating the system
does move farther away from the passive state in this region. Also, on
increasing the temperature of the left bath, the ergotropy increases as the
system moves even further away from the passive state due to the increase in
the temperature gradient.

All the above results seem to suggest that whenever two energy levels
intersect, the physical observables of our system behave in an atypical
manner. To further probe this, we study the variation of energy levels and
currents for a larger range of the asymmetric parameter $\Delta$ in Fig. 12.
We see in Fig. 12(a) that only in the region where $\Delta$ is negative with
small magnitude $(-2.5<\Delta<0)$, all the energy levels have distinct values.
In Fig. 12(b), we see that this region also corresponds to the only parameter
range where we get parallel currents in the branches. Moving on either side of
this region leads to the energy levels coming closer to each other and the
transition from parallel to circulating currents. Specifically in Fig. 12(a)
for even more negative values of $\Delta(\sim-10)$, the energy levels
$E_{2},E_{3}$ and $E_{1}$ move towards each other. On the other side, as we
move to positive values of $\Delta(\sim 5)$, initially $E_{4}$, $E_{5}$ move
towards each other followed by $E_{6}$ coming closer to them for even larger
positive values of $\Delta(\sim 11)$. Correspondingly in 12(b), we see that
there is a transition from parallel current to clockwise circulating current
as we move away from $\Delta=0$ on either side, though the regions are not
symmetric with respect to $\Delta$. Also, on moving more towards the positive
side of $\Delta(\sim 6)$, there is a point of transition from clockwise
circulating current to anticlockwise circulating current. This point closely
corresponds to the intersection of energy levels $E_{2}$ and $E_{3}$.
Considering the above observations, it appears that whenever two energy levels
come towards each other, it can lead to a transition from parallel to
circulating currents in our system.

However, it is important to note that the current magnification discussed
above is only possible for temperature gradients of the same order as the
interaction strengths. In the high temperature regime current magnification is
absent as the inter-spin correlations are small and the currents just flow
parallel in the branches (see Fig. 10(d)). Though the above discussion does
help us to predict suitable conditions for current circulation namely ergodic
constraints and additional degeneracies. We note that analysing the phonon
transfer mechanism between different energy levels might help us to get better
insights into current magnification in quantum systems but this is beyond the
scope of the present manuscript Upadhyay  _et al._ (2021); Karg ı  _et al._
(2019).

##  IV Conclusion

To conclude, we study heat current magnification due to branch spin number
asymmetry in classical and quantum spin systems. We find that current
magnification is absent in the classical Q2R model for symmetric branch
interaction strengths. However, we encounter an interesting feature that the
heat current flowing through a spin branch is inversely proportional to the
number of spin-spin bonds in that branch. We employ different spin-spin
interaction strengths in the upper and lower branch and show that this
inequality is enough to generate as well as manipulate the current
magnification. This happens even if we have the same number of spins in both
the branches. We then provide a possible physical mechanism responsible for
current magnification in such systems. Similar features and heat current
values are found if we use the CCA model instead of the Q2R model. This shows
that the presence of additional momentum energy in the system does not
contribute appreciably to the steady state currents. This is in accordance
with the earlier studies with such dynamics. We then study a five spin Quantum
system with modified Heisenberg XXZ interaction and preserved magnetisation
using the Redfield master equation. We, with detailed numerical analysis, show
that it is possible to generate current magnification just by the branch spin
number asymmetry in this model for a suitable range of asymmetry parameter
$\Delta$. Our results indicate that the onset of current magnification is
accompanied by a sudden dip in total current which may be triggered due to the
intersection of two energy levels for certain values of $\Delta$. It is seen
that this value of $\Delta$ also corresponds to population inversion and
atypical trends for other physical observables e.g. magnetisation. Since the
two factors namely the ergodic constraint due to fixed magnetisation and
additional degeneracies in the system push the system away from the passive
state, we deduce that they may be the main causes for current magnification.
This deduction is also supported by the ergotropy evaluations for our system.
It is also seen that the total current is immune to the inversion of
individual magnetisation of the spins while the branch currents are not. In
both the classical as well as quantum case, a point exists where the current
in one of the branch is zero and the total current flows only through one of
the branches. This point corresponds to the transition between parallel and
circular currents. We also find that for both the classical and quantum
models, current magnification is only observed when temperature gradient and
intra-system interaction strength have similar order of energy. This also
points to the importance of correlation between spins required for current
magnification. In further studies, the relation of current magnification with
ergodicity can be studied as both the classical and quantum systems considered
here have non ergodic dynamics Stauffer (2000); Creutz (1986), the effect of
different spectral densities can also be explored in the quantum case.

##  V acknowledgement

RM gratefully acknowledges financial support from Science and Engineering
Research Board (SERB), India, under the Core Research Grant (Project No.
CRG/2020/000620). We also thank Özgür E. Müstecaplıoğlu and M. Tahir Naseem of
the Koç University, Istanbul for useful discussions and insights in this
project.

## References

  * Kaushik  _et al._ (2018) S. Kaushik, S. Kaushik, and R. Marathe, The European Physical Journal B 91, 87 (2018).
  * Upadhyay  _et al._ (2021) V. Upadhyay, M. T. Naseem, R. Marathe, and O. E. Müstecaplıoğlu, Phys. Rev. E 104, 054137 (2021).
  * Landi  _et al._ (2014) G. T. Landi, E. Novais, M. J. de Oliveira, and D. Karevski, Phys. Rev. E 90, 042142 (2014).
  * Jiang  _et al._ (2010) J. W. Jiang, J. S. Wang, and B. Li, Europhys. Lett. 89, 46005 (2010).
  * Ordonez-Miranda  _et al._ (2017) J. Ordonez-Miranda, Y. Ezzahri, and K. Joulain, Phys. Rev. E 95, 022128 (2017).
  * Segal and Nitzan (2005) D. Segal and A. Nitzan, Phys. Rev. Lett. 94, 034301 (2005).
  * Karg ı  _et al._ (2019) C. Karg ı, M. T. Naseem, T. c. v. Opatrný, O. E. Müstecaplıoğlu, and G. Kurizki, Phys. Rev. E 99, 042121 (2019).
  * Segal (2008) D. Segal, Phys. Rev. Lett. 100, 105901 (2008).
  * Zhang  _et al._ (2013) L. Zhang, J. Thingna, D. He, J.-S. Wang, and B. Li, EPL (Europhysics Letters) 103, 64002 (2013).
  * Terraneo  _et al._ (2002) M. Terraneo, M. Peyrard, and G. Casati, Phys. Rev. Lett. 88, 094302 (2002).
  * Saira  _et al._ (2007) O.-P. Saira, M. Meschke, F. Giazotto, A. M. Savin, M. Möttönen, and J. P. Pekola, Phys. Rev. Lett. 99, 027203 (2007).
  * Chung Lo  _et al._ (2008) W. Chung Lo, L. Wang, and B. Li, J. Phys. Soc. Jpn. 77, 054402 (2008).
  * Ruan and Wang (2020) Q. Ruan and L. Wang, Phys. Rev. Research 2, 023087 (2020).
  * Naseem  _et al._ (2020) M. T. Naseem, A. Misra, O. E. Müstecaplioğlu, and G. Kurizki, Phys. Rev. Research 2, 033285 (2020).
  * Jiang  _et al._ (2015) J.-H. Jiang, M. Kulkarni, D. Segal, and Y. Imry, Phys. Rev. B 92, 045309 (2015).
  * Zhang  _et al._ (2018) Y. Zhang, Z. Yang, X. Zhang, B. Lin, G. Lin, and J. Chen, EPL (Europhysics Letters) 122, 17002 (2018).
  * Sánchez  _et al._ (2017) R. Sánchez, H. Thierschmann, and L. W. Molenkamp, Phys. Rev. B 95, 241401 (2017).
  * Messina  _et al._ (2012) R. Messina, M. Antezza, and P. Ben-Abdallah, Phys. Rev. Lett. 109, 244302 (2012).
  * Hwang  _et al._ (2009) J. Hwang, M. Pototschnig, R. Lettow, G. Zumofen, A. Renn, S. Götzinger, and V. Sandoghdar, Nature 460 (2009).
  * Dhar (2008) A. Dhar, Adv. Phys. 57, 457 (2008).
  * Poulsen  _et al._ (2022a) K. Poulsen, A. C. Santos, L. B. Kristensen, and N. T. Zinner, Phys. Rev. A 105, 052605 (2022a).
  * Poulsen  _et al._ (2022b) K. Poulsen, A. C. Santos, and N. T. Zinner, Phys. Rev. Lett. 128, 240401 (2022b).
  * Chioquetta  _et al._ (2021) A. Chioquetta, E. Pereira, G. T. Landi, and R. C. Drumond, Phys. Rev. E 103, 032108 (2021).
  * Li  _et al._ (2004) B. Li, L. Wang, and G. Casati, Phys. Rev. Lett. 93, 184301 (2004).
  * Cho and McKenzie (2005a) S. Y. Cho and R. H. McKenzie, Phys. Rev. B 71, 045317 (2005a).
  * Lepri  _et al._ (2003) S. Lepri, R. Livi, and A. Politi, Physics Reports 377, 1 (2003).
  * Li  _et al._ (2012) N. Li, J. Ren, L. Wang, G. Zhang, P. Hänggi, and B. Li, Rev. Mod. Phys. 84, 1045 (2012).
  * Pereira (2017) E. Pereira, Phys. Rev. E 96, 012114 (2017).
  * Kosloff (2013) R. Kosloff, Entropy 15, 2100 (2013).
  * Kirchberg and Nitzan (2022) H. Kirchberg and A. Nitzan, The Journal of Chemical Physics 156, 094306 (2022).
  * Xu  _et al._ (2019) X. Xu, K. Choo, V. Balachandran, and D. Poletti, Entropy 21, 228 (2019).
  * Dugar and Chien (2022) P. Dugar and C.-C. Chien, Phys. Rev. E 105, 064111 (2022).
  * Roy (2007) D. Roy, Journal of Physics: Condensed Matter 20, 025206 (2007).
  * Rai  _et al._ (2010) D. Rai, O. Hod, and A. Nitzan, The Journal of Physical Chemistry C 114, 20583 (2010).
  * Marathe  _et al._ (2010) R. Marathe, A. Dhar, and A. M. Jayannavar, Phys. Rev. E 82, 031117 (2010).
  * Bandopadhyay  _et al._ (2004) S. Bandopadhyay, P. S. Deo, and A. M. Jayannavar, Phys. Rev. B 70, 075315 (2004).
  * Jayannavar  _et al._ (1995) A. Jayannavar, P. Singha Deo, and T. Pareek, Physica B: Condensed Matter 212, 261 (1995), proceedings of the Workshop on Novel Physics in Low-Dimensional Electron Systems.
  * Cho and McKenzie (2005b) S. Y. Cho and R. H. McKenzie, Phys. Rev. B 71, 045317 (2005b).
  * Xiao  _et al._ (2013) S. Xiao, X. Yang, L. Cai, and T. C. Liu, The European Physical Journal B 86, 77 (2013).
  * Nakanishi and Tsukada (2001) S. Nakanishi and M. Tsukada, Phys. Rev. Lett. 87, 126801 (2001).
  * Eckmann and Zabey (2004) J. P. Eckmann and E. Zabey, Journal of Statistical Physics 114, 515 (2004).
  * Jayannavar and Singha Deo (1995) A. M. Jayannavar and P. Singha Deo, Phys. Rev. B 51, 10175 (1995).
  * Vichniac (1984) G. Y. Vichniac, Physica D: Nonlinear Phenomena 10, 96 (1984).
  * Saito  _et al._ (1999) K. Saito, S. Takesue, and S. Miyashita, Phys. Rev. E 59, 2783 (1999).
  * Creutz (1986) M. Creutz, Annals of Physics 167, 62 (1986).
  * Schulte  _et al._ (1987) M. Schulte, W. Stiefelhagen, and E. S. Demme, Journal of Physics A: Mathematical and General 20, L1023 (1987).
  * Stauffer (2000) D. Stauffer, Computer Physics Communications 127, 113 (2000).
  * Coldea  _et al._ (2010) R. Coldea, D. A. Tennant, E. M. Wheeler, E. Wawrzynska, D. Prabhakaran, M. Telling, K. Habicht, P. Smeibidl, and K. Kiefer, Science 327, 177 (2010).
  * Wolf (2000) W. P. Wolf, Brazilian Journal of Physics 30 (2000).
  * Micadei  _et al._ (2019) K. Micadei, J. P. S. Peterson, A. M. Souza, R. S. Sarthour, I. S. Oliveira, G. T. Landi, T. B. Batalhão, R. M. Serra, and E. Lutz, Nat. Commun. 10, 2456 (2019).
  * Zhang  _et al._ (2009) J. Zhang, F. M. Cucchietti, C. M. Chandrashekar, M. Laforest, C. A. Ryan, M. Ditty, A. Hubbard, J. K. Gamble, and R. Laflamme, Phys. Rev. A 79, 012305 (2009).
  * Seo  _et al._ (2013) M. Seo, H. K. Choi, S.-Y. Lee, N. Kim, Y. Chung, H.-S. Sim, V. Umansky, and D. Mahalu, Phys. Rev. Lett. 110, 046803 (2013).
  * Rogge and Haug (2008) M. C. Rogge and R. J. Haug, Phys. Rev. B 77, 193306 (2008).
  * Bermudez  _et al._ (2013) A. Bermudez, M. Bruderer, and M. B. Plenio, Phys. Rev. Lett. 111, 040601 (2013).
  * Porras and Cirac (2004) D. Porras and J. I. Cirac, Phys. Rev. Lett. 92, 207901 (2004).
  * Pollard and Friesner (1994) W. T. Pollard and R. A. Friesner, The Journal of Chemical Physics 100, 5054 (1994).
  * Wolfram (1983) S. Wolfram, Rev. Mod. Phys. 55, 601 (1983).
  * Sandvik (2010) A. W. Sandvik, AIP Conference Proceedings 1297, 135 (2010).
  * Coleman (2015) P. Coleman, “Simple examples of second quantization,” in  _Introduction to Many-Body Physics_ (Cambridge University Press, 2015) p. 71–94.
  * Breuer and Petruccione (2002) H. P. Breuer and F. Petruccione,  _The theory of open quantum systems_ (Oxford university press, 2002).
  * Cattaneo  _et al._ (2019) M. Cattaneo, G. L. Giorgi, S. Maniscalco, and R. Zambrini, New Journal of Physics 21, 113045 (2019).
  * Tupkary  _et al._ (2022) D. Tupkary, A. Dhar, M. Kulkarni, and A. Purkayastha, Phys. Rev. A 105, 032208 (2022).
  * Wootters (1998) W. K. Wootters, Phys. Rev. Lett. 80, 2245 (1998).
  * Allahverdyan  _et al._ (2004) A. E. Allahverdyan, R. Balian, and T. M. Nieuwenhuizen, Europhysics Letters 67, 565 (2004).