# Entanglement and thermokinetic uncertainty relations in coherent mesoscopic
transport

Kacper Prech<EMAIL_ADDRESS>Department of Physics, University of
Basel, Klingelbergstrasse 82, 4056 Basel, Switzerland  Philip Johansson  Elias
Nyholm  Physics Department, Lund University, Box 118, 22100 Lund, Sweden
Gabriel T. Landi  Instituto de Física da Universidade de S̃ao Paulo, 05314-970
S̃ao Paulo, Brazil  Claudio Verdozzi  Physics Department and ETSF, Lund
University, PO Box 118, 221 00 Lund, Sweden  Peter Samuelsson  Physics
Department and NanoLund, Lund University, Box 118, 22100 Lund, Sweden  Patrick
P. Potts  Department of Physics, University of Basel, Klingelbergstrasse 82,
4056 Basel, Switzerland

###### Abstract

A deeper understanding of the differences between quantum and classical
dynamics promises great potential for emerging technologies. Nevertheless,
some aspects remain poorly understood, particularly concerning the role of
quantum coherence in open quantum systems. On the one hand, coherence leads to
entanglement and even nonlocality. On the other, it may lead to a suppression
of fluctuations, causing violations of thermokinetic uncertainty relations
(TUR and KUR) that are valid for classical processes. These represent two
different manifestations of coherence, one depending only on the state of the
system (static) and one depending on two-time correlation functions
(dynamical). Here we employ these manifestations of coherence to determine
when mesoscopic quantum transport can be captured by a classical model based
on stochastic jumps, and when such a model breaks down, implying nonclassical
behavior. To this end, we focus on a minimal model of a double quantum dot
coupled to two thermal reservoirs. In this system, quantum tunneling induces
Rabi oscillations and results in both entanglement and nonlocality, as well as
TUR and KUR violations. These effects, which describe the breakdown of a
classical description, are accompanied by a peak in coherence. Our results
provide guiding principles for the design of out-of-equilibrium devices that
exhibit nonclassical behavior.

##  I Introduction

The inability of classical physics to account for experimental observations on
the atomic scale resulted in the emergence of quantum mechanics in the last
century [1]. Quantum physics exhibits behavior that may be strikingly
different from that of classical physics. One fundamentally distinguishing
feature is quantum coherence, i.e., the ability to maintain a superposition of
different quantum states. Coherence can be understood as a resource in quantum
information processing [2], to realize tasks such as quantum computation or
quantum metrology. Coherence also plays an important role in the
thermodynamics of quantum systems and influences fluctuations of work and
entropy production [3, 4, 5].

Despite the enormous interest in quantum advantages for technological
applications, the difference between quantum and classical behavior has not
been fully understood yet. One particular instance is that of quantum
transport and thermodynamics in open quantum systems. Such systems may usually
be described by multiple coupled sites, in contact with thermal reservoirs at
different temperatures and/or chemical potentials (Fig. 1). After a long time
elapses, this system will reach a nonequilibrium steady state (NESS),
characterized by currents flowing between the reservoirs. This setup forms the
basis for a variety of devices and applications, including thermoelectrics [6,
7], autonomous engines [8], and nanoscale sensors [9].

Figure 1: Manifestations of coherent tunneling in mesoscopic transport. The
buildup of quantum coherence, caused by Rabi oscillations, may result in
entanglement and even nonlocality. In addition, coherent tunneling may also
result in a suppression of fluctuations, leading to violations of
thermokinetic uncertainty relations.

Transport through such a system is enabled by quantum tunneling between sites.
However, under certain circumstances, a classical model based on stochastic
jumps may fully captures the nonequilibrium behavior. To tap into the
potential provided by quantum coherence, it is thus of crucial importance to
determine when a classical model is adequate to describe an out-of-equilibrium
quantum system and when it breaks down. To determine the breakdown of a
classical description, we consider two manifestations of coherence:
entanglement in the NESS [10, 11, 12, 13, 14] and violations of thermokinetic
uncertainty relations [15, 16, 17, 18, 19, 20, 21, 22].

Entanglement [23] is an exceptionally strong form of correlation that may
result in nonlocal behavior [24], which may have an important technological
impact as it can be exploited for, e.g., quantum cryptography [25, 26] or to
provide advantages in sensing applications [27, 28, 29]. Entanglement, though,
depends solely on the quantum state of the system at a given time; i.e., it is
a static quantity.

The second manifestation is provided by violations of the so-called
thermodynamic uncertainty relations (TURs) [30, 31, 32, 33] and kinetic
uncertainty relations (KURs) [33], as well as their recently proposed
unification, the TKUR [34]. These are bounds on the signal-to-noise ratio of
currents in classical systems. In quantum-coherent systems, these bounds no
longer need to be respected, so that TUR or KUR violations can be interpreted
as manifestations of quantum coherence. In contrast to entanglement, the TUR
and KUR bounds depend on two-time correlations. These two types of
manifestations of coherence are therefore complementary: while entanglement
can be regarded as a static manifestation of coherence, TUR and KUR violations
are dynamical in nature.

Here we focus on a minimal model, containing just two sites and two reservoirs
(see Fig. 1). In this setup, both NESS entanglement [10, 11, 12] as well as
TUR violations [16, 17] are present. Investigating these manifestations of
coherence together, and comparing to a stochastic jump model, allows us
determine when quantum coherence results in a behavior that cannot be captured
by a classical model. We find that the buildup of coherence due to quantum
tunneling can become large enough to generate not only entanglement but even
nonlocality, in contrast to previous studies [11]. Furthermore, the Rabi
oscillations induced by quantum tunneling result in a reduction of
fluctuations that may overcome both the TUR as well as the KUR. To the best of
our knowledge, such KUR (and TKUR) violations have not been reported in
quantum systems before. Interestingly, both manifestations of coherence occur
when the tunnel coupling between the sites is comparable to the system-bath
coupling, where the coherence in the NESS shows a local maximum. In addition,
entanglement is found for strong tunnel couplings, where coherence exhibits
its global maximum. In this regime, however, the dynamics is well captured by
a classical jump model that involves the entangled eigenstates of the two-site
system. Thermokinetic uncertainty relations may therefore not be violated in
this strong-coupling regime.

For concreteness, we consider a double quantum dot, where electrons can hop
between two levels. The basic model parameters are outlined in Fig. 2. In
order to contrast classical models and quantum-coherent models, we employ
different quantum master equations [35, 36] (local vs global) and investigate
the effect of interactions. For a noninteracting system, we benchmark local
and global quantum master equations with nonequilibrium Green’s functions
(NEGFs). We note that we benchmark not only average quantities, depending on
the steady state, but also current fluctuations which depend on two-time
correlation functions. This further clarifies the regimes of validity of the
local and global master equation which have been debated recently [35, 37, 38,
39]. In the presence of interactions, we employ a thermodynamically consistent
semilocal master equation [36, 40], as well as NEGFs. While our results are
valid for the two-site model sketched in Figs. 1 and 2, the insights gained in
the breakdown of a classical description carry over to more complicated
transport scenarios, serving as guiding principles in designing out-of-
equilibrium devices that can behave nonclassically. We illustrate this by
investigating a chain of three quantum dots.

This paper is organized as follows. The system we consider and the models we
employ for its description are discussed in Sec. II. The manifestations of
coherence are introduced in Sec. III. In Sec. IV, we present the results for
noninteracting electrons and in Sec. V we discuss the role of interactions.
Section VI contains a discussion on more complicated transport scenarios,
including a chain of three quantum dots. Conclusions are provided in Sec. VII.

##  II The system and models

Figure 2: Sketch of the system. Two spin-polarized quantum dots (QDs)
($\ell=L,R$) with on-site energies $\epsilon$ are tunnel-coupled with strength
$g$ and exhibit an interdot Coulomb interaction $U$. Each QD is weakly coupled
to a fermionic reservoir with temperature $T_{\ell}$ and chemical potential
$\mu_{\ell}$. The coupling between system and reservoir is denoted by
$\gamma_{\ell}$.

###  II.1 Double quantum dot

We consider two quantum dots (QDs) connected in series, labeled $\ell=L,R$
[41, 42, 43]. Each QD is weakly coupled to a fermionic reservoir with
temperature $T_{\ell}$ and chemical potential $\mu_{\ell}$, respectively (Fig.
2). As shown previously, this system exhibits both manifestations of coherence
we are interested in: First, by driving an electron current through the system
it is possible to generate an entangled state [10]. Second, fluctuations of
the same current may be suppressed below the classical TUR bound [16]. We
consider spin-polarized electrons, such that at most one electron may reside
on each QD. The Hamiltonian of the system reads

$\hat{H}=\sum_{\ell=L,R}\epsilon\hat{c}_{\ell}^{\dagger}\hat{c}_{\ell}+g\left(\hat{c}_{L}^{\dagger}\hat{c}_{R}+\hat{c}_{R}^{\dagger}\hat{c}_{L}\right)+U\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R},$
(1)

where $\epsilon$ is the onsite energy of the QDs, $g$ the interdot tunnel
coupling, $U$ the Coulomb repulsion between electrons, and $\hat{c}_{\ell}$
($\hat{c}_{\ell}^{\dagger}$) the annihilation (creation) operator of the
electron on QD $\ell$, respectively. The effect of unequal onsite energies is
considered in Appendix A.

The effect of the reservoirs is described by multiple methods, including local
and global Lindblad master equations as well as NEGFs. We now introduce these
methods for noninteracting electrons ($U=0$); the interacting case is treated
in Sec. V. The strength of the coupling between QD $\ell$ and the
corresponding reservoir (assumed to be energy-independent) will be denoted by
$\gamma_{\ell}$ for all methods. Throughout, we consider a voltage bias that
is applied symmetrically such that $\mu_{L}=\epsilon+eV/2$ and
$\mu_{R}=\epsilon-eV/2$.

###  II.2 Lindblad master equations

The main method we use to investigate the dynamics of the double QD is
provided by Lindblad master equations. These rely on Born-Markov
approximations that for our system are justified when (we set $k_{B}=\hbar=1$
throughout)

$\gamma_{\ell}\ll\max\\{T_{\ell^{{}^{\prime}}},|\epsilon\pm
g-\mu_{\ell^{{}^{\prime}}}|\\},\qquad\text{(Born-Markov)}$ (2)

for any $\ell$ and $\ell^{{}^{\prime}}$ and for both signs. This condition
ensures that the relevant bath properties are flat around the Bohr frequencies
$\epsilon\pm g$. To ensure thermodynamic consistency, we follow a recently
derived framework [36, 40]. Both the local as well as the global master
equation that we employ can be written in the form

$\frac{d\hat{\rho}}{dt}=-i\left[\hat{H},\hat{\rho}\right]+\mathcal{L}_{L}\hat{\rho}+\mathcal{L}_{R}\hat{\rho}:=\mathcal{L}\hat{\rho},$
(3)

where $\mathcal{L}_{\ell}$ is the Lindblad superoperator that accounts for the
dissipation to the bath $\ell$. For the local master equation, the Lindblad
superoperators are given by [36]

$\mathcal{L}_{\ell}=n_{\ell}\gamma_{\ell}\mathcal{D}[\hat{c}^{\dagger}_{\ell}]+(1-n_{\ell})\gamma_{\ell}\mathcal{D}[\hat{c}_{\ell}],$
(4)

where

$\mathcal{D}[\hat{O}]\hat{\rho}=\hat{O}\hat{\rho}\hat{O}^{\dagger}-\frac{1}{2}\\{\hat{O}^{\dagger}\hat{O},\hat{\rho}\\},$
(5)

and $n_{\ell}:=n_{\ell}(\epsilon)$, with the Fermi-Dirac distribution

$n_{\ell}(\omega)=\frac{1}{\exp{\left(\frac{\omega-\mu_{\ell}}{T_{\ell}}\right)}+1}.$
(6)

When the Born-Markov approximations are justified, the local master equation
(ME) may be employed when [36]

$g\ll\max\\{T_{\ell},|\epsilon-\mu_{\ell}|\\}\qquad\text{(local ME)}.$ (7)

In analogy to Eq. (2), this condition ensures that the relevant bath
properties are flat across both Bohr frequencies. In this work, we use the
local master equation extensively, as its regime of validity coincides with
the regime where both manifestations of coherence appear.

In addition, we make use of the global master equation which relies on the
secular approximation [44]. The Lindblad superoperators in the global approach
read [36]

$\mathcal{L}_{\ell}=\frac{1}{2}\sum_{s=\pm}\left\\{n_{\ell}^{s}\gamma_{\ell}\mathcal{D}[\hat{c}_{s}^{\dagger}]+(1-n_{\ell}^{s})\gamma_{\ell}\mathcal{D}[\hat{c}_{s}]\right\\},$
(8)

where

$\epsilon_{\pm}=\epsilon\pm g,\hskip
42.67912pt\hat{c}_{\pm}=\frac{1}{\sqrt{2}}\left(\hat{c}_{R}\pm\hat{c}_{L}\right),$
(9)

and $n_{\ell}^{s}:=n_{\ell}(\epsilon_{s})$. The annihilation operators
$\hat{c}_{\pm}$ correspond to the eigenmodes of Hamiltonian (1). The secular
approximation is justified when [44]

$g\gg\max\\{\gamma_{L},\gamma_{R}\\}\qquad\text{(global ME).}$ (10)

In both Lindblad master equations, the steady-state average current flowing
across the system can be computed with the expression

$\langle I\rangle=\text{Tr}\left[\hat{N}\mathcal{L}_{L}\hat{\rho}\right],$
(11)

where

$\hat{N}=\hat{c}_{L}^{\dagger}\hat{c}_{L}+\hat{c}_{R}^{\dagger}\hat{c}_{R}$
(12)

is the operator of the number of electrons in the double quantum dot. To
determine the current fluctuations, we resort to the method of full counting
statistics discussed in Appendix B.

We note that in both master equations, we neglect a Lamb shift that results in
a shift of the on-site energies.

###  II.3 Nonequilibrium Green’s functions

For noninteracting electrons, all relevant quantities can be solved exactly
using the method of NEGFs. The transmission function $\mathcal{T}(\omega)$,
which specifies the transition rate of electrons with energy $\omega$, is a
fundamental quantity to investigate transport across the system. For our
system it is given by [45, 17]

$\mathcal{T}(\omega)=\frac{\gamma_{L}\gamma_{R}g^{2}}{|\left(\omega-\epsilon+i\frac{\gamma_{L}}{2}\right)\left(\omega-\epsilon+i\frac{\gamma_{R}}{2}\right)-g^{2}|^{2}},$
(13)

where we assumed the wideband approximation [46], i.e., energy independence of
the coupling strength $\gamma_{\ell}$ and omitting any Lamb shift. The
transmission function allows us to compute the average of the current and its
fluctuations, given by [47, 48, 49]

$\langle
I\rangle=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi}\mathcal{T}(\omega)\left[n_{L}(\omega)-n_{R}(\omega)\right],$
(14)

and

$\begin{split}\langle\\!\langle
I^{2}\rangle\\!\rangle&=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi}\mathcal{T}(\omega)\\{n_{L}(\omega)+n_{R}(\omega)-2n_{L}(\omega)n_{R}(\omega)\\\
&-\mathcal{T}(\omega)[n_{L}(\omega)-n_{R}(\omega)]^{2}\\}.\end{split}$ (15)

The elements of the density matrix of the system can be computed from the
lesser Green’s functions together with Wick’s theorem. More information on the
NEGF formalism can be found in Appendix A.3.

###  II.4 Classical model

The last method that we employ is a classical, stochastic model for the double
QD system [50, 15]. We consider the time evolution of the vector of
probabilities $\vec{p}=[p_{0},p_{L},p_{R},p_{D}]^{T}$, where entries denote,
respectively, the probability that the system is empty, only the left QD is
occupied, only the right QD is occupied, and both QDs are occupied. In the
classical model, the time evolution of $\vec{p}$ is governed by the rate
equation

$\frac{d}{dt}\vec{p}=W\vec{p},$ (16)

where the elements $W_{i\neq j}$ of the $4\times 4$ matrix $W$ are the
transition rates from state $j$ to state $i$, and $W_{ii}=-\sum_{j\neq
i}W_{ji}$. The transition rates between states with a different total number
of electrons are taken from the local master equation, i.e.,
$W_{L0}=W_{DR}=n_{L}\gamma_{L}$, $W_{0L}=W_{RD}=(1-n_{L})\gamma_{L}$, and
similarly for $L\leftrightarrow R$. The transition rate describing hopping
between the dots can be obtained from perturbation theory (cf. Appendix A.4)
and reads

$W_{LR}=W_{RL}=\frac{4g^{2}}{\gamma_{L}+\gamma_{R}}.$ (17)

The average current can be evaluated as

$\langle I\rangle=W_{L0}p_{0}-W_{0L}p_{L}+W_{DR}p_{R}-W_{RD}p_{D}.$ (18)

To compute the fluctuations of the current we resort to full counting
statistics (cf. Appendix B), as in the case of the Lindblad master equations.

##  III Manifestations of coherence

###  III.1 Coherence

We study coherence in the quantum state $\hat{\rho}$ of the system with
respect to the basis of the particle occupation number of each quantum dot. A
local Fock basis
$\\{|n_{L},n_{R}\rangle:=(\hat{c}_{L}^{\dagger})^{n_{L}}(\hat{c}_{R}^{\dagger})^{n_{R}}|0,0\rangle\\}$
is constructed such that the first (second) number denotes the occupation of
the left (right) QD. Then, the density matrix of the system can be written in
the basis $\\{|0,0\rangle,|1,0\rangle,|0,1\rangle,|1,1\rangle\\}$ as a matrix

$\hat{\rho}=\begin{pmatrix}p_{0}&0&0&0\\\ 0&p_{L}&\alpha&0\\\
0&\alpha^{*}&p_{R}&0\\\ 0&0&0&p_{D}\end{pmatrix}.$ (19)

It has a single off-diagonal element $\alpha$, corresponding to a single
electron being in a superposition between the left and the right dots. Any
other off-diagonal element is strictly zero due to the charge superselection
rule [51, 52]. Throughout, we will quantify coherence with the absolute value
$|\alpha|$, which is equivalent to one half of the $l_{1}$-norm of coherence
[53].

###  III.2 Entanglement

The system sketched in Fig. 2 and its modifications have attracted a
significant amount of attention because they may feature entanglement in the
steady state [10, 11, 12, 13, 14]. This entanglement may arise from two
different mechanisms. For tunnel couplings comparable to the system-bath
coupling, the entanglement is generated by a current passing through the
double dot [12]. This scheme to produce nonseparable states can be extended to
multiqubit machines [54]. It can also be improved by local filtering [55], and
it can be mediated by a cavity mode [56].

For tunnel couplings much larger than the temperature, only the lowest energy
eigenstate is occupied. This state is an entangled state such that the thermal
state becomes entangled. In this paper, we are mainly interested in the
generation of entanglement by passing a current through the system because
this is the regime where we also observe other manifestations of coherence.

We emphasize that the introduced scheme generates mode entanglement between
two QDs, where the number of particles plays the role of the degree of freedom
[57, 58, 52, 59]. This is qualitatively different from entanglement between
two particles, such as spin entanglement in electronic systems. We note that
the entanglement of a single electron [60] has been debated in the literature
due to the charge superselection rule. However, following [58], having two
copies of the system allows to overcome this issue and makes the entanglement
of a single electron useful.

To quantify the amount of entanglement we use the entanglement measure of
concurrence [61, 62]. For the present system, the concurrence is given by

$\mathcal{C}=\text{max}\left\\{2|\alpha|-2\sqrt{p_{0}p_{D}},0\right\\}.$ (20)

From this expression, we can already anticipate that Coulomb interactions have
a favorable effect on the entanglement, because they reduce $p_{D}$. This has
a twofold effect as it allows for a larger occupation in the single-electron
subspace, which may result in a larger $\alpha$.

###  III.3 Thermodynamic and kinematic uncertainty relations

The TUR and KUR are two different inequalities that bound the signal-to-noise
ratio of a current in classical systems. The TUR has been developed in the
field of classical stochastic thermodynamics [63, 64, 65, 66, 67], which is a
framework to describe fluctuations of thermodynamic quantities, such as work
or heat, in nanoscale systems. It provides the bound in terms of the average
entropy production rate, which quantifies the irreversibility of a process.
The TUR can be applied to obtain a trade-off between power production,
constancy, and efficiency in stochastic heat engines [68] and to infer the
efficiency of molecular motors [69].

The KUR was derived using the fluctuation-response inequality for out-of-
equilibrium stochastic dynamics [70] and finds applications outside of the
thermodynamic settings. In the case of the KUR, the signal-to-noise ratio is
bounded by the average dynamical activity, which is a measure of the total
rate of transitions between the states of the system. The question whether the
KUR can be violated in a quantum system has, to the best of our knowledge, not
been explored. We note that in contrast to entropy production, which can be
derived using any model, subtleties arise when computing the dynamical
activity that features in the KUR with different models. These are discussed
below.

While the TUR provides a tight bound when dissipation is small, the KUR is
more relevant for systems that are far from equilibrium. Either of these two
different inequalities may thus provide a stronger bound on the fluctuations
depending on the situation. As discussed below, these inequalities have also
recently been combined into the TKUR [34].

The TUR quantifies a trade-off between an average of a current, $\langle
I\rangle$, its fluctuations $\langle\\!\langle I^{2}\rangle\\!\rangle$, and
the total entropy production rate $\langle\sigma\rangle$. It was initially
proven in the long-time limit for nonequilibrium systems that follow a time-
homogeneous Markovian dynamics and obey local detailed balance. It is given by
[30, 31]

$\mathcal{Q}_{T}:=\frac{2\langle I\rangle^{2}}{\langle\\!\langle
I^{2}\rangle\\!\rangle\langle\sigma\rangle}\leq 1.$ (21)

In addition, several modified versions of the TUR have been derived, which
include a finite time generalization [71, 72], measurement and feedback
scenarios [73], as well as time-dependent periodic [74] and nonperiodic [75]
driving. In this work, however, we restrict our attention to the TUR given in
Eq. (21). In the steady-state regime, entropy production comes solely from the
dissipation of heat within the reservoirs [6],

$\langle\sigma\rangle=-\frac{J_{L}}{T_{L}}-\frac{J_{R}}{T_{R}},$ (22)

where $J_{\ell}$ is the steady state heat current from the reservoir $\ell$.
How these heat currents can be computed from the different models is discussed
in Appendix A. For reservoirs with equal temperature $T$, the entropy
production simplifies to

$\langle\sigma\rangle=\frac{(\mu_{L}-\mu_{R})\langle
I\rangle}{T}=\frac{eV\langle I\rangle}{T}.$ (23)

Moreover, Eq. 21 reduces to

$\mathcal{Q}_{T}:=\frac{2T\langle I\rangle}{eV\langle\\!\langle
I^{2}\rangle\\!\rangle}\leq 1,$ (24)

which only depends on the Fano factor $\langle\\!\langle
I^{2}\rangle\\!\rangle/\langle I\rangle$.

The KUR provides an alternative bound on the signal-to-noise ratio [33],

$\mathcal{Q}_{K}:=\frac{\langle I\rangle^{2}}{\langle\\!\langle
I^{2}\rangle\\!\rangle\langle\mathcal{A}\rangle}\leq 1,$ (25)

where $\langle\mathcal{A}\rangle$ is the dynamical activity that quantifies
the average rate of jumps the system undergoes. For a system described by a
classical rate equation, such as Eq. (16), the dynamical activity is defined
as [33]

$\langle\mathcal{A}\rangle=\sum_{i\neq j}W_{ij}p_{j}.$ (26)

It quantifies the average total rate of jumps between all states of the
system. In the quantum model, the process of electrons moving between the left
and right QDs is modeled by a coherent evolution, not classical jumps. For
this reason, the dynamical activity for Markovian quantum master equations is
usually defined such that it only quantifies the rate of jumps induced by
Lindblad jump operators [76, 77]. For the local master equation, this would
read

$\langle\mathcal{A}_{q}\rangle=\langle\mathcal{A}\rangle-\frac{4g^{2}}{\gamma_{L}+\gamma_{R}}(p_{L}+p_{R}),$
(27)

resulting in a strictly smaller dynamical activity. Employing
$\langle\mathcal{A}_{q}\rangle$ in the KUR could result in the following
problem: In a regime where the classical model and the master equation provide
the same probabilities and fluctuations (i.e., the system essentially behaves
classically), we could get a violation of the KUR simply because
$\langle\mathcal{A}_{q}\rangle$ does not take into account the transitions
between the left and right dots. To avoid such problems, we employ Eq. (26)
throughout for the dynamical activity, where the rates $W_{ij}$ are always
taken from the classical model (16), and the probabilities $p_{j}$ are
obtained with the respective method that is employed.

Figure 3: Coherence, current, and fluctuations for noninteracting electrons.
Results are shown for the local master equation (red), global master equation
(green), classical model (blue dashed), and NEGFs (black dashed). (a)
Coherence $|\alpha|$, (b) average current $\langle I\rangle$, and (c) current
fluctuations $\langle\\!\langle I^{2}\rangle\\!\rangle$, all as a function of
$g/\gamma$. Parameters: $\gamma_{L}=\gamma_{R}=\gamma$, $T/\gamma=100$,
$\mu_{L}=eV/2$, $\mu_{R}=-eV/2$, $eV/T=6.5$, $\epsilon=0$, $U=0$. Background
shadings mark approximate regions where the local master equation (red), the
global master equation (green), and the classical model (blue) are valid,
respectively.

Finally, we consider the TKUR that combines the above bounds as [34]

$\mathcal{Q}_{TK}:=\frac{\langle\\!\langle
I\rangle\\!\rangle^{2}}{\langle\\!\langle
I^{2}\rangle\\!\rangle}\frac{4\langle\mathcal{A}\rangle}{\langle\sigma\rangle^{2}}f\left(\frac{\langle\sigma\rangle}{2\langle\mathcal{A}\rangle}\right)^{2}\leq
1,$ (28)

where $f(x)$ is the inverse function of $x\tanh{x}$. This inequality optimizes
between the TUR and the KUR, and thus gives a stronger bound than both of
them.

When investigating violations of uncertainty relations, we use the quantifier

$\mathcal{V}_{j}=\text{max}\\{\mathcal{Q}_{j}-1,0\\},$ (29)

where $j=T,K,TK$ correspond to TUR, KUR, and TKUR, respectively. We refer to
this quantity as $j$ violation. If it is positive, the respective uncertainty
relation is violated. We note that numerically comparing the values of these
quantifiers between different types of violations is not meaningful.

##  IV Noninteracting electrons

In this section, we provide our results for noninteracting electrons. In this
case, there exist exact solutions obtained with NEGFs (cf. Appendix A.3),
which serve as a benchmark. Before illustrating the manifestations of
coherence, we consider the influence of the tunneling strength $g$ on
coherence, the current, and its fluctuations. This allows for illustrating the
coherent nature of the dynamics and to compare the different models. All
analytical expressions provided in this section are derived with the local
Lindblad master equation. For analytical results obtained with the global
approach we refer the reader to Appendix A.2.

###  IV.1 Coherent dynamics: Comparing the models

We first discuss how the coherence in the system changes with the interdot
tunnel coupling in the steady state. We find

$|\alpha|=\frac{2g|n_{L}-n_{R}|\gamma_{L}\gamma_{R}}{(\gamma_{L}+\gamma_{R})(4g^{2}+\gamma_{L}\gamma_{R})}.$
(30)

Interestingly, this expression exhibits a peak at $g\simeq\gamma$ (with
$\gamma_{L}=\gamma_{R}=\gamma$) [see Fig. 3(a)]. This peak in coherence can be
understood as follows: The coherent tunneling induces Rabi oscillations
between the left and the right dots, which may result in a buildup of
coherence. For $g\ll\gamma$, the decoherence induced by the bath suppresses
any buildup of coherence akin to the Zeno effect. For $g\gg\gamma$, the Rabi
oscillations become very fast, suppressing coherence by phase averaging. For
intermediate $g$, coherence can build up, resulting in a peak. We note that
the location of this peak does not depend on temperature or chemical potential
(which only enter in the Fermi-Dirac distributions $n_{\ell}$).

While the peak in coherence is captured well by the local master equation,
this approach breaks down when $g$ becomes of the order of $T$ [cf. Eq. (7)].
In this regime, where the global master equation is justified, coherence grows
again. As discussed above, this is a result of the entangled ground state
being more populated than the excited state. In the large $g$ limit, the
steady state reduces to the ground state, a pure singlet state with maximal
coherence. As illustrated in Fig. 3(a), we find excellent agreement between
the exact NEGF solution and the Lindblad master equation solutions in their
respective regime of validity.

Additionally, the current and its fluctuations provide insight into the effect
of coherence on the dynamics. For the current, the local master equation
yields

$\langle
I\rangle=\frac{4g^{2}(n_{L}-n_{R})\gamma_{L}\gamma_{R}(\gamma_{L}+\gamma_{R})}{(\gamma_{L}+\gamma_{R})^{2}(4g^{2}+\gamma_{L}\gamma_{R})}.$
(31)

The average current is illustrated in Fig. 3(b). Starting from small $g$, its
value grows as $g^{2}$ and saturates when the tunnel coupling $g$ becomes
larger than the system-bath coupling $\gamma$, which becomes the bottleneck
for transport. Similarly to the coherence, when $g$ is of the order of $T$ we
observe a breakdown of the local approach. For larger values of $g$, where the
global approach describes the dynamics well, the average current decreases, as
the eigenenergies of the Hamiltonian leave the bias window, and transport of
electrons is impeded. In addition, the plot includes the average current
predicted by the classical model (cf. Appendix A.4). It is in complete
agreement with the local master equation, and as such is accurate in the
regime of validity of the local approach.

The current fluctuations are given by

$\begin{split}\langle\\!\langle I^{2}\rangle\\!\rangle&=\langle
I\rangle\frac{n_{L}+n_{R}-2n_{L}n_{R}}{n_{L}-n_{R}}\\\ &-2\langle
I\rangle^{2}\left(\frac{1}{\gamma_{L}+\gamma_{R}}+\frac{\gamma_{L}+\gamma_{R}}{4g^{2}+\gamma_{L}\gamma_{R}}\right),\end{split}$
(32)

which is illustrated in Fig. 3(c). Comparison with the NEGF solution shows
that the local and global master equations reproduce the exact current
fluctuations in their respective regime of validity. To the best of our
knowledge, this is the first time that these regimes of validity have been
benchmarked by quantities that do not only depend on the quantum state alone
but instead depend on two-time correlation functions. As a function of
$g/\gamma$, the fluctuations display analogous features to the average.
However, when $g$ is of the same order as $\gamma$, the fluctuations are
suppressed below the value of the classical model (see also Ref. [15] for a
similar observation). The regime characterized by the peak in coherence thus
coincides with the regime where the fluctuations of the current are suppressed
relative to the classical model. These features underlie the manifestations of
coherence that we focus on below.

Figure 4: Comparison of coherence $|\alpha|$, concurrence $\mathcal{C}$, TUR
violation $\mathcal{V}_{T}$, KUR violation $\mathcal{V}_{K}$, and TKUR
Violation $\mathcal{V}_{TK}$ for noninteracting electrons. (a) Manifestations
of coherence as a function of $g/\gamma$ and $eV/T$. Results are obtained with
the local master equation. Parameters: $\gamma_{L}=\gamma_{R}=\gamma$,
$T/\gamma=100$, $\mu_{L}=eV/2$, $\mu_{R}=-eV/2$, $\epsilon=0$, $U=0$. For
clarity we only show the region, enveloped by the black dashed lines, where
$\mathcal{V}_{TK}>0$. [(b), (c)] Cross sections of panel (a) (gold and black
vertical lines) corresponding to $eV/T=6.5$ and $eV/T=11.5$, respectively.

###  IV.2 Manifestations of coherence

The expression for the concurrence [cf. Eq. (20)] is given by

$\displaystyle\mathcal{C}=\max\Bigg{\\{}0,~{}$ $\displaystyle
2|\alpha|-\frac{2\gamma_{L}\gamma_{R}}{4g^{2}+\gamma_{L}\gamma_{R}}\sqrt{\left[n_{L}n_{R}+\frac{4g^{2}}{\gamma_{L}\gamma_{R}}\bar{n}^{2}\right]}$
$\displaystyle\times\sqrt{\left[(1-n_{L})(1-n_{R})+\frac{4g^{2}}{\gamma_{L}\gamma_{R}}(1-\bar{n})^{2}\right]}\Bigg{\\}},$
(33)

where

$\bar{n}=\frac{n_{L}\gamma_{L}+n_{R}\gamma_{R}}{\gamma_{L}+\gamma_{R}},$ (34)

and $|\alpha|$ is given in Eq. (30). Maximal concurrence occurs for infinite
bias, $(\epsilon-\mu_{L})/T_{L}\to-\infty$ and
$(\epsilon-\mu_{R})/T_{R}\to\infty$; i.e., $eV\rightarrow\infty$. In this
regime the maximal value is $\mathcal{C}=(\sqrt{5}-1)/4\approx 0.31$ and can
be reached by setting $g/\gamma=(\sqrt{5}-1)/4$. The TUR violation, on the
other hand, is computed from Eq. (29) together with Eq. (24). It becomes
saturated at $g/\gamma=\sqrt{15}/6$ [16] irrespective of the applied bias,
similarly as it was found for the coherence. The corresponding maximum is
$\mathcal{V}_{T}\approx 0.141$, which occurs at $eV/T\approx 4.8$. Similarly,
we find KUR violations, i.e., $\mathcal{V}_{K}>0$. Its optimal value
$\mathcal{V}_{K}\approx 0.216$ occurs when $(\epsilon-\mu_{L})/T_{L}\to\infty$
and $(\epsilon-\mu_{R})/T_{R}\to-\infty$, which is the same regime saturating
concurrence, and when $g/\gamma\approx 0.59$.

All three manifestations are presented in Fig. 4 as a function of $g$ and
$eV$. The range of values of $g$ coincides with the window of the tunnel
couplings that exhibit both a peak of coherence and reduction of current
fluctuations. This is directly demonstrated in Figs. 4(b) and 4(c), where we
compare $|\alpha|$ with $\mathcal{C}$, $\mathcal{V}_{T}$, $\mathcal{V}_{K}$,
and $\mathcal{V}_{TK}$ on two cross sections of Fig. 4(a). Both plots show
excellent agreement between the local master equation (solid color lines) and
NEGFs (corresponding dashed black lines). The regime where $g$ and $\gamma$
are comparable thus features both entanglement, a static manifestation of
coherence describing the steady state, as well as TUR and KUR violations,
dynamical manifestations of coherence which arise from the nonclassical
dynamics of the system. Crucially, we see how these two manifestations of
coherence do not always appear together. As discussed above, in the regime of
strong $g$, entanglement is found. However, the dynamics of the system are
well described by the global master equation, a classical rate equation
involving the entangled eigenstates of the Hamiltonian. Consequently, we do
not find TUR and KUR violations.

There are also regimes where entanglement and TUR violations are found to
partially overlap, such that a crossover can be achieved by increasing the
voltage bias. Specifically, the TUR is only violated for small voltage bias,
whereas entanglement requires large voltage bias. Increasing $eV$ gives rise
to increasing entropy production, which causes the TUR to become less tight.
Entanglement grows with increasing current, which itself increases with the
voltage. For large voltages, we also observe KUR violations due to the
suppression of current fluctuations. In agreement with previous findings
$\mathcal{V}_{K}\approx 0.216$, TUR violations are observed close to
equilibrium while KUR violations appear far from equilibrium. As the TKUR is
tighter than both the TUR and the KUR, its violations encompass both the TUR
as well as the KUR violations.

We end this section by considering the usefulness of the entangled quantum
state to achieve nonlocality [24] and to perform quantum teleportation [78]
(see Appendix D). This has been done previously in a similar system [11]. For
noninteracting electrons, we find a teleportation fidelity of
$f=(7+\sqrt{5})/12\approx 0.77$ (cf. Appendix D.2), which is above the
classical limit of $2/3$. As expected, this corresponds to the value found in
Ref. [11] in the noninteracting case. Concerning nonlocality, in our system
with noninteracting electrons we do not find any violation of the Clauser-
Horne-Shimony-Holt (CHSH) inequality [79] (cf. Appendix D.1), in agreement
with the findings of Ref. [11]. As discussed in Sec. V.3, inclusion of Coulomb
interactions allows for a higher teleportation fidelity and even CHSH
violations, implying Bell nonlocality. This is in contrast to the results of
Ref. [11], where no chemical potential was present and population inversion
was obtained using negative temperatures.

##  V Interacting electrons: Enhanced coherence, stronger violations, and
nonlocal states

Figure 5: The effect of Coulomb interactions on coherence and concurrence. (a)
Coherence $|\alpha|$ as a function of $g/\gamma$ for different strength of
Coulomb interaction: $U/T=0$ (black), $U/T=3$ (blue), $U/T=5$ (green), and
$U\to\infty$ (red). The plot was produced with the results of the semilocal
master equation. Parameters: $T/\gamma=100$, $eV/T=6.5$, $\epsilon=0$. [(b),
(c)] Coherence $|\alpha|$ and concurrence $\mathcal{C}$ as a function of
$U/g$, respectively, obtained with the semilocal master equation (black) and
NEGFs (red). Parameters: $g/\gamma=0.3$, $T/\gamma=10$, $eV/T=3.6$, and
$\epsilon=0$.

Coulomb interactions between the electrons change the dynamics of the double
QD system by incorporating an additional energetic cost for occupying both
QDs. In this section we show how the presence of interactions affects the
manifestations of coherence.

###  V.1 Methods

We investigate the dynamics of the system using methods analogous to the
noninteracting case, adjusted to take into account the effect of Coulomb
interactions between electrons.

####  V.1.1 Master equation

In the previous section we have identified the regime of interest to be where
the peak of coherence appears together with all the manifestations of
coherence. Since the same holds for interacting electrons, we restrict our
attention to the regime where tunnel couplings $g$ are comparable to
$\gamma_{\ell}$. A Lindblad master equation suitable for treating this regime
in the presence of Coulomb interactions (termed semilocal) has recently been
microscopically derived [36]. The corresponding superoperators read

$\begin{split}\mathcal{L}_{\ell}&=\gamma_{\ell}\left[n_{\ell}\mathcal{D}[(1-\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}})\hat{c}_{\ell}^{\dagger}]+(1-n_{\ell})\mathcal{D}[(1-\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}})\hat{c}_{\ell}]\right]\\\
&+\gamma_{\ell}\left[n_{\ell}^{U}\mathcal{D}[\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}}\hat{c}_{\ell}^{\dagger}]+(1-n_{\ell}^{U})\mathcal{D}[\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}}\hat{c}_{\ell}]\right],\end{split}$
(35)

where

$n_{\ell}:=n_{\ell}(\epsilon),$ (36) $n_{\ell}^{U}:=n_{\ell}(\epsilon+U),$
(37)

and $\bar{\ell}\neq\ell$. It is valid in the regime

$g\ll\max\\{T_{\ell},|\epsilon+U-\mu_{\ell}|\\}\qquad\text{(semilocal ME)}.$
(38)

The key difference between the semilocal master equation and the local one,
which we used for the noninteracting case, is that now the jump rates from
each QD to the adjacent reservoir depend on the occupation of the other QD.
Importantly, the results from the semilocal master equation reduce to the
results from the local master equation in the limit $U\rightarrow 0$ [36].

In order to compute the average current and its fluctuations we proceed
analogously to the case of noninteracting particles. The average current in
the steady state is given by Eq. (11), and the fluctuations are obtained with
full counting statistics (Appendix B).

####  V.1.2 Classical model

As for noninteracting electrons, we derived a classical model of the form of
Eq. (16) using the perturbative approach described in Appendix A.4. The matrix
elements describing jumps of electrons between the reservoirs and the system
read $W_{L0}=\gamma_{L}n_{L}$, $W_{0L}=\gamma_{L}(1-n_{L})$,
$W_{DR}=\gamma_{L}n_{L}^{U}$, $W_{RD}=\gamma_{L}(1-n_{L}^{U})$, and similarly
for $L\leftrightarrow R$. These are the same rates that appear in the
semilocal master equation [cf. Eq. (35)]. The interdot tunnel rate is given by

$W_{LR}=\frac{4g^{2}}{\gamma_{L}(1-n_{L}+n_{L}^{U})+\gamma_{R}(1-n_{R}+n_{R}^{U})},$
(39)

and $W_{RL}=W_{LR}$.

####  V.1.3 NEGFs

The effect of Coulomb interactions in the NEGF can be introduced via a many-
body self-energy. This is an in-principle exact procedure. In practice, one
must, however resort to approximate self-energies, calculated, e.g., via many-
body approximations. Here, we consider a rather simple approximation, namely
the second Born perturbative scheme [80, 81, 82] (for details we refer the
reader to Appendix C.2). The range of $U/g$ where the second Born
approximation gives a reliable description of electronic correlations is
limited and much smaller than the range we consider with the Lindblad master
equation in the remainder of the paper. Nonetheless, it may still be useful to
compare the two methods in a restricted domain where both approaches produce
meaningful results.

We thus notice that in this limited $U/g$ range the two treatments give
results that are in good mutual agreement for the coherences, while
discernible discrepancies occur for the concurrence, due to well-known
possible shortcomings of the second Born approximation when determining double
occupancies (see Appendix C.2).

Figure 6: Evolution of concurrence $\mathcal{C}$, TUR violation
$\mathcal{V}_{T}$, KUR violation $\mathcal{V}_{K}$, and TKUR violation
$\mathcal{V}_{TK}$ for different strengths of Coulomb interaction: (a)
$U/T=0$, (b) $U/T=3$, (c) $U/T=5$, and (d) $U/T=10$. Plots as a function of
$g/\gamma$ and $eV/T$ were obtained with the local master equation.
Parameters: $\gamma_{L}=\gamma_{R}=\gamma$, $T/\gamma=100$, $\mu_{L}=eV/2$,
$\mu_{R}=-eV/2$, $\epsilon=0$. Dashed black lines envelop the regions of TKUR
violations.

###  V.2 Enhanced manifestations of coherence

Since Coulomb interactions introduce an energetic cost of having both QDs
simultaneously filled, the probability of double occupation, $p_{D}$, is
suppressed and vanishes in the limit $U\to\infty$. This reduced occupation of
the doubly filled subspace generally has a beneficial effect on the amount of
coherence, as it increases the probability for the subspace with a single
electron on the dot. However, as explained below, Coulomb interactions do not
always result in an increased amount of coherence. The expression for
coherence reads

$\begin{split}|\alpha|=&\frac{2g\gamma_{L}\gamma_{R}}{\gamma_{L}+\gamma_{R}}\\\
&\times\frac{n_{L}(1-n_{R}^{U})-n_{R}(1-n_{L}^{U})}{4g^{2}(1+\bar{n}-\bar{n}_{U})+(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})},\end{split}$
(40)

where we have introduced

$\bar{n}_{U}=\frac{\gamma_{L}n_{L}^{U}+\gamma_{R}n_{R}^{U}}{\gamma_{L}+\gamma_{R}},$
(41)

and $\Delta_{l}=n_{l}-n_{l}^{U}$. In the limit $U\to\infty$ the coherence
reduces to

$|\alpha|=\frac{2g\gamma_{L}\gamma_{R}\left(n_{L}-n_{R}\right)}{(\gamma_{L}+\gamma_{R})\left[4g^{2}(1+\bar{n})+(1-n_{L}n_{R})\gamma_{L}\gamma_{R}(1-\bar{n})\right]}.$
(42)

As illustrated in Fig. 5(a), increasing $U$ results in a larger peak in
coherence that is shifted to smaller values of $g/\gamma$. Because of this
shift, there are values of $g/\gamma$ where coherence decreases when
increasing $U$.

From the definition of concurrence, given by Eq. (20), we anticipate that
Coulomb interactions increase the amount of entanglement in the system for two
reasons: First, because $p_{D}$ is smaller, and second, because this may
increase the off-diagonal element $\alpha$. This is illustrated in Figs. 5(b)
and 5(c), where we survey the influence of $U/g$ on the coherence and
concurrence, respectively. In Fig. 5(b) coherence shows moderate growth with
the strength of the Coulomb interaction. For comparison, we feature the
results of the semilocal master equation (black) and of NEGFs (red) on the
restricted interval, demonstrating a good agreement between the two methods. A
similar plot showing the concurrence is presented in Fig. 5(c). For
concurrence we observe a strong impact of Coulomb interactions, which
indicates that the growth of concurrence comes predominantly from the
suppression of $p_{D}$. We note that the agreement between the semilocal
master equation and NEGFs is better for coherence than for concurrence. We
attribute this to the fact that the concurrence relies not only on single
particle quantities but requires the computation of $p_{D}$. As discussed in
Appendix C.2, this quantity may be more sensitive to the employed
approximations.

In addition to entanglement, Coulomb interactions impact the dynamical
manifestations, TUR and KUR violations, which is illustrated in Fig. 6. In the
series of panels we show how $\mathcal{C}$, $\mathcal{V}_{T}$, and
$\mathcal{V}_{K}$ evolve as $U/T$ increases. Concurrence and TUR violations
are considerably enhanced by interactions, while the enhancement of KUR
violations is less prominent. While the range of $g/\gamma$ where entanglement
is present is extended, the same is not true for the remaining manifestations.
This is a consequence of the fact that larger $\mathcal{C}$ is caused
predominantly by the reduction of $p_{D}$ as opposed to an increase in
coherence. In addition, the values of $g/\gamma$ that maximize
$\mathcal{V}_{T}$ and $\mathcal{V}_{K}$ become smaller and align with the
optimal value for $\mathcal{C}$. This is similar to the observed shift in the
peak of coherence.

###  V.3 Nonlocality and entanglement for an interacting double uqantum dot

The main drawback of looking at the entanglement solely through the prism of
any entanglement measure, such as concurrence, is that it does not provide a
definite measure of the usefulness to achieve nonlocality [24] or perform
quantum information tasks, such as quantum teleportation [78]. The problem of
generating operational nonclassicality in the class of systems akin to the
double quantum dot was systematically analyzed by Bohr Brask et al. [11],
where no nonlocality was found. In their approach, however, the system was
brought out of equilibrium via a temperature gradient alone but allowing for
negative temperatures to describe population inversion. Here we show that, by
considering a bias voltage, Coulomb interactions may result in nonlocality.

The optimal regime to maximize concurrence and observe nonlocality is provided
by the following conditions: First, jumps from the right reservoir into the
system are completely suppressed, i.e., $n_{R}=n_{R}^{U}=0$. Second, jumps
from the left reservoir are suppressed only when there is already an electron
in the system, i.e., $n_{L}^{U}=0$ and $n_{L}=1$. This is achieved in the
limits $U\to\infty$ and $eV\to\infty$, with $U/eV\to\infty$.

Under these conditions, the probability of the double occupation $p_{D}$
vanishes and the coherence reduces to

$|\alpha|=\frac{2g\gamma_{L}\gamma_{R}}{\gamma_{L}\gamma_{R}^{2}+4g^{2}(2\gamma_{L}+\gamma_{R})}.$
(43)

The expression for concurrence simplifies to

$\mathcal{C}=\frac{4g\gamma_{L}\gamma_{R}}{\gamma_{L}\gamma_{R}^{2}+4g^{2}(2\gamma_{L}+\gamma_{R})}.$
(44)

Upon examining $\frac{\partial{\mathcal{C}}}{{\partial{g}}}$, we find the
maximal value $\mathcal{C}=\sqrt{2}/2\approx 0.71$ when
$\gamma_{R}/g=2\sqrt{2}$ and $\gamma_{L}/g\to\infty$. This is significantly
larger than $(\sqrt{5}-1)/4\approx 0.31$, which we found to be the maximum in
the noninteracting case. If the couplings to the baths are assumed to be
equal, the largest concurrence occurs at $\gamma/g=2\sqrt{3}$, resulting in
$\mathcal{C}=\sqrt{3}/3\approx 0.58$. We note that in the absence of a
chemical potential (even allowing for negative temperatures), Coulomb
interactions do not alter the maximal value for the concurrence [11].

Motivated by the substantial increase in concurrence due to Coulomb
interactions, we revisit the question of Bell nonlocality (cf. Appendix D.1),
focusing on the CHSH scenario [79], which was also considered in Ref. [11]. In
this scenario, if the statistics of the outcomes of measurements on a quantum
state admit a local hidden-variable model, then the CHSH quantity obeys
$\text{CHSH}\leq 2$ [24]. Any violation of this bound is a signature of
nonlocality, and the maximal $\text{CHSH}=2\sqrt{2}$ can be achieved with
maximally entangled Bell states. For the double quantum dot, in the regime
discussed above, we find the maximal value $\text{CHSH}=2\sqrt{3/2}$, which is
achieved for the same set of parameters $\gamma_{L},\gamma_{R}$, and $g$ that
saturate concurrence. When $\gamma_{R}/g=2\sqrt{2}$, CHSH displays nonlocality
for $\gamma_{L}/g>5.66$.

For the optimal fidelity of quantum teleportation [83] we find
$f=(4+\sqrt{2})/6\approx 0.9$, which outperforms the system without Coulomb
interactions, where $f=(7+\sqrt{5})/12\approx 0.77$, which coincides with the
maximum found in Ref. [11]. The details of the calculations are outlined in
Appendix D.2.

##  VI More general models

Figure 7: TUR violations and coherence in the chain of three QDs: (a) TUR
violation, $\mathcal{V}_{T}$, and (b) one-half of the $\ell_{1}$-norm of
coherence, $\frac{1}{2}\mathcal{C}_{\ell_{1}}$, as a function of two coherent
couplings $g_{1}/\gamma$ and $g_{2}/\gamma$; c) diagonal cross-sections
corresponding to $g_{1}=g_{2}=g$. Parameters: $eV/T=5$, $T/\gamma=100$,
$\epsilon=0$.

Up to this point we have been concerned with the minimal model consisting of
only two sites, which allows to clearly identify the regimes of nonclassical
behavior. The fact that similar conclusions hold for both interacting and
noninteracting electrons implies that our results may serve to anticipate the
behavior of more general mesoscopic systems. To illustrate this, we
investigate TUR violations in a chain of three QDs, where electrons in the
system may coherently traverse between the neighboring QDs. The corresponding
Hamiltonian is given by

$\hat{H}=\sum_{\ell=L,M,R}\epsilon\hat{c}_{\ell}^{\dagger}c_{\ell}+g_{1}(\hat{c}_{L}^{\dagger}c_{M}+\hat{c}_{M}^{\dagger}c_{L})+g_{2}(\hat{c}_{M}^{\dagger}c_{R}+\hat{c}_{R}^{\dagger}c_{M}),$
(45)

where the subscript $M$ refers to the middle QD, and $g_{1}$ ($g_{2}$) denotes
the coherent coupling between the left and the middle (the right and the
middle) QDs. We investigate the chain with the NEGFs using analogous methods
as those presented in Appendix A.3. As anticipated from our simple model, we
find nonclassical behavior where both coherent couplings, $g_{1}$ and $g_{2}$,
are of the order of the coupling to the baths $\gamma$. This is illustrated in
Figs. 7(a) and 7(b), which show TUR violation $\mathcal{V}_{T}$ and the amount
of coherence in the system, respectively. As a quantifier of the amount of
coherence we use one-half of the $\ell_{1}$-norm of coherence, here denoted by
$C_{\ell_{1}}$, which reduces to $|\alpha|$ in the two-site model. The peaks
of $\mathcal{V}_{T}$ or $C_{\ell_{1}}$ are visualized in Fig. 7(c), where we
show the diagonal cross-sections with equal coherent couplings,
$g_{1}=g_{2}=g$, of Figs. 7(a) (black) and 7(b) (red), respectively. These
graphs further demonstrate that the TUR violation and coherence are maximized
for $g_{1}=g_{2}$ both being comparable to $\gamma$. Moreover, in the chain of
three QDs we observe stronger TUR violation than in a double quantum dot (DQD)
operating under the same set of parameters (see Fig. 4). This indicates that
nonclassical behavior is more easily accessible in longer chains. The fact
that our conclusions, drawn from the investigation of the DQD, hold in the
chain of three QDs provides evidence for the applicability of these results in
more general systems. We note that these results are also consistent with Ref.
[84, 85], where TUR violations in longer chains of QDs were investigated.

##  VII Conclusions and outlook

In this paper we investigated the nonclassical behavior of a mesoscopic open
quantum system by analyzing two complementary manifestations of coherence:
Entanglement, a static manifestation, and violations of thermokinetic
uncertainty relations, dynamical manifestations of coherence. We considered a
double quantum dot weakly coupled to fermionic reservoirs, where the transport
of electrons between the dots can be modeled either by coherent tunneling, or
by stochastic jumps.

We identified the regime where the system exhibits truly nonclassical
behavior: On the one hand, coherence exhibits a peak, which may result in an
entangled state and, in the presence of interactions, even nonlocality. On the
other hand, coherent tunneling suppresses the fluctuations of the current
below the classical model, allowing for violations of the TUR and the KUR.
This behavior occurs in the regime where the coherent tunnel-coupling is of
the order of the coupling to the baths, which is captured by the local master
equation. This is different from the regime of large tunnel couplings, which
is captured by the global master equation, where the state of the system may
be strongly entangled, but the dynamics of the system can be expressed through
a classical rate equation, and no TUR or KUR violations are present.

Our systematic characterization of the electron transport across the system
performed with different master equations and NEGFs illustrates the need to go
beyond the steady state in order to fully capture the nonclassical behavior in
mesoscopic transport. While the steady state captures the entanglement, the
violations of thermokinetic uncertainty relations cannot be captured by the
steady state alone. Those different manifestations of coherence are not
equivalent and do not need to coincide. We note that in the regime of
nonclassical behavior, not only the quantum state of the DQD but also the
current fluctuations, which depend on two-time correlations, are well captured
by the local master equation. While our results are obtained for a specific
model, we believe that these conclusions are relevant for a broad range of
open quantum systems and provide guiding principles for the design of out-of-
equilibrium devices that exhibit nonclassical transport. This is evidenced by
the fact that our predictions about the nonclassical behavior, drawn from the
noninteracting two-site model, also hold for interacting electrons and in the
chain of three QDs.

Our results may be probed experimentally in a number of different platforms,
including electronic QDs, which can be implemented, e.g., using two-
dimensional electron gases [86, 87] or nanowires [88, 89], as well as
superconducting qubits [90, 91, 92]. While measurements of the electric
current and its fluctuations across a DQD is a standard experimental technique
[93], accessing the quantum state of the DQD is more challenging.
Alternatively, the density matrix and, therefore, the entanglement may be
probed in an implementation of our model using two coupled superconducting
qubits [90].

###### Acknowledgements.

We acknowledge contributions by Meeri Mölsä during her bachelor thesis in the
early stages of this work. P.P.P. and K.P. acknowledge funding from the Swiss
National Science Foundation (Eccellenza Professorial Fellowship
PCEFP2_194268). P.S. and C.V. acknowledge support from the Swedish Research
Council, Grants No. 2018-03921 and No. 2017-03945, respectively.

##  Appendix A Noninteracting electrons

In the main text of the paper we restrict our attention to the case with the
same onsite energies of both QDs. In this appendix, however, we provide a
derivation of the results in the general case of noninteracting electrons. The
Hamiltonian operator of the system is given by

$\hat{H}=\sum_{\ell=L,R}\epsilon_{\ell}\hat{c}_{\ell}^{\dagger}\hat{c}_{\ell}+g\left(\hat{c}_{L}^{\dagger}\hat{c}_{R}+\hat{c}_{R}^{\dagger}\hat{c}_{L}\right),$
(46)

where $\epsilon_{\ell}$ is the on-site energy of the QD $\ell$. In the
noninteracting case, the coefficients of density matrix (19) can be expressed
as

$\begin{split}&p_{0}=\left(1-\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle\right)\left(1-\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle\right)-\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle,\\\
&p_{L}=\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle\left(1-\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle\right)+\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle,\\\
&p_{R}=\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle\left(1-\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle\right)+\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle,\\\
&p_{D}=\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle-\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle,\\\
&\alpha=\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle,\end{split}$ (47)

where we have used a corollary of Wick’s theorem [94],

$\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle=\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle-\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle.$
(48)

###  A.1 Local master equation

In the local master equation satisfying thermodynamic consistency the Lindblad
superoperators appearing in Eq. (3) are given by [36]

$\mathcal{L}_{\ell}=n_{\ell}\gamma_{\ell}\mathcal{D}[\hat{c}^{\dagger}_{\ell}]+(1-n_{\ell})\gamma_{\ell}\mathcal{D}[\hat{c}_{\ell}],$
(49)

where now $n_{\ell}:=n_{\ell}(\bar{\epsilon})$, and
$\bar{\epsilon}=(\epsilon_{L}+\epsilon_{R})/2$. In order to find the density
matrix in the steady state, we consider the vector
$\vec{v}=[\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle,\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle,\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle,\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle]^{T}$.
Its equation of motion reads

$\frac{d}{dt}\vec{v}=A\vec{v}+\vec{b},$ (50)

where

$A=\begin{pmatrix}-\gamma_{L}&0&-ig&ig\\\ 0&-\gamma_{R}&ig&-ig\\\
-ig&ig&2i\delta-\frac{\gamma_{L}+\gamma_{R}}{2}&0\\\
ig&-ig&0&-2i\delta-\frac{\gamma_{L}+\gamma_{R}}{2}\end{pmatrix},$ (51)

$\vec{b}=[\gamma_{L}n_{L},\gamma_{R}n_{R},0,0]^{T}$, and
$\delta=(\epsilon_{L}-\epsilon_{R})/2$. The steady-state solution
$-A^{-1}\vec{b}$ then yields

$\begin{split}&\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle=\frac{4g^{2}(\gamma_{L}+\gamma_{R})(n_{L}\gamma_{L}+n_{R}\gamma_{R})+n_{L}\gamma_{L}\gamma_{R}\left((\gamma_{L}+\gamma_{R})^{2}+16\delta^{2}\right)}{(\gamma_{L}+\gamma_{R})^{2}(\gamma_{L}\gamma_{R}+4g^{2})+16\gamma_{L}\gamma_{R}\delta^{2}},\\\
&\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle=\frac{4g^{2}(\gamma_{L}+\gamma_{R})(n_{L}\gamma_{L}+n_{R}\gamma_{R})+n_{R}\gamma_{L}\gamma_{R}\left((\gamma_{L}+\gamma_{R})^{2}+16\delta^{2}\right)}{(\gamma_{L}+\gamma_{R})^{2}(\gamma_{L}\gamma_{R}+4g^{2})+16\gamma_{L}\gamma_{R}\delta^{2}},\\\
&\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle=\frac{-2ig(n_{L}-n_{R})\gamma_{L}\gamma_{R}(\gamma_{L}+\gamma_{R}+4i\delta)}{(\gamma_{L}+\gamma_{R})^{2}(\gamma_{L}\gamma_{R}+4g^{2})+16\gamma_{L}\gamma_{R}\delta^{2}},\\\
&\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle=\frac{2ig(n_{L}-n_{R})\gamma_{L}\gamma_{R}(\gamma_{L}+\gamma_{R}-4i\delta)}{(\gamma_{L}+\gamma_{R})^{2}(\gamma_{L}\gamma_{R}+4g^{2})+16\gamma_{L}\gamma_{R}\delta^{2}}.\end{split}$
(52)

Using Eq. (47) we obtain

$\begin{split}&p_{0}=\frac{4g^{2}(1-\bar{n})^{2}+(1-n_{L})(1-n_{R})\gamma_{L}\gamma_{R}\left(1+\xi^{2}\right)}{4g^{2}+\gamma_{L}\gamma_{R}(1+\xi^{2})},\\\
&p_{D}=\frac{4g^{2}\bar{n}^{2}+n_{L}n_{R}\gamma_{L}\gamma_{R}\left(1+\xi^{2}\right)}{4g^{2}+\gamma_{L}\gamma_{R}(1+\xi^{2})},\\\
&p_{L}=\frac{4g^{2}\bar{n}(1-\bar{n})+n_{L}(1-n_{R})\gamma_{L}\gamma_{R}\left(1+\xi^{2}\right)}{4g^{2}+\gamma_{L}\gamma_{R}(1+\xi^{2})},\\\
&p_{R}=\frac{4g^{2}\bar{n}(1-\bar{n})+(1-n_{L})n_{R}\gamma_{L}\gamma_{R}\left(1+\xi^{2}\right)}{4g^{2}+\gamma_{L}\gamma_{R}(1+\xi^{2})},\\\
&\alpha=\frac{2ig(n_{L}-n_{R})\gamma_{L}\gamma_{R}(1+i\xi)}{(\gamma_{L}+\gamma_{R})(4g^{2}+\gamma_{L}\gamma_{R}(1+\xi^{2}))},\\\
\end{split}$ (53)

where

$\bar{n}=\frac{n_{L}\gamma_{L}+n_{R}\gamma_{R}}{\gamma_{L}+\gamma_{R}}\hskip
19.91684pt\text{and}\hskip
14.22636pt\xi=\frac{4\delta}{\gamma_{L}+\gamma_{R}}.$ (54)

Using the formula for the concurrence, given by Eq. (20), we obtain the
following expression,

$\mathcal{C}=2|\alpha|-\frac{2\gamma_{L}\gamma_{R}}{4g^{2}+\gamma_{L}\gamma_{R}(1+\xi^{2})}\sqrt{\left[(1-n_{L})(1-n_{R})(1+\xi^{2})+\frac{4g^{2}}{\gamma_{L}\gamma_{R}}(1-\bar{n})^{2}\right]\left[n_{L}n_{R}(1+\xi^{2})+\frac{4g^{2}}{\gamma_{L}\gamma_{R}}\bar{n}^{2}\right]},$
(55)

if this quantity is positive and zero otherwise.

####  A.1.1 Average current

In the steady state, the average current passing through the system from the
left reservoir can be obtained from the expression given in Eq. (11).
Equivalently, one may compute the current as an expectation value of the
current operator,

$\hat{I}=ig(\hat{c}_{L}^{\dagger}\hat{c}_{R}-\hat{c}_{R}^{\dagger}\hat{c}_{L}).$
(56)

In the steady state, we find

$\langle
I\rangle=\frac{4g^{2}(n_{L}-n_{R})\gamma_{L}\gamma_{R}(\gamma_{L}+\gamma_{R})}{(\gamma_{L}+\gamma_{R})^{2}(4g^{2}+\gamma_{L}\gamma_{R})+16\gamma_{L}\gamma_{R}\delta^{2}}.$
(57)

####  A.1.2 Heat current

We define the heat currents in the thermodynamically consistent framework [36,
40], which guarantees the thermodynamically consistent heat flow and maintains
the accuracy of the local form of the master equation. In this framework the
average heat current from the reservoir $\ell$ to the double quantum dot is
given by [36]

$J_{\ell}=\text{Tr}\left[(\hat{H}_{TD}-\mu_{\ell}\hat{N})\mathcal{L}_{\ell}\hat{\rho}\right],$
(58)

where we have introduced a thermodynamic bookkeeping Hamiltonian

$\hat{H}_{TD}=\bar{\epsilon}(\hat{c}_{L}^{\dagger}\hat{c}_{L}+\hat{c}_{R}^{\dagger}\hat{c}_{R}).$
(59)

We find

$J_{L}=\text{Tr}\left[\left(\bar{\epsilon}-\mu_{L}\right)\hat{N}\mathcal{L}_{L}\hat{\rho}\right]=\left(\bar{\epsilon}-\mu_{L}\right)\langle
I\rangle,$ (60)

and

$J_{R}=\text{Tr}\left[\left(\bar{\epsilon}-\mu_{R}\right)\hat{N}\mathcal{L}_{R}\hat{\rho}\right]=-\left(\bar{\epsilon}-\mu_{R}\right)\langle
I\rangle.$ (61)

We note that if $\epsilon_{L}=\epsilon_{R}$, which is assumed in the main text
of this paper, the steady-state heat currents reduce to the ones obtained from
the conventional definition, where in Eq. (58) the Hamiltonian $\hat{H}$ is
used instead of $\hat{H}_{TD}$.

####  A.1.3 Current fluctuations

In order to investigate the fluctuations of the current, we resort to the
method of full counting statistics. All the details about this procedure are
explained in Appendix B.

The counting-field-dependent Liouvillian $\mathcal{L}(\chi)$, introduced in
Eq. (112), is obtained by inserting the terms $e^{\pm i\chi}$ into the
Lindblad superoperator [Eq. (49)] of the local master equation,

$\mathcal{L}(\chi)\hat{\rho}(\chi,t)=\mathcal{L}\hat{\rho}(\chi,t)+\gamma_{L}n_{L}(e^{i\chi}-1)\hat{c}_{L}^{\dagger}\hat{\rho}(\chi,t)\hat{c}_{L}+\gamma_{L}(1-n_{L})(e^{-i\chi}-1)\hat{c}_{L}\hat{\rho}(\chi,t)\hat{c}_{L}^{\dagger},$
(62)

where $\mathcal{L}$ is the Liouvillian in the absence of a counting field
given in Eq. (3). Writing $\mathcal{L}(\chi)$ in the basis
$\\{p_{0},p_{L},p_{R},p_{D},\text{Re}[\alpha],\text{Im}[\alpha]\\}$
corresponding to the elements of the density matrix [cf. Eq. (19)] yields

$\mathbb{L}(\chi)=\footnotesize\begin{pmatrix}-n_{L}\gamma_{L}-n_{R}\gamma_{R}&e^{-i\chi}(1-n_{L})\gamma_{L}&(1-n_{R})\gamma_{R}&0&0&0\\\
e^{i\chi}n_{L}\gamma_{L}&-(1-n_{L})\gamma_{L}-n_{R}\gamma_{R}&0&(1-n_{R})\gamma_{R}&0&-2g\\\
n_{R}\gamma_{R}&0&-n_{L}\gamma_{L}-(1-n_{R})\gamma_{R}&e^{-i\chi}(1-n_{L})\gamma_{L}&0&2g\\\
0&n_{R}\gamma_{R}&e^{i\chi}n_{L}\gamma_{L}&-(1-n_{L})\gamma_{L}-(1-n_{R})\gamma_{R}&0&0\\\
0&0&0&0&-\frac{\gamma_{L}+\gamma_{R}}{2}&2\delta\\\
0&g&-g&0&-2\delta&-\frac{\gamma_{L}+\gamma_{R}}{2}\end{pmatrix}\normalsize.$
(63)

The expression for the current fluctuations is given by

$\langle\\!\langle I^{2}\rangle\\!\rangle=\langle
I\rangle\frac{n_{L}+n_{R}-2n_{L}n_{R}}{n_{L}-n_{R}}-\frac{2\langle
I\rangle^{2}}{\gamma_{L}+\gamma_{R}}\left(1+\frac{(\gamma_{L}+\gamma_{R})^{2}+16\delta^{2}(\gamma_{L}-\gamma_{R})^{2}/(\gamma_{L}+\gamma_{R})^{2}}{4g^{2}+\gamma_{L}\gamma_{R}+16\delta^{2}\gamma_{L}\gamma_{R}/(\gamma_{L}+\gamma_{R})^{2}}\right).$
(64)

###  A.2 Global master equation

In the global master equation electrons jump into the eigenmodes of the
Hamiltonian [cf. (46)], which can be written in diagonal form,

$\hat{H}=\epsilon_{+}\hat{c}_{+}^{\dagger}\hat{c}_{+}+\epsilon_{-}\hat{c}_{-}^{\dagger}\hat{c}_{-},$
(65)

where we have introduced the annihilation operators

$\begin{split}&\hat{c}_{+}=\cos{(\theta/2)}\hat{c}_{L}+\sin{(\theta/2)}\hat{c}_{R},\\\
&\hat{c}_{-}=-\sin{(\theta/2)}\hat{c}_{L}+\cos{(\theta/2)}\hat{c}_{R},\end{split}$
(66)

and the corresponding eigenenergies

$\epsilon_{\pm}=\bar{\epsilon}\pm\sqrt{\delta^{2}+g^{2}}.$ (67)

Here $\bar{\epsilon}=(\epsilon_{L}+\epsilon_{R})/2$, and the angle $\theta$
($0\leq\theta\leq\pi$) is defined by
$\cos{(\theta)}=\delta/\sqrt{g^{2}+\delta^{2}}$. The Lindblad superoperators
in the global master equation are given by [36]

$\mathcal{L}_{\ell}=\sum_{s=\pm}n_{\ell}^{s}\gamma_{\ell}^{s}\mathcal{D}[\hat{c}_{s}^{\dagger}]+(1-n_{\ell}^{s})\gamma_{\ell}^{s}\mathcal{D}[\hat{c}_{s}],$
(68)

where we have introduced the rescaled coupling rates

$\begin{split}&\gamma_{L}^{+}=\gamma_{L}\cos^{2}{(\theta/2)}\text{,
}\gamma_{L}^{-}=\gamma_{L}\sin^{2}{(\theta/2)},\\\
&\gamma_{R}^{+}=\gamma_{R}\sin^{2}{(\theta/2)}\text{,
}\gamma_{R}^{-}=\gamma_{R}\cos^{2}{(\theta/2)},\end{split}$ (69)

and

$n_{\ell}^{s}:=n_{\ell}(\epsilon_{s}).$ (70)

From the global master equation, we find the equations of motion

$\frac{d}{dt}\langle\hat{c}_{s}^{\dagger}\hat{c}_{s}\rangle=(\gamma_{L}^{s}n_{L}^{s}+\gamma_{R}^{s}n_{R}^{s})-(\gamma_{L}^{s}+\gamma_{R}^{s})\langle\hat{c}_{s}^{\dagger}\hat{c}_{s}\rangle.$
(71)

The corresponding steady-state solution is given by

$\langle\hat{c}_{s}^{\dagger}\hat{c}_{s}\rangle=\frac{\gamma_{L}^{s}n_{L}^{s}+\gamma_{R}^{s}n_{R}^{s}}{\gamma_{L}^{s}+\gamma_{R}^{s}}=:\bar{n}_{s},$
(72)

and the remaining expectation values vanish, i.e.,
$\langle\hat{c}_{+}^{\dagger}\hat{c}_{-}\rangle=\langle\hat{c}_{-}^{\dagger}\hat{c}_{+}\rangle=0$.

The density matrix elements in the local basis can be obtained from Eq. (47),
using the relations

$\begin{split}&\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle=\sin^{2}(\theta/2)\langle\hat{c}_{-}^{\dagger}\hat{c}_{-}\rangle+\cos^{2}(\theta/2)\langle\hat{c}_{+}^{\dagger}\hat{c}_{+}\rangle,\\\
&\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle=\cos^{2}(\theta/2)\langle\hat{c}_{-}^{\dagger}\hat{c}_{-}\rangle+\sin^{2}(\theta/2)\langle\hat{c}_{+}^{\dagger}\hat{c}_{+}\rangle,\\\
&\langle\hat{c}_{L}^{\dagger}\hat{c}_{R}\rangle=\sin{(\theta/2)}\cos{(\theta/2)}\left(\langle\hat{c}_{+}^{\dagger}\hat{c}_{+}\rangle-\langle\hat{c}_{-}^{\dagger}\hat{c}_{-}\rangle\right),\\\
&\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle=\sin{(\theta/2)}\cos{(\theta/2)}\left(\langle\hat{c}_{+}^{\dagger}\hat{c}_{+}\rangle-\langle\hat{c}_{-}^{\dagger}\hat{c}_{-}\rangle\right).\end{split}$
(73)

The coherence then reads

$|\alpha|=\frac{g}{2\sqrt{g^{2}+\delta^{2}}}|\bar{n}_{+}-\bar{n}_{-}|,$ (74)

whereas the concurrence is given by

$\begin{split}\mathcal{C}=2|\alpha|-2\sqrt{\bar{n}_{+}\bar{n}_{-}(1-\bar{n}_{+})(1-\bar{n}_{-})}.\end{split}$
(75)

####  A.2.1 Average current

As for the local master equation, the average may be obtained from Eq. (11).
We note that for the global master equation, we may not use a current operator
[35]. We find

$\langle
I\rangle=\frac{(n_{L}^{+}-n_{R}^{+})\gamma_{L}^{+}\gamma_{R}^{+}}{\gamma_{L}^{+}+\gamma_{R}^{+}}+\frac{(n_{L}^{-}-n_{R}^{-})\gamma_{L}^{-}\gamma_{R}^{-}}{\gamma_{L}^{-}+\gamma_{R}^{-}}.$
(76)

####  A.2.2 Heat current

The heat current from reservoir $\ell$ may be written as [36]

$J_{\ell}=\text{Tr}\left[\left(H-\mu_{\ell}\hat{N}\right)\mathcal{L}_{\ell}\hat{\rho}\right],$
(77)

resulting in

$J_{\ell}=\sum_{s=\pm}\frac{(n_{\ell}^{s}-n_{\bar{\ell}}^{s})\gamma_{\ell}^{s}\gamma_{\bar{\ell}}^{s}}{\gamma_{\ell}^{s}+\gamma_{\bar{\ell}}^{s}}(\epsilon_{s}-\mu_{\ell}),$
(78)

where $\bar{\ell}\neq\ell$.

####  A.2.3 Current fluctuations

Similarly to the local master equation, we find a $\chi$-dependent
Liouvillian:

$\mathcal{L}(\chi)\hat{\rho}(\chi,t)=\mathcal{L}\hat{\rho}(\chi,t)+\sum_{s=\pm}(e^{i\chi}-1)\gamma_{L}^{s}n_{L}^{s}\hat{c}_{s}^{\dagger}\hat{\rho}(\chi,t)\hat{c}_{s}+\sum_{s=\pm}(e^{-i\chi}-1)\gamma_{L}^{s}(1-n_{L}^{s})\hat{c}_{s}\hat{\rho}(\chi,t)\hat{c}_{s}^{\dagger}.$
(79)

In the global master equation, the natural choice of the basis is given by
$\\{p_{0},p_{+},p_{-},p_{D}\\}$, where $p_{\pm}=\langle
0,0|\hat{c}_{\pm}\hat{\rho}\hat{c}_{\pm}^{\dagger}|0,0\rangle$. In this basis
the Liouvillian reads

$\mathbb{L}(\chi)=\footnotesize\begin{pmatrix}-\sum_{s,\ell}n_{l}^{s}\gamma_{\ell}^{s}&e^{-i\chi}(1-n_{L}^{+})\gamma_{L}^{+}+(1-n_{R}^{+})\gamma_{R}^{+}&e^{-i\chi}(1-n_{L}^{-})\gamma_{L}^{-}+(1-n_{R}^{-})\gamma_{R}^{-}&0\\\
e^{i\chi}n_{L}^{+}\gamma_{L}^{+}+n_{R}^{+}\gamma_{R}^{+}&-\sum_{\ell}\left[(1-n_{\ell}^{+})\gamma_{\ell}^{+}+n_{\ell}^{-}\gamma_{\ell}^{-}\right]&0&e^{-i\chi}(1-n_{\ell}^{-})\gamma_{L}^{-}+(1-n_{R}^{-})\gamma_{R}^{-}\\\
e^{i\chi}n_{L}^{-}\gamma_{L}^{-}+n_{R}^{-}\gamma_{R}^{-}&0&-\sum_{\ell}\left[(1-n_{\ell}^{-})\gamma_{\ell}^{-}+n_{\ell}^{+}\gamma_{\ell}^{+}\right]&e^{-i\chi}(1-n_{L}^{+})\gamma_{L}^{+}+(1-n_{R}^{+})\gamma_{R}^{+}\\\
0&e^{i\chi}n_{L}^{-}\gamma_{L}^{-}+n_{R}^{-}\gamma_{R}^{-}&e^{i\chi}n_{L}^{+}\gamma_{L}^{+}+n_{R}^{+}\gamma_{R}^{+}&-\sum_{s,\ell}(1-n_{\ell}^{s})\gamma_{\ell}^{s}\\\
\end{pmatrix}\normalsize.$ (80)

Using the method of full counting statistics outlined in Appendix B, we find

$\langle\\!\langle
I^{2}\rangle\\!\rangle=\sum_{s=\pm}\frac{\gamma_{L}^{s}\gamma_{R}^{s}}{\gamma_{L}^{s}+\gamma_{L}^{s}}\left(n_{L}^{s}+n_{R}^{s}-2n_{L}^{s}n_{R}^{s}-\frac{2\gamma_{L}^{s}\gamma_{R}^{s}\left(n_{L}^{s}-n_{R}^{s}\right)^{2}}{(\gamma_{L}^{s}+\gamma_{L}^{s})^{2}}\right).$
(81)

###  A.3 Nonequilibrium Green’s functions

The key building blocks in the steady-state description with NEGFs are the
retarded ${\bm{G}}^{r}(\omega)$ and the lesser ${\bm{G}}^{<}(\omega)$ Green’s
functions [95, 96]:

$\displaystyle{\bm{G}}^{r}(\omega)=1/\left(\omega{\bm{1}}-\mathbf{H}-{\bm{\Sigma}}^{r}(\omega)\right),$
(82)
$\displaystyle{\bm{G}}^{<}(\omega)={\bm{G}}^{r}(\omega){\bm{\Sigma}}^{<}(\omega){\bm{G}}^{a}(\omega).$
(83)

We use boldface symbols to represent $2\times 2$ matrices in the site indices
(L, R) of the QD system. ${\bm{1}}$ is the identity matrix, the Hamiltonian of
the system [cf. Eq. (46)] is represented by

$\mathbf{H}=\begin{pmatrix}\epsilon_{L}&g\\\ g&\epsilon_{R}\\\ \end{pmatrix},$
(84)

and ${\bm{G}}^{a}=({\bm{G}}^{r})^{\dagger}$. In the absence of Coulomb
interaction the retarded or lesser self-energy ${\bm{\Sigma}}^{r/<}(\omega)$
describes the coupling between the double quantum dot and the reservoirs. In
the wideband limit we have [46, 96, 97]

${\bm{\Sigma}}^{<}(\omega)=i\begin{pmatrix}n_{L}(\omega)\gamma_{L}&0\\\
0&n_{R}(\omega)\gamma_{R}\\\ \end{pmatrix},$ (85)

and

${\bm{\Sigma}}^{r}(\omega)=-\frac{i}{2}\begin{pmatrix}\gamma_{L}&0\\\
0&\gamma_{R}\\\ \end{pmatrix}.$ (86)

From Eqs. (82) and (83) we find the following four elements of
${\bm{G}}^{<}(\omega)$, given by

$\begin{split}&G_{LL}^{<}(\omega)=i\frac{|\omega-\epsilon_{R}+\frac{i\gamma_{R}}{2}^{2}\gamma_{L}n_{L}(\omega)+g^{2}\gamma_{R}n_{R}(\omega)}{|(\omega-\epsilon_{R}+\frac{i\gamma_{R}}{2})(\omega-\epsilon_{L}+\frac{i\gamma_{L}}{2})-g^{2}|^{2}},\\\
&G_{RR}^{<}(\omega)=i\frac{|\omega-\epsilon_{L}+\frac{i\gamma_{L}}{2}^{2}\gamma_{R}n_{R}(\omega)+g^{2}\gamma_{L}n_{L}(\omega)}{|(\omega-\epsilon_{R}+\frac{i\gamma_{R}}{2})(\omega-\epsilon_{L}+\frac{i\gamma_{L}}{2})-g^{2}|^{2}},\\\
&G_{LR}^{<}(\omega)=ig\frac{(\omega-\epsilon_{R}+\frac{i\gamma_{R}}{2})\gamma_{L}n_{L}(\omega)+(\omega-\epsilon_{L}-\frac{i\gamma_{L}}{2})\gamma_{R}n_{R}(\omega)}{|(\omega-\epsilon_{L}+\frac{i\gamma_{L}}{2})(\omega-\epsilon_{R}+\frac{i\gamma_{R}}{2}))-g^{2}|^{2}},\\\
&G_{RL}^{<}(\omega)=ig\frac{(\omega-\epsilon_{R}-\frac{i\gamma_{R}}{2})\gamma_{L}n_{L}(\omega)+(\omega-\epsilon_{L}+\frac{i\gamma_{L}}{2})\gamma_{R}n_{R}(\omega)}{|(\omega-\epsilon_{L}+\frac{i\gamma_{L}}{2})(\omega-\epsilon_{R}+\frac{i\gamma_{R}}{2}))-g^{2}|^{2}}.\\\
\end{split}$ (87)

The expectation values
$\langle\hat{c}_{\ell}^{\dagger}\hat{c}_{\ell^{{}^{\prime}}}\rangle$
($\ell,\ell^{{}^{\prime}}=L,R$) are obtained by taking the inverse Fourier
transform of the corresponding lesser Green’s functions [95, 96]

$\langle\hat{c}_{\ell}^{\dagger}\hat{c}_{\ell^{{}^{\prime}}}\rangle=-i\mathcal{F}\left[G_{\ell^{{}^{\prime}}\ell}^{<}(\omega)\right],$
(88)

where we use the convention

$\mathcal{F}\left[G(\omega)\right]=\lim_{\tau\to
0}\frac{1}{2\pi}\int_{-\infty}^{\infty}d\omega e^{-i\omega\tau}G(\omega).$
(89)

Cumulants of the electron current through the system can be obtained from the
transmission function $\mathcal{T}(\omega)$. The expression for the average
current reads [47, 48, 49]

$\langle
I\rangle=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi}\mathcal{T}(\omega)\left[n_{L}(\omega)-n_{R}(\omega)\right],$
(90)

and the fluctuations are given by

$\begin{split}\langle\\!\langle
I^{2}&\rangle\\!\rangle=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi}\mathcal{T}(\omega)\\{n_{L}(\omega)+n_{R}(\omega)\\\
&-2n_{L}(\omega)n_{R}(\omega)-\mathcal{T}(\omega)[n_{L}(\omega)-n_{R}(\omega)]^{2}\\}.\end{split}$
(91)

The transmission function itself can be expressed with the Green’s functions
[96]

$\mathcal{T}(\omega)=\text{Tr}\left[{\bm{G}}^{a}(\omega){\bm{\Gamma}}_{L}{\bm{G}}^{r}(\omega){\bm{\Gamma}}_{R}\right],$
(92)

where

${\bm{\Gamma}}_{L}=\text{Diag}\left(\gamma_{L},0\right)\hskip
14.22636pt\text{and}\hskip
14.22636pt{\bm{\Gamma}}_{R}=\text{Diag}\left(0,\gamma_{R}\right).$ (93)

For the double quantum dot system it is given by [45, 17]

$\mathcal{T}(\omega)=\frac{\gamma_{L}\gamma_{R}g^{2}}{|\left(\omega-\epsilon_{L}+i\frac{\gamma_{L}}{2}\right)\left(\omega-\epsilon_{R}+i\frac{\gamma_{R}}{2}\right)-g^{2}|^{2}}.$
(94)

In the same framework, the heat current from reservoir $\ell$ reads

$J_{\ell}=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi}(\omega-\mu_{\ell})\mathcal{T}(\omega)\left[n_{\ell}(\omega)-n_{\bar{\ell}}(\omega)\right],$
(95)

where $\bar{\ell}\neq\ell$.

###  A.4 Stochastic model of the system

To motivate the classical model, given by the rate equation in Eq. (16), we
consider the local master equation

$\frac{d\hat{\rho}(t)}{dt}=\mathcal{L}\hat{\rho}(t),$ (96)

with $\mathcal{L}$ given in Eqs. (3) and (4).

We now introduce Nakajima-Zwanzig superoperators $\mathcal{P}$ and
$\mathcal{Q}=\mathcal{I}-\mathcal{P}$, such that $\mathcal{P}$ projects onto
the population subspace of $\hat{\rho}$. Here $\mathcal{I}$ is the identity
superoperator. We note that $\mathcal{P}$ and $\mathcal{Q}$ are projectors,
i.e., $\mathcal{P}^{2}=\mathcal{P}$, $\mathcal{Q}^{2}=\mathcal{Q}$, and
$\mathcal{Q}\mathcal{P}=0$. We also decompose the Liouvillian of the system
into $\mathcal{L}=\mathcal{L}_{0}+\mathcal{V}$, where $\mathcal{V}$ is the
contribution due to the inter-dot tunneling, i.e.,
$\mathcal{V}\hat{\rho}=-ig[\hat{c}_{L}^{\dagger}\hat{c}_{R}+\hat{c}_{R}^{\dagger}\hat{c}_{L},\hat{\rho}]$.
For small $g$ we can treat $\mathcal{V}$ as a perturbation. In the following
derivation, we will use the relations

$\mathcal{P}\mathcal{V}\mathcal{P}=[\mathcal{L}_{0},\mathcal{P}]=[\mathcal{L}_{0},\mathcal{Q}]=0.$
(97)

By acting with $\mathcal{P}$ and $\mathcal{Q}$ on Eq. (96) we obtain the
equations

$\frac{d\mathcal{P}\hat{\rho}(t)}{dt}=\mathcal{P}\mathcal{L}\mathcal{P}\hat{\rho}(t)+\mathcal{P}\mathcal{L}\mathcal{Q}\hat{\rho}(t)$
(98)

and

$\frac{d\mathcal{Q}\hat{\rho}(t)}{dt}=\mathcal{Q}\mathcal{L}\mathcal{P}\hat{\rho}(t)+\mathcal{Q}\mathcal{L}\mathcal{Q}\hat{\rho}(t).$
(99)

The formal solution to the second equation reads

$\mathcal{Q}\hat{\rho}(t)=\mathcal{G}(t,0)\mathcal{Q}\hat{\rho}(0)+\int_{0}^{t}ds\mathcal{G}(t,s)\mathcal{Q}\mathcal{L}\mathcal{P}\hat{\rho}(s),$
(100)

where we have introduced the propagator

$\mathcal{G}(t,s)=e^{\int_{s}^{t}d\tau\mathcal{Q}\mathcal{L}}.$ (101)

After inserting Eq. (100) into Eq. (98) and assuming that
$\mathcal{Q}\hat{\rho}(0)=0$, we find

$\frac{d\mathcal{P}\hat{\rho}(t)}{dt}=\mathcal{P}\mathcal{L}\mathcal{P}\hat{\rho}(t)+\mathcal{P}\mathcal{L}\int_{0}^{t}ds\mathcal{G}(t,s)\mathcal{Q}\mathcal{L}\mathcal{P}\hat{\rho}(s).$
(102)

By inserting $\mathcal{L}=\mathcal{L}_{0}+\mathcal{V}$, applying Eq. (97), and
using the properties of Nakajima-Zwanzig projectors, we obtain the formal
solution

$\begin{split}&\frac{d\mathcal{P}\hat{\rho}(t)}{dt}=\mathcal{L}_{0}\mathcal{P}\hat{\rho}(t)\\\
&+(\mathcal{P}\mathcal{V}\mathcal{Q})\int_{0}^{t}dse^{(\mathcal{Q}\mathcal{L}_{0}\mathcal{Q}+\mathcal{Q}\mathcal{V}\mathcal{Q})s}(\mathcal{Q}\mathcal{V}\mathcal{P})\mathcal{P}\hat{\rho}(t-s),\end{split}$
(103)

where we have also made a substitution $s\to t-s$. The eigenvalues of
$\mathcal{V}$ and $\mathcal{L}_{0}$ are proportional to $g$ and
$\gamma_{\ell}$, respectively. For $g\ll\gamma_{\ell}$, we can treat
$\mathcal{V}$ as a perturbation and neglect its contribution to the
exponential
$e^{(\mathcal{Q}\mathcal{L}_{0}\mathcal{Q}+\mathcal{Q}\mathcal{V}\mathcal{Q})s}$.
This exponential thus decays on the timescale $1/\gamma_{\ell}$. We now assume
that the diagonal elements of the density matrix remain approximately constant
during this timescale, an assumption which is not justified (see below for
more information). This allows us to substitute $\mathcal{P}\hat{\rho}(t)$ for
$\mathcal{P}\hat{\rho}(t-s)$ in Eq. (103) and to extend the integral to
infinity. Equation (103) then reduces to a Markovian master equation for
$\mathcal{P}\hat{\rho}$,

$\frac{d\mathcal{P}\hat{\rho}(t)}{dt}=(\mathcal{L}_{0}-\mathcal{P}\mathcal{V}\mathcal{Q}\mathcal{L}_{0}^{-1}\mathcal{Q}\mathcal{V}\mathcal{P})\mathcal{P}\hat{\rho}(t),$
(104)

where

$\mathcal{L}_{0}^{-1}=-\int_{0}^{\infty}dse^{\mathcal{L}_{0}s}\mathcal{Q}$
(105)

is the Drazin inverse of $\mathcal{L}_{0}$.

Equation (104) can be cast into the form given in Eq. (16) with the stochastic
matrix

$W=\footnotesize\begin{pmatrix}-n_{L}\gamma_{L}-n_{R}\gamma_{R}&(1-n_{L})\gamma_{L}&(1-n_{R})\gamma_{R}&0\\\
n_{L}\gamma_{L}&-(1-n_{L})\gamma_{L}-n_{R}\gamma_{R}-\frac{4g^{2}}{\gamma_{L}+\gamma_{R}}&\frac{4g^{2}}{\gamma_{L}+\gamma_{R}}&(1-n_{R})\gamma_{R}\\\
n_{R}\gamma_{R}&\frac{4g^{2}}{\gamma_{L}+\gamma_{R}}&-n_{L}\gamma_{L}-(1-n_{R})\gamma_{R}-\frac{4g^{2}}{\gamma_{L}+\gamma_{R}}\gamma_{R}&(1-n_{L})\gamma_{L}\\\
0&n_{R}\gamma_{R}&n_{L}\gamma_{L}&-(1-n_{L})\gamma_{L}-(1-n_{R})\gamma_{R}\end{pmatrix}\normalsize.$
(106)

All entries in $W$ are identical to the ones in the local master equation [cf.
(63)] except the element $W_{LR}=W_{RL}$, which comes from the term
$\mathcal{P}\mathcal{V}\mathcal{Q}\mathcal{L}_{0}^{-1}\mathcal{Q}\mathcal{V}\mathcal{P}$
(here we assumed $\delta=0$ for simplicity).

We stress that the approximation giving rise to Eq. (104) is not justified. It
assumes that the diagonal elements remain constant over the timescale
$1/\gamma_{\ell}$. This is, however, exactly the timescale over which the
diagonal elements do change. For our purposes, this is not problematic since
our aim is to derive an analogous classical description for the system. We do
not claim that this classical analog is exact in any limit.

Interestingly, we find the steady-state populations in the classical model to
be the same as in the local master equation given in Eq. (52) for all choices
of parameters. The same holds for the average current, which agrees with
expression (31) obtained with the local master equation. The fact that the
derived rate equation reproduces the average behavior of the local master
equation further justifies its use as an analogous classical model.

To obtain the current fluctuations, counting fields may be introduced in
analogy to Eq. (63). Using full counting statistics, we then obtain

$\begin{split}&\langle\\!\langle
I^{2}\rangle\\!\rangle=\frac{4g^{2}\gamma_{L}\gamma_{R}(n_{L}(1-n_{R})+n_{R}(1-n_{L}))}{(\gamma_{L}+\gamma_{R})(4g^{2}+\gamma_{L}\gamma_{R})}\\\
&-\frac{32[g^{2}\gamma_{L}\gamma_{R}(n_{L}-n_{R})]^{2}[4g^{2}-\gamma_{L}\gamma_{R}+(\gamma_{L}+\gamma_{R})^{2}]}{[(\gamma_{L}+\gamma_{R})(4g^{2}+\gamma_{L}\gamma_{R})]^{3}}.\end{split}$
(107)

##  Appendix B Full counting statistics

In order to obtain the expressions for the current fluctuations, we employ the
technique of full counting statistics. In this framework the density matrix of
the system is resolved into $\hat{\rho}(n,t)$, the (unnormalized) density
matrices given that a net number of $n$ electrons tunneled from the left
reservoir into the DQD during the time interval $[0,t]$ [98, 99, 100]. The
probability that $n$ electrons have tunneled can then be obtained by

$P[n]=\text{Tr}[\hat{\rho}(n,t)].$ (108)

The original density matrix can be recovered by taking the sum

$\hat{\rho}(t)=\sum_{n}\hat{\rho}(n,t).$ (109)

The cumulants of the net number of transported electrons can be obtained from
the cumulant generating function $S(\chi,t)$, given by

$S(\chi,t)=\log{\sum_{n}P[n]e^{i\chi
n}}=\log{\text{Tr}\left[\hat{\rho}(\chi,t)\right]},$ (110)

where

$\hat{\rho}(\chi,t)=\sum_{n}\hat{\rho}(n,t)e^{i\chi n},$ (111)

and $\chi$ is known as the counting field.

The time evolution of $\hat{\rho}(\chi,t)$ can be found by inserting the
$\chi$-dependent term, which tracks electron jumps, into the original Lindblad
master equation [100]. This results in a counting-field-dependent master
equation,

$\frac{d\hat{\rho}(\chi,t)}{dt}=\mathcal{L}(\chi)\hat{\rho}(\chi,t).$ (112)

The specific form of $\mathcal{L}(\chi)$ depends on the choice of the master
equation, and the original Liouvillian can be recovered with
$\mathcal{L}=\mathcal{L}(0)$. The formal solution to Eq. (112) is given by

$\hat{\rho}(\chi,t)=e^{\mathcal{L}(\chi)t}\hat{\rho}(0,0).$ (113)

The cumulants of the electron current in the steady state can be obtained from
the scaled cumulant generating function of $n$ in the long-time limit,

$\langle\\!\langle
I^{k}\rangle\\!\rangle=\frac{d^{k}}{d(i\chi)^{k}}\left(\lim_{t\rightarrow\infty}\frac{S(\chi,t)}{t}\right)|_{\chi=0}=\frac{d^{k}}{d(i\chi)^{k}}\lambda(\chi)|_{\chi=0},$
(114)

where $\lambda(\chi)$ is the eigenvalue of $\mathcal{L}(\chi)$ with the
largest real part. The average of the current, $\langle
I\rangle=\langle\\!\langle I\rangle\\!\rangle$, and its fluctuations,
$\langle\\!\langle I^{2}\rangle\\!\rangle$, can be evaluated by taking
$k=1,2$, respectively.

Since a direct of computation of $\lambda(\chi)$ is not straightforward for
our system, we make use of the scheme described in [101]. We express the
Liouvillian $\mathcal{L}(\chi)$ as a matrix $\mathbb{L}(\chi)$ in a basis that
can be freely chosen and may be different depending on whether we use the
local or global master equation. Let $P_{\chi}(x)$ be the characteristic
polynomial

$P_{\chi}(x)=\det{\left[\mathbb{L}(\chi)-x\mathbb{I}\right]},$ (115)

where $\mathbb{I}$ is the identity matrix. Upon expanding in $x$ we obtain

$P_{\chi}(x)=\sum_{j=0}^{M}a_{j}x^{j}=\sum_{j=0}^{M}\sum_{k=0}^{\infty}\frac{a_{j}^{(k)}(i\chi)^{k}}{k!}x^{j},$
(116)

where

$a_{j}=\frac{\partial^{j}}{\partial x^{j}}\frac{P_{\chi}(x)}{j!}|_{x=0}\hskip
11.38092pt\text{and}\hskip
11.38092pta_{j}^{(k)}=\frac{d^{k}}{d(i\chi)^{k}}a_{j}|_{\chi=0}.$ (117)

Here $M$ is the dimension of the matrix $\mathbb{L}(\chi)$. By definition, the
eigenvalue $\lambda(\chi)$ of $\mathbb{L}(\chi)$ (with the largest real part)
fulfills

$P_{\chi}(\lambda(\chi))=0.$ (118)

It follows that

$0=\frac{d}{d(i\chi)}P_{\chi}(\lambda(\chi))|_{\chi=0}=a_{0}^{(1)}+a_{1}^{(0)}\langle
I\rangle.$ (119)

In order to obtain the right-hand side, we first expanded
$P_{\chi}(\lambda(\chi))$ using Eq. (116), used $\lambda(0)=0$, as well as Eq.
(114). We thus find

$\langle I\rangle=-\frac{a_{0}^{(1)}}{a_{1}^{(0)}}.$ (120)

The second cumulant of the current may be obtained in a similar fashion from

$0=\frac{d^{2}}{d(i\chi)^{2}}P_{\chi}(\lambda(\chi))|_{\chi=0},$ (121)

resulting in

$\langle\\!\langle
I^{2}\rangle\\!\rangle=-\frac{a_{0}^{(2)}+2a_{1}^{(1)}\langle
I\rangle+2a_{2}^{(0)}\langle I\rangle^{2}}{a_{1}^{(0)}}.$ (122)

##  Appendix C Interacting electrons

When investigating the system interacting electrons we restrict ourselves to
the case with equal onsite energies of quantum dots, where the Hamiltonian of
the system is given by Eq. (1). In the presence of Coulomb repulsion, the
coefficients of the density matrix in the basis already introduced in Eq. (47)
are given by

$\begin{split}&p_{0}=1-\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle-\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle+\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle,\\\
&p_{L}=\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\rangle-\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle,\\\
&p_{R}=\langle\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle-\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle,\\\
&p_{D}=\langle\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R}\rangle,\\\
&\alpha=\langle\hat{c}_{R}^{\dagger}\hat{c}_{L}\rangle.\\\ \end{split}$ (123)

###  C.1 Semilocal master equation

In the semilocal master equation, the Lindblad superoperator for bath $\ell$
is given by [36]

$\mathcal{L}_{\ell}=\gamma_{\ell}\left[n_{\ell}\mathcal{D}[(1-\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}})\hat{c}_{\ell}^{\dagger}]+n_{\ell}^{U}\mathcal{D}[\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}}\hat{c}_{\ell}^{\dagger}]+(1-n_{\ell})\mathcal{D}[(1-\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}})\hat{c}_{\ell}]+(1-n_{\ell}^{U})\mathcal{D}[\hat{c}_{\bar{\ell}}^{\dagger}\hat{c}_{\bar{\ell}}\hat{c}_{\ell}]\right],$
(124)

where

$n_{\ell}:=n_{\ell}(\epsilon),$ (125)

and

$n_{\ell}^{U}:=n_{\ell}(\epsilon+U).$ (126)

The steady state of the system is given by the eigenvector corresponding to
the zero eigenvalue of the Liouvillian. TheLiouvillian can be computed from
Eq. (134) evaluated at zero counting field $\chi=0$. We find

$\begin{split}&p_{0}=\frac{4g^{2}(1-\bar{n})(1-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})\left((1-n_{L})(1-n_{R})-n_{L}^{U}n_{R}^{U}+\sum_{l}\frac{\gamma_{l}n_{l}^{U}(n_{\bar{l}}-n_{l}\Delta_{\bar{l}})}{\gamma_{L}+\gamma_{R}}\right)}{4g^{2}(1+\bar{n}-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\Delta_{L}\Delta_{R})(1-\bar{n}+\bar{n}_{U})}\\\
&p_{L}=\frac{4g^{2}\bar{n}(1-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})\left(n_{L}(1-n_{R})+\sum_{l}\frac{\gamma_{l}n_{l}n_{l}^{U}\Delta_{\bar{l}}}{\gamma_{L}+\gamma_{R}}+\frac{\gamma_{R}(n_{L}^{U}n_{R}-n_{L}n_{R}^{U})}{\gamma_{L}+\gamma_{R}}\right)}{4g^{2}(1+\bar{n}-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\Delta_{L}\Delta_{R})(1-\bar{n}+\bar{n}_{U})}\\\
&p_{R}=\frac{4g^{2}\bar{n}(1-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})\left((1-n_{L})n_{R}+\sum_{l}\frac{\gamma_{l}n_{l}n_{l}^{U}\Delta_{\bar{l}}}{\gamma_{L}+\gamma_{R}}+\frac{\gamma_{L}(n_{R}^{U}n_{L}-n_{R}n_{L}^{U})}{\gamma_{L}+\gamma_{R}}\right)}{4g^{2}(1+\bar{n}-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\Delta_{L}\Delta_{R})(1-\bar{n}+\bar{n}_{U})}\\\
&p_{D}=\frac{4g^{2}\bar{n}\bar{n}_{U}+\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})\sum_{l}\frac{\gamma_{l}n_{l}^{U}(n_{\bar{l}}-n_{l}\Delta_{\bar{l}})}{\gamma_{L}+\gamma_{R}}}{4g^{2}(1+\bar{n}-\bar{n}_{U})+\gamma_{L}\gamma_{R}(1-\Delta_{L}\Delta_{R})(1-\bar{n}+\bar{n}_{U})}\\\
&\alpha=\frac{2ig\gamma_{L}\gamma_{R}\left[n_{L}(1-n_{R}^{U})-n_{R}(1-n_{L}^{U})\right]}{(\gamma_{L}+\gamma_{R})\left[4g^{2}(1+\bar{n}-\bar{n}_{U})+(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})\right]},\\\
\end{split}$ (127)

where

$\bar{n}_{U}=\frac{\gamma_{L}n_{L}^{U}+\gamma_{R}n_{R}^{U}}{\gamma_{L}+\gamma_{R}},$
(128)

and $\Delta_{l}=n_{l}-n_{l}^{U}$. Expressions for coherence $|\alpha|$ as well
as the concurrence $\mathcal{C}$, which is defined in Eq. (20), can be
obtained from the coefficients shown above.

####  C.1.1 Average electron current

For the average electron current through the system, given by Eq. (11), we
find

$\langle
I\rangle=\frac{4g^{2}\gamma_{L}\gamma_{R}\left[n_{L}(1-n_{R}^{U})-n_{R}(1-n_{L}^{U})\right]}{(\gamma_{L}+\gamma_{R})\left[4g^{2}(1+\bar{n}-\bar{n}_{U})+(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})\right]}.$
(129)

####  C.1.2 Heat currents

The heat current from bath $\ell$ is defined in complete analogy to the local
master equation, i.e., Eq. (58) [36]. In the presence of Coulomb interactions,
the bookkeeping Hamiltonian reads

$\hat{H}_{TD}=\epsilon(\hat{c}_{L}^{\dagger}\hat{c}_{L}+\hat{c}_{R}^{\dagger}\hat{c}_{R})+U\hat{c}_{L}^{\dagger}\hat{c}_{L}\hat{c}_{R}^{\dagger}\hat{c}_{R},$
(130)

resulting in

$\begin{split}\langle Q_{L}\rangle&=(\epsilon-\mu_{L})\langle I\rangle\\\
&+\frac{U\gamma_{L}\gamma_{R}}{\gamma_{L}+\gamma_{R}}\frac{4g^{2}(n_{L}^{U}-n_{R}^{U})\bar{n}+\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})(n_{L}^{U}n_{R}-n_{R}^{U}n_{L}-n_{L}n_{L}^{U}\Delta_{R}+n_{R}n_{R}^{U}\Delta_{L})}{4g^{2}(1+\bar{n}-\bar{n}_{U})+(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})},\end{split}$
(131)

and

$\begin{split}\langle Q_{R}\rangle&=-(\epsilon-\mu_{R})\langle I\rangle\\\
&-\frac{U\gamma_{L}\gamma_{R}}{\gamma_{L}+\gamma_{R}}\frac{4g^{2}(n_{L}^{U}-n_{R}^{U})\bar{n}+\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})(n_{L}^{U}n_{R}-n_{R}^{U}n_{L}-n_{L}n_{L}^{U}\Delta_{R}+n_{R}n_{R}^{U}\Delta_{L})}{4g^{2}(1+\bar{n}-\bar{n}_{U})+(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})},\end{split}$
(132)

respectively.

####  C.1.3 Current fluctuations

The current fluctuations are obtained in complete analogy to the method used
in the noninteracting case, i.e., full counting statistics (cf. Appendix B).
The $\chi$-dependent Liouvillian reads

$\begin{split}\frac{d\hat{\rho}}{dt}&=\mathcal{L}\hat{\rho}+\gamma_{L}(e^{i\chi}-1)\left[n_{L}(1-c_{R}^{\dagger}c_{R})\hat{c}_{L}^{\dagger}\hat{\rho}\hat{c}_{L}(1-c_{R}^{\dagger}c_{R})+n_{L}^{U}(\hat{c}_{R}^{\dagger}\hat{c}_{R})\hat{c}_{L}^{\dagger}\hat{\rho}\hat{c}_{L}(\hat{c}_{R}^{\dagger}\hat{c}_{R})\right]\\\
&+\gamma_{L}(e^{-i\chi}-1)\left[n_{L}(1-c_{R}^{\dagger}c_{R})\hat{c}_{L}\hat{\rho}\hat{c}_{L}^{\dagger}(1-c_{R}^{\dagger}c_{R})+n_{L}^{U}(\hat{c}_{R}^{\dagger}\hat{c}_{R})\hat{c}_{L}\hat{\rho}\hat{c}_{L}^{\dagger}(\hat{c}_{R}^{\dagger}\hat{c}_{R})\right]\end{split}$
(133)

In the basis
$\\{p_{0},p_{L},p_{R},p_{D},\text{Re}[\alpha],\text{Im}[\alpha]\\}$ it is
given by

$\mathbb{L}(\chi)=\begin{pmatrix}-n_{L}\gamma_{L}-n_{R}\gamma_{R}&e^{-i\chi}(1-n_{L})n_{L}&(1-n_{R})\gamma_{R}&0&0&0\\\
e^{i\chi}n_{L}\gamma_{L}&-(1-n_{L})\gamma_{L}-n_{R}^{U}\gamma_{R}&0&(1-n_{R}^{U})\gamma_{R}&0&-2g\\\
n_{R}\gamma_{R}&0&-n_{L}^{U}\gamma_{L}-(1-n_{R})\gamma_{R}&e^{-i\chi}(1-n_{L}^{U})\gamma_{L}&0&2g\\\
0&n_{R}\gamma_{R}&e^{i\chi}n_{L}\gamma_{L}&-(1-n_{L})\gamma_{L}-(1-n_{R})\gamma_{R}&0&0\\\
0&0&0&0&\gamma_{U}&2\delta\\\ 0&g&-g&0&-2\delta&\gamma_{U}\end{pmatrix},$
(134)

where $\gamma_{U}=-\sum_{l}\frac{1-n_{l}+n_{l}^{U}}{2}\gamma_{l}$. We thereby
find the current fluctuations

$\begin{split}\langle\langle
I^{2}\rangle\rangle&=\frac{n_{L}(1-n_{R})+n_{R}(1-n_{L})+2n_{L}n_{R}\bar{n}_{U}-2n_{L}^{U}n_{R}^{U}\bar{n}+\frac{\gamma_{L}(n_{L}n_{R}^{U}-n_{R}n_{L}^{U})+\gamma_{R}(n_{R}n_{L}^{U}-n_{L}n_{R}^{U})}{\gamma_{L}+\gamma_{R}}}{(n_{L}(1-n_{R}^{U})-n_{R}(1-n_{L}^{U}))}\langle
I\rangle\\\ &+\frac{2\langle
I\rangle^{2}}{\gamma_{L}+\gamma_{R}}\left(\frac{n_{L}^{U}(1-n_{R})-n_{R}^{U}(1-n_{L})}{n_{L}(1-n_{R}^{U})-n_{R}(1-n_{L}^{U})}-\frac{2(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}+8g^{2}+(\gamma_{L}+\gamma_{R})^{2}(1-\bar{n}+\bar{n}_{U})}{4g^{2}(1+\bar{n}-\bar{n}_{U})+(1-\Delta_{L}\Delta_{R})\gamma_{L}\gamma_{R}(1-\bar{n}+\bar{n}_{U})}\right).\end{split}$
(135)

###  C.2 Nonequilibrium Green’s functions

To determine the densities of the system via NEGFs, we can still proceed as in
Appendix A.3. However, we have to account, within NEGF, for the Coulomb
interaction via a many-body contribution to the self-energy [80, 81, 82]. As
before, the key building blocks in the steady-state description with NEGF in
the presence of interactions are the retarded ${\bm{G}}^{r}(\omega)$ and
lesser ${\bm{G}}^{<}(\omega)$ Green’s functions, defined in Eqs. (82) and (83)
and for convenience repeated here:

$\displaystyle{\bm{G}}^{r}(\omega)=1/\left(~{}\omega{\bm{1}}-{\bm{H}}-{\bm{\Sigma}}^{r}(\omega)\right),$
(136)
$\displaystyle{\bm{G}}^{<}(\omega)={\bm{G}}^{r}(\omega){\bm{\Sigma}}^{<}(\omega){\bm{G}}^{a}(\omega).$
(137)

The self-energy ${\bm{\Sigma}}$ now has two parts, coming from the embedding
(emb) and many-body (MB) contributions:
${\bm{\Sigma}}^{r/<}={\bm{\Sigma}}^{r/<}_{\text{emb}}+{\bm{\Sigma}}^{r/<}_{\text{MB}}$,
where the embedding term has been discussed before in Appendix. A.3. The
correlation effects due to interactions are introduced via the many-body self-
energy
${\bm{\Sigma}}_{\text{MB}}^{r/<}={\bm{\Sigma}}_{\text{MB}}^{r/<}[{\bm{G}}^{r},{\bm{G}}^{<}]$.
For ${\bm{\Sigma}}_{\text{MB}}$ we use the second Born approximation [97, 102,
103], keeping all Feynman diagrams up to second order. These equations are
then solved self-consistently, using the Pulay scheme [104] to improve
convergence. Self-consistency ensures that general conservation laws [81, 80,
82] are obeyed (e.g., if $I_{L}$ is the current from lead $L$, then
$I_{R}=-I_{L}$). The embedding self-energy is independent of the self-
consistency cycle and is calculated separately only once. For the particle
density in the device region and the average current $\langle I\rangle$ coming
from the left reservoir, the expressions are [96], respectively,

$\displaystyle\langle\hat{N}_{\ell}\rangle=\langle\hat{c}_{\ell}^{\dagger}\hat{c}_{\ell}\rangle=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi
i}G_{\ell\ell}^{<}(\omega),$ (138) $\displaystyle\langle
I\rangle=\\!\int_{-\infty}^{\infty}\frac{d\omega}{2\pi
i}\text{Tr}\left[{\bm{\Gamma}}_{L}\left({\bm{G}}^{<}(\omega)-2\pi
in_{L}(\omega){\bm{A}}(\omega)\right)\right].$ (139)

Here $2\pi{\bm{A}(\omega)}=i({\bm{G}}^{r}(\omega)-{\bm{G}}^{a}(\omega))$.
Finally, the concurrence $\mathcal{C}$ is computed starting from the
expressions for $\langle\hat{N}_{\ell}\rangle$ ($\ell=L,R$), directly obtained
from Eq. (138), and from $\langle\hat{N}_{L}\hat{N}_{R}\rangle$. The latter
quantity is obtained according to a prescription introduced in Ref. [105].

We start with the defining expression for $\Sigma$. This in standard NEGF
notation on the Keldysh contour reads

$\displaystyle\int d3\ \Sigma(1,3)G(3,2)=-i\int d3\
\nu(1,3)G_{2}(1,3^{+};2,3^{++}),$ (140)

in terms of the two-particle Green’s function $G_{2}$:

$G_{2}(1,2;3,4)=-\langle\Psi|\mathcal{T}\left[\hat{c}(1)\hat{c}(2)\hat{c}^{\dagger}(4)\hat{c}^{\dagger}(3)\right]|\Psi\rangle.$
(141)

Here and in Eq. (140), $1\equiv(\ell_{1},z_{1})$ (and similarly for 2, 3, 4),
with $\ell_{1}$ a site index and $z_{1}$ a time variable on the Keldysh
contour.

We now make the site indices explicit in Eq. (140), i.e.,
$1\rightarrow(\ell,z),\ 2\rightarrow(\ell^{\prime},z^{\prime}),\
3\rightarrow(\bar{\ell},\bar{z})$, and set
$(\ell^{\prime},z^{\prime})=(\ell,z^{+})=(L,z^{+})$. Furthermore, with
nonlocal interactions, the interaction term becomes
$\nu_{\ell\bar{\ell}}(z,\bar{z})=U(\delta_{\ell
L}\delta_{\bar{\ell}R}+\delta_{\ell R}\delta_{\bar{\ell}L})\delta(z-\bar{z})$.
Equation 140 then becomes

$\displaystyle\sum_{\bar{\ell}}\int d\bar{z}\
\Sigma_{L\bar{\ell}}(z,\bar{z})G_{\bar{\ell}L}(\bar{z},z^{+})$
$\displaystyle=-iUG_{LRLR}(z,z^{+};z^{+},z^{++}),$ (142)

where, here and in the following, $\bar{\ell}=L,R$ in the sums. Then, using
the definition of $G_{2}$ in Eq. (141),

$\displaystyle\sum_{\bar{\ell}}\int d\bar{z}\
\Sigma_{L\bar{\ell}}(z,\bar{z})G_{\bar{\ell}L}(\bar{z},z^{+})$
$\displaystyle=iU\langle\Psi|\mathcal{T}\left[\hat{c}_{L}(z)\hat{c}_{R}(z^{+})\hat{c}_{R}^{\dagger}(z^{++})\hat{c}^{\dagger}_{L}(z^{+})\right]|\Psi\rangle,$
(143)

and by applying the Langreth rules, we arrive at an expression for physical
times:

$\displaystyle\frac{-i}{U}\sum_{\bar{\ell}}\int
d\bar{t}\big{[}\Sigma^{<}_{L\bar{\ell}}(t,\bar{t})G^{a}_{\bar{\ell}L}(\bar{t},t^{+})+\Sigma^{r}_{L\bar{\ell}}(t,\bar{t})G^{<}_{\bar{\ell}L}(\bar{t},t^{+})\big{]}$
$\displaystyle=\langle\Psi|\left[\hat{c}_{R}^{\dagger}(t^{++})\hat{c}_{R}(t^{+})\hat{c}^{\dagger}_{L}(t^{+})\hat{c}_{L}(t)\right]|\Psi\rangle=\langle\hat{N}_{R}\hat{N}_{L}\rangle.$
(144)

Finally, moving to frequency space, we have

$\displaystyle\langle\hat{N}_{L}\hat{N}_{R}\rangle=$
$\displaystyle\frac{1}{U}\sum_{\bar{\ell}}\int\frac{d\omega}{2\pi i}\
(\Sigma^{<}_{L\bar{\ell}}(\omega)G^{a}_{\bar{\ell}L}(\omega)$ (145)
$\displaystyle+\Sigma^{r}_{L\bar{\ell}}(\omega)G^{<}_{\bar{\ell}L}(\omega)).$

For computational convenience we separate the self-energy into a correlation
part and one coming from the Hartree-Fock interaction:

$\displaystyle\langle\hat{N}_{R}\hat{N}_{L}\rangle=(\langle\hat{N}_{R}\rangle\langle\hat{N}_{L}\rangle-\langle\hat{c}_{R}\hat{c}_{L}\rangle\langle\hat{c}_{L}\hat{c}_{R}\rangle)\hskip
56.9055pt$ (146)
$\displaystyle+\frac{1}{U}\sum_{\bar{\ell}}\int\frac{d\omega}{2\pi i}\
(\Sigma^{<}_{c,L\bar{\ell}}(\omega)G^{a}_{\bar{\ell}L}(\omega)+\Sigma^{r}_{c,L\bar{\ell}}(\omega)G^{<}_{\bar{\ell}L}(\omega)),$

where the subscript c on the self-energy indicates that it contains only the
contribution from correlations.

The applicability of this expression has been examined in Ref. [105]. There it
was shown that some approximate but conserving many-body self-energies (among
them, those obtained in the second Born approximation) do not always guarantee
a non-negative value for local double occupancies
$\langle\hat{N}_{\ell\uparrow}\hat{N}_{\ell\downarrow}\rangle$ in a spinful
model with Hubbard interactions. Here, by contrast, we are considering a
spinless model, with nonlocal interactions and correlations and in all
applications of Eq. (145) in the present paper, no negative correlations were
encountered.

##  Appendix D Operational nonclassicality

Following Ref. [10], we consider the usefulness of the DQD steady state for
Bell nonlocality and quantum teleportation.

###  D.1 Nonlocality

In the CHSH scenario [79] two players that share a bipartite quantum state,
make measurements $x$ and $y$, and obtain outcomes $a$ and $b$. The
conditional probability $p(a,b|x,y)$ is nonlocal if it does not admit a local
hidden-variable model

$p(a,b|x,y)=\sum_{\lambda}p(\lambda)p(a|x,\lambda)p(b|y,\lambda)$ (147)

for any $\lambda$ with an underlying distribution $p(\lambda)$ [106]. In the
CHSH scenario $a,b,x,y\in\\{0,1\\}$, and $p(a,b|x,y)$ can be explained by a
local hidden-variable model whenever [79]

$\text{CHSH}=\sum_{a,b,x,y}(-1)^{a+b+xy}p(a,b|x,y)\leq 2.$ (148)

For the density matrix of the double quantum dot, given in Eq. (19), the
maximal value for the CHSH quantity is given by [11]

$\text{CHSH}=2\sqrt{8|\alpha|^{2}+(1-\Delta)^{2}-\min\\{4|\alpha|^{2},(1-\Delta)^{2}\\}},$
(149)

where $\Delta=2p_{0}+2p_{D}$.

For noninteracting electrons, we examine the maximum possible CHSH in the
regime $eV\to\infty$, where the corresponding Fermi distributions reduce to
$n_{L}\to 1$ and $n_{R}\to 0$. Using the expressions for the coefficients of
the density matrix, given in Eq. (53), we do not find any parameters of the
system that give rise to a violation of CHSH inequality.

In the presence of Coulomb interactions, we proceed similarly, and consider
the regime described in Sec. V.3. Under this set of conditions, the
probability of the double occupation $p_{D}$ vanishes, and the maximal value
of the CHSH quantity reduces to

$\text{CHSH}=2\sqrt{8|\alpha|^{2}+(1-2p_{0})^{2}-\min\\{4|\alpha|^{2},(1-2p_{0})^{2}\\}},$
(150)

where the corresponding $|\alpha|$ and $p_{0}$ are given by Eqs. (43) and
(127), respectively. We find that the maximal value of the CHSH quantity is
$\text{CHSH}=2\sqrt{3/2}$. It is achieved when $\gamma_{R}/g=2\sqrt{2}$ and
$\gamma_{L}/g\to\infty$.

###  D.2 Quantum teleportation

Quantum teleportation [78] is a protocol where an unknown quantum state is
sent between two parties by the means of a shared entangled state and
classical communication. In the original implementation, a maximally entangled
Bell state $|\psi^{-}\rangle=(|01\rangle-|10\rangle)/\sqrt{2}$ was used,
giving the average fidelity of teleportation $f=1$. However, when the parties
implementing the protocol do not share a Bell state, but any other less
entangled state $\hat{\rho}$, the fidelity of the teleportation is lower with
the threshold value $f=2/3$ for separable states. The average fidelity
$f=(1+2F)/3$ can be computed from the so-called singlet fraction [83]

$F=\max_{U}\langle\psi^{-}|(U\otimes I)~{}\hat{\rho}~{}(U\otimes
I)|\psi^{-}\rangle,$ (151)

where the optimization is performed over the unitary transformations $U$. For
the density matrix of the double quantum dot, given in Eq. (19), we have that
if its singlet fraction exceeds that of any separable state, i.e., $F>1/2$,
then it is given by [11]

$F=|\alpha|+\frac{1-p_{0}-p_{D}}{2}.$ (152)

In the absence of Coulomb interactions $F$ is saturated at the maximal bias,
i.e., when $eV\to\infty$. Using the expressions for the steady state of the
density matrix, given by Eq. (53), we find the corresponding singlet fraction,

$F=\frac{\gamma_{L}\gamma_{R}\left(8g^{2}+4g(\gamma_{L}+\gamma_{R})+(\gamma_{L}+\gamma_{R})^{2}\right)}{2(\gamma_{L}+\gamma_{R})^{2}(4g^{2}+\gamma_{L}\gamma_{R})}.$
(153)

It is maximized for the symmetric coupling strength to the baths that fulfills
$g/\gamma=(\sqrt{5}-1)/4$, resulting in $F=(3+\sqrt{5})/8\approx 0.65$.

Under the conditions described in Sec. V.3 for interacting electrons,
$p_{D}=0$, whereas $|\alpha|$ and $p_{0}$ are given by Eqs. (43) and (127),
respectively. The singlet fraction is then given by

$F=\frac{(8g^{2}+4g\gamma_{R}+\gamma_{R}^{2})\gamma_{L}}{2\gamma_{R}^{2}\gamma_{L}+8g^{2}(\gamma_{R}+2\gamma_{L})}.$
(154)

$F$ is maximized when $\gamma_{R}/g=2\sqrt{2}$ and $\gamma_{L}\to\infty$,
which coincides with the parameters that saturate concurrence, and amounts to
$F=(2+\sqrt{2})/4\approx 0.85$. If the coupling strengths to the baths are
equal, the maximal singlet fraction occurs at $\gamma/g=1+\sqrt{13}$ and is
equal to $F=(5+\sqrt{13})/12\approx 0.72$.

## References

  * Jammer [1966] M. Jammer,  _The Conceptual Development of Quantum Mechanics_ (McGraw-Hill, New York, 1966).
  * Streltsov  _et al._ [2017] A. Streltsov, G. Adesso, and M. B. Plenio, Colloquium: Quantum coherence as a resource, Rev. Mod. Phys. 89, 041003 (2017).
  * Scandi  _et al._ [2020] M. Scandi, H. J. D. Miller, J. Anders, and M. Perarnau-Llobet, Quantum work statistics close to equilibrium, Phys. Rev. Research 2, 023377 (2020).
  * Latune  _et al._ [2020] C. L. Latune, I. Sinayskiy, and F. Petruccione, Negative contributions to entropy production induced by quantum coherences, Phys. Rev. A 102, 042220 (2020).
  * Tajima and Funo [2021] H. Tajima and K. Funo, Superconducting-like heat current: Effective cancellation of current-dissipation trade-off by quantum coherence, Phys. Rev. Lett. 127, 190604 (2021).
  * Benenti  _et al._ [2017] G. Benenti, G. Casati, K. Saito, and R. Whitney, Fundamental aspects of steady-state conversion of heat to work at the nanoscale, Phys. Rep. 694, 1 (2017).
  * Pekola and Karimi [2021] J. P. Pekola and B. Karimi, Colloquium: Quantum heat transport in condensed matter systems, Rev. Mod. Phys. 93, 041001 (2021).
  * Mitchison [2019] M. T. Mitchison, Quantum thermal absorption machines: refrigerators, engines and clocks, Contemp. Phys. 60, 164 (2019).
  * Degen  _et al._ [2017] C. L. Degen, F. Reinhard, and P. Cappellaro, Quantum sensing, Rev. Mod. Phys. 89, 035002 (2017).
  * Brask  _et al._ [2015] J. B. Brask, G. Haack, N. Brunner, and M. Huber, Autonomous quantum thermal machine for generating steady-state entanglement, New J. Phys. 17, 113029 (2015).
  * Bohr Brask  _et al._ [2022] J. Bohr Brask, F. Clivaz, G. Haack, and A. Tavakoli, Operational nonclassicality in minimal autonomous thermal machines, Quantum 6, 672 (2022).
  * Khandelwal  _et al._ [2020] S. Khandelwal, N. Palazzo, N. Brunner, and G. Haack, Critical heat current for operating an entanglement engine, New J. Phys. 22, 073039 (2020).
  * Man  _et al._ [2019] Z.-X. Man, A. Tavakoli, J. B. Brask, and Y.-J. Xia, Improving autonomous thermal entanglement generation using a common reservoir, Phys. Scr. 94, 075101 (2019).
  * Wang  _et al._ [2019] Z. Wang, W. Wu, and J. Wang, Steady-state entanglement and coherence of two coupled qubits in equilibrium and nonequilibrium environments, Phys. Rev. A 99, 042320 (2019).
  * Kießlich  _et al._ [2006] G. Kießlich, P. Samuelsson, A. Wacker, and E. Schöll, Counting statistics and decoherence in coupled quantum dots, Phys. Rev. B 73, 033312 (2006).
  * Ptaszyński [2018] K. Ptaszyński, Coherence-enhanced constancy of a quantum thermoelectric generator, Phys. Rev. B 98, 085425 (2018).
  * Agarwalla and Segal [2018] B. K. Agarwalla and D. Segal, Assessing the validity of the thermodynamic uncertainty relation in quantum systems, Phys. Rev. B 98, 155438 (2018).
  * Kalaee  _et al._ [2021] A. A. S. Kalaee, A. Wacker, and P. P. Potts, Violating the thermodynamic uncertainty relation in the three-level maser, Phys. Rev. E 104, L012103 (2021).
  * Cangemi  _et al._ [2020] L. M. Cangemi, V. Cataudella, G. Benenti, M. Sassetti, and G. De Filippis, Violation of thermodynamics uncertainty relations in a periodically driven work-to-work converter from weak to strong dissipation, Phys. Rev. B 102, 165418 (2020).
  * Rignon-Bret  _et al._ [2021] A. Rignon-Bret, G. Guarnieri, J. Goold, and M. T. Mitchison, Thermodynamics of precision in quantum nanomachines, Phys. Rev. E 103, 012133 (2021).
  * Menczel  _et al._ [2021] P. Menczel, E. Loisa, K. Brandner, and C. Flindt, Thermodynamic uncertainty relations for coherently driven open quantum systems, J. Phys. A: Math. Theor. 54, 314002 (2021).
  * Liu and Segal [2019] J. Liu and D. Segal, Thermodynamic uncertainty relation in quantum thermoelectric junctions, Phys. Rev. E 99, 062141 (2019).
  * Horodecki  _et al._ [2009] R. Horodecki, P. Horodecki, M. Horodecki, and K. Horodecki, Quantum entanglement, Rev. Mod. Phys. 81, 865 (2009).
  * Brunner  _et al._ [2014] N. Brunner, D. Cavalcanti, S. Pironio, V. Scarani, and S. Wehner, Bell nonlocality, Rev. Mod. Phys. 86, 419 (2014).
  * Nielsen and Chuang [2000] M. A. Nielsen and I. L. Chuang,  _Quantum Computation and Quantum Information_ (Cambridge University Press, Cambridge, 2000).
  * Renner and Wolf [0] R. Renner and R. Wolf, Quantum advantage in cryptography, AIAA Journal 0, 1 (0).
  * Giovannetti  _et al._ [2004] V. Giovannetti, S. Lloyd, and L. Maccone, Quantum-enhanced measurements: Beating the standard quantum limit, Science 306, 1330 (2004).
  * Giovannetti  _et al._ [2006] V. Giovannetti, S. Lloyd, and L. Maccone, Quantum metrology, Phys. Rev. Lett. 96, 010401 (2006).
  * Escher  _et al._ [2011] B. M. Escher, R. L. de Matos Filho, and L. Davidovich, General framework for estimating the ultimate precision limit in noisy quantum-enhanced metrology, Nat. Phys. 7, 406 (2011).
  * Barato and Seifert [2015] A. C. Barato and U. Seifert, Thermodynamic uncertainty relation for biomolecular processes, Phys. Rev. Lett. 114, 158101 (2015).
  * Gingrich  _et al._ [2016] T. R. Gingrich, J. M. Horowitz, N. Perunov, and J. L. England, Dissipation bounds all steady-state current fluctuations, Phys. Rev. Lett. 116, 120601 (2016).
  * Horowitz and Gingrich [2019] J. Horowitz and T. Gingrich, Thermodynamic uncertainty relations constrain non-equilibrium fluctuations, Nat. Phys. 16, 1 (2019).
  * Terlizzi and Baiesi [2018] I. D. Terlizzi and M. Baiesi, Kinetic uncertainty relation, J. Phys. A: Math. Theor. 52, 02LT03 (2018).
  * Vo  _et al._ [2022] V. T. Vo, T. V. Vu, and Y. Hasegawa, Unified thermodynamic–kinetic uncertainty relation, J. Phys. A: Math. Theor. 55, 405004 (2022).
  * Hofer  _et al._ [2017] P. P. Hofer, M. Perarnau-Llobet, L. D. M. Miranda, G. Haack, R. Silva, J. B. Brask, and N. Brunner, Markovian master equations for quantum thermal machines: local versus global approach, New J. Phys. 19, 123037 (2017).
  * Potts  _et al._ [2021] P. P. Potts, A. A. S. Kalaee, and A. Wacker, A thermodynamically consistent markovian master equation beyond the secular approximation, New J. Phys. 23, 123013 (2021).
  * González  _et al._ [2017] J. O. González, L. A. Correa, G. Nocerino, J. P. Palao, D. Alonso, and G. Adesso, Testing the validity of the ‘local’ and ‘global’ gkls master equations on an exactly solvable model, Open Systems & Information Dynamics 24, 1740010 (2017).
  * Ángel Rivas  _et al._ [2010] Ángel Rivas, A. D. K. Plato, S. F. Huelga, and M. B. Plenio, Markovian master equations: a critical study, New Journal of Physics 12, 113032 (2010).
  * Scali  _et al._ [2021] S. Scali, J. Anders, and L. A. Correa, Local master equations bypass the secular approximation, Quantum 5, 451 (2021).
  * Trushechkin [2021] A. Trushechkin, Unified gorini-kossakowski-lindblad-sudarshan quantum master equation beyond the secular approximation, Phys. Rev. A 103, 062226 (2021).
  * van der Wiel  _et al._ [2002] W. G. van der Wiel, S. D. Franceschi, J. M. Elzerman, T. Fujisawa, S. Tarucha, and L. P. Kouwenhoven, Electron transport through double quantum dots, Rev. Mod. Phys. 75, 1 (2002).
  * Sukhorukov  _et al._ [2001] E. V. Sukhorukov, G. Burkard, and D. Loss, Noise of a quantum dot system in the cotunneling regime, Phys. Rev. B 63, 125315 (2001).
  * Golovach and Loss [2004] V. N. Golovach and D. Loss, Transport through a double quantum dot in the sequential tunneling and cotunneling regimes, Phys. Rev. B 69, 245327 (2004).
  * Breuer and Petruccione [2003] H.-P. Breuer and F. Petruccione,  _The theory of open quantum systems_ (Oxford University Press, Oxford, 2003).
  * Sumetskii and Fel’shtyn [1991] M. Y. Sumetskii and M. L. Fel’shtyn, Absolute transparency of an inelastic channel and the photovoltaic effect in the resonance tunneling through the two-well heterojunction, JETP Lett. 53, 24 (1991).
  * Jauho  _et al._ [1994] A.-P. Jauho, N. S. Wingreen, and Y. Meir, Time-dependent transport in interacting and noninteracting resonant-tunneling systems, Phys. Rev. B 50, 5528 (1994).
  * Levitov and Lesovik [1993] L. Levitov and G. Lesovik, Charge distribution in quantum shot noise, JETP Lett. 58, 225 (1993).
  * Levitov  _et al._ [1996] L. Levitov, H. Lee, and G. Lesovik, Electron counting statistics and coherent states of electric current, J. Math. Phys. 37 (1996).
  * Levitov and Reznikov [2004] L. S. Levitov and M. Reznikov, Counting statistics of tunneling current, Phys. Rev. B 70, 115305 (2004).
  * Sprekeler  _et al._ [2004] H. Sprekeler, G. Kießlich, A. Wacker, and E. Schöll, Coulomb effects in tunneling through a quantum dot stack, Phys. Rev. B 69, 125328 (2004).
  * Wick  _et al._ [1952] G. C. Wick, A. S. Wightman, and E. P. Wigner, The intrinsic parity of elementary particles, Phys. Rev. 88, 101 (1952).
  * Bartlett  _et al._ [2007] S. D. Bartlett, T. Rudolph, and R. W. Spekkens, Reference frames, superselection rules, and quantum information, Rev. Mod. Phys. 79, 555 (2007).
  * Baumgratz  _et al._ [2014] T. Baumgratz, M. Cramer, and M. B. Plenio, Quantifying coherence, Phys. Rev. Lett. 113, 140401 (2014).
  * Tavakoli  _et al._ [2020] A. Tavakoli, G. Haack, N. Brunner, and J. B. Brask, Autonomous multipartite entanglement engines, Phys. Rev. A 101, 012315 (2020).
  * Tavakoli  _et al._ [2018] A. Tavakoli, G. Haack, M. Huber, N. Brunner, and J. B. Brask, Heralded generation of maximal entanglement in any dimension via incoherent coupling to thermal baths, Quantum 2, 73 (2018).
  * Tacchino  _et al._ [2018] F. Tacchino, A. Auffèves, M. F. Santos, and D. Gerace, Steady state entanglement beyond thermal limits, Phys. Rev. Lett. 120, 063604 (2018).
  * Friis [2016] N. Friis, Unlocking fermionic mode entanglement, New J. Phys. 18, 061001 (2016).
  * Dasenbrook  _et al._ [2016] D. Dasenbrook, J. Bowles, J. B. Brask, P. P. Hofer, C. Flindt, and N. Brunner, Single-electron entanglement and nonlocality, New J. Phys. 18, 043036 (2016).
  * Tan  _et al._ [1991] S. M. Tan, D. F. Walls, and M. J. Collett, Nonlocality of a single photon, Phys. Rev. Lett. 66, 252 (1991).
  * Potts  _et al._ [2016] P. Potts, D. Dasenbrook, and C. Flindt, On-demand entanglement generation using dynamic single-electron sources, Phys. Status Solidi B 254 (2016).
  * Hill and Wootters [1997] S. A. Hill and W. K. Wootters, Entanglement of a pair of quantum bits, Phys. Rev. Lett. 78, 5022 (1997).
  * Wootters [1998] W. K. Wootters, Entanglement of formation of an arbitrary state of two qubits, Phys. Rev. Lett. 80, 2245 (1998).
  * Seifert [2012] U. Seifert, Stochastic thermodynamics, fluctuation theorems and molecular machines, Rep. Prog. Phys. 75, 126001 (2012).
  * Jarzynski [2011] C. Jarzynski, Equalities and inequalities: Irreversibility and the second law of thermodynamics at the nanoscale, Annu. Rev. Condens. Matter Phys. 2, 329 (2011).
  * Bochkov and Kuzovlev [2013] G. N. Bochkov and Y. E. Kuzovlev, Fluctuation–dissipation relations. achievements and misunderstandings, Phys.-Uspekhi 56, 590 (2013).
  * Malek Mansour and Baras [2017] M. Malek Mansour and F. Baras, Fluctuation theorem: A critical review, Chaos 27, 104609 (2017).
  * Harris and Schütz [2007] R. J. Harris and G. M. Schütz, Fluctuation theorems for stochastic dynamics, J. Stat. Mech.: Theory Exp. 2007 (07), P07020.
  * Pietzonka and Seifert [2018] P. Pietzonka and U. Seifert, Universal trade-off between power, efficiency, and constancy in steady-state heat engines, Phys. Rev. Lett. 120, 190602 (2018).
  * Pietzonka  _et al._ [2016] P. Pietzonka, A. C. Barato, and U. Seifert, Universal bound on the efficiency of molecular motors, J. Stat. Mech.: Theory Exp. 2016 (12), 124004.
  * Dechant and ichi Sasa [2020] A. Dechant and S. ichi Sasa, Fluctuation-response inequality out of equilibrium, PNAS 117, 6430 (2020).
  * Pietzonka  _et al._ [2017] P. Pietzonka, F. Ritort, and U. Seifert, Finite-time generalization of the thermodynamic uncertainty relation, Phys. Rev. E 96, 012101 (2017).
  * Horowitz and Gingrich [2017] J. M. Horowitz and T. R. Gingrich, Proof of the finite-time thermodynamic uncertainty relation for steady-state currents, Phys. Rev. E 96, 020103(R) (2017).
  * Potts and Samuelsson [2019] P. P. Potts and P. Samuelsson, Thermodynamic uncertainty relations including measurement and feedback, Phys. Rev. E 100, 052137 (2019).
  * Proesmans and den Broeck [2017] K. Proesmans and C. V. den Broeck, Discrete-time thermodynamic uncertainty relation, EPL 119, 20001 (2017).
  * Koyuk and Seifert [2020] T. Koyuk and U. Seifert, Thermodynamic uncertainty relation for time-dependent driving, Phys. Rev. Lett. 125, 260604 (2020).
  * Hasegawa [2020] Y. Hasegawa, Quantum thermodynamic uncertainty relation for continuous measurement, Phys. Rev. Lett. 125, 050601 (2020).
  * Van Vu and Saito [2022] T. Van Vu and K. Saito, Thermodynamics of precision in markovian open quantum dynamics, Phys. Rev. Lett. 128, 140602 (2022).
  * Bennett  _et al._ [1993] C. H. Bennett, G. Brassard, C. Crépeau, R. Jozsa, A. Peres, and W. K. Wootters, Teleporting an unknown quantum state via dual classical and einstein-podolsky-rosen channels, Phys. Rev. Lett. 70, 1895 (1993).
  * Clauser  _et al._ [1969] J. F. Clauser, M. A. Horne, A. Shimony, and R. A. Holt, Proposed experiment to test local hidden-variable theories, Phys. Rev. Lett. 23, 880 (1969).
  * Kadanoff and Baym [1962] L. P. Kadanoff and G. Baym,  _Quantum Statistical Mechanics_ (Benjamin, New York, 1962).
  * Baym [1962] G. Baym, Self-Consistent Approximations in Many-Body Systems, Phys. Rev. 127, 1391 (1962).
  * Keldysh [1965] L. Keldysh, Diagram technique for nonequilibrium processes, Sov. Phys. JETP 20, 1018 (1965).
  * Horodecki  _et al._ [1999] M. Horodecki, P. Horodecki, and R. Horodecki, General teleportation channel, singlet fraction, and quasidistillation, Phys. Rev. A 60, 1888 (1999).
  * Saryal  _et al._ [2019] S. Saryal, H. M. Friedman, D. Segal, and B. K. Agarwalla, Thermodynamic uncertainty relation in thermal transport, Phys. Rev. E 100, 042101 (2019).
  * Ehrlich and Schaller [2021] T. Ehrlich and G. Schaller, Broadband frequency filters with quantum dot chains, Phys. Rev. B 104, 045424 (2021).
  * Petersson  _et al._ [2010] K. D. Petersson, J. R. Petta, H. Lu, and A. C. Gossard, Quantum coherence in a one-electron semiconductor charge qubit, Phys. Rev. Lett. 105, 246804 (2010).
  * Frey  _et al._ [2012] T. Frey, P. J. Leek, M. Beck, A. Blais, T. Ihn, K. Ensslin, and A. Wallraff, Dipole coupling of a double quantum dot to a microwave resonator, Phys. Rev. Lett. 108, 046807 (2012).
  * Barker  _et al._ [2019] D. Barker, S. Lehmann, L. Namazi, M. Nilsson, C. Thelander, K. A. Dick, and V. F. Maisi, Individually addressable double quantum dots formed with nanowire polytypes and identified by epitaxial markers, Applied Physics Letters 114, 183502 (2019).
  * Nilsson  _et al._ [2016] M. Nilsson, L. Namazi, S. Lehmann, M. Leijnse, K. A. Dick, and C. Thelander, Single-electron transport in inas nanowire quantum dots formed by crystal phase engineering, Phys. Rev. B 93, 195422 (2016).
  * Aamir  _et al._ [2022] M. A. Aamir, C. C. Moreno, S. Sundelin, J. Biznárová, M. Scigliuzzo, K. E. Patel, A. Osman, D. P. Lozano, I. Strandberg, and S. Gasparinetti, Engineering symmetry-selective couplings of a superconducting artificial molecule to microwave waveguides, Phys. Rev. Lett. 129, 123604 (2022).
  * Izmalkov  _et al._ [2004] A. Izmalkov, M. Grajcar, E. Il’ichev, T. Wagner, H.-G. Meyer, A. Y. Smirnov, M. H. S. Amin, A. M. van den Brink, and A. M. Zagoskin, Evidence for entangled states of two coupled flux qubits, Phys. Rev. Lett. 93, 037003 (2004).
  * Majer  _et al._ [2005] J. B. Majer, F. G. Paauw, A. C. J. ter Haar, C. J. P. M. Harmans, and J. E. Mooij, Spectroscopy on two coupled superconducting flux qubits, Phys. Rev. Lett. 94, 090501 (2005).
  * Okazaki  _et al._ [2013] Y. Okazaki, S. Sasaki, and K. Muraki, Shot noise spectroscopy on a semiconductor quantum dot in the elastic and inelastic cotunneling regimes, Phys. Rev. B 87, 041302(R) (2013).
  * [94] L. G. Molinari, arXiv:1710.09248 [math-ph] .
  * Bruus and Flensberg [2004] H. Bruus and K. Flensberg,  _Many-Body Quantum Theory in Condensed Matter Physics_ (Oxford University Press, Oxford, 2004).
  * Meir and Wingreen [1992] Y. Meir and N. S. Wingreen, Landauer formula for the current through an interacting electron region, Phys. Rev. Lett. 68, 2512 (1992).
  * Myöhänen  _et al._ [2008] P. Myöhänen, A. Stan, G. Stefanucci, and R. van Leeuwen, A many-body approach to quantum transport dynamics: Initial correlations and memory effects, EPL 84, 67001 (2008).
  * Esposito  _et al._ [2009] M. Esposito, U. Harbola, and S. Mukamel, Nonequilibrium fluctuations, fluctuation theorems, and counting statistics in quantum systems, Rev. Mod. Phys. 81, 1665 (2009).
  * Plenio and Knight [1998] M. B. Plenio and P. L. Knight, The quantum-jump approach to dissipative dynamics in quantum optics, Rev. Mod. Phys. 70, 101 (1998).
  * Schaller [2013] G. Schaller,  _Open Quantum Systems Far from Equilibrium_, Vol. 881 (Springer Cham, Cham, 2013).
  * Bruderer  _et al._ [2014] M. Bruderer, L. D. Contreras-Pulido, M. Thaller, L. Sironi, D. Obreschkow, and M. B. Plenio, Inverse counting statistics for stochastic and open quantum systems: the characteristic polynomial approach, New J. Phys. 16, 033030 (2014).
  * Puig von Friesen  _et al._ [2009] M. Puig von Friesen, C. Verdozzi, and C.-O. Almbladh, Successes and Failures of Kadanoff-Baym Dynamics in Hubbard Nanoclusters, Phys. Rev. Lett. 103, 176404 (2009).
  * Karlsson and Verdozzi [2014] D. Karlsson and C. Verdozzi, Transport of correlated electrons through disordered chains: A perspective on entanglement, conductance, and disorder averaging, Phys. Rev. B 90, 201109(R) (2014).
  * Pulay [1980] P. Pulay, Convergence acceleration of iterative sequences. the case of scf iteration, Chem. Phys. Lett. 73, 393 (1980).
  * von Friesen  _et al._ [2011] M. P. von Friesen, C. Verdozzi, and C.-O. Almbladh, Can we always get the entanglement entropy from the Kadanoff-Baym equations? The case of the T-matrix approximation, EPL 95, 27005 (2011).
  * Bell [1964] J. S. Bell, On the einstein podolsky rosen paradox, Phys. Phys. Fiz. 1, 195 (1964).