Quantum coherence-control of thermal energy transport: The V model as a case study

Here, we study a minimal model, the three-level V system coupled to two heat baths, and investigate the role of quantum coherences in heat transport in both the transient regime and in the nonequilibrium steady-state. In our model, energy is exchanged between the baths through two parallel pathways, which can be made distinct through the nondegeneracy of excited levels (energy splitting $\Delta$) and a control parameter $\alpha$, which adjusts the strength of one of the arms. Using a nonsecular quantum master equation of Redfield form, we succeed in deriving closed-form expressions for the quantum coherences and the heat current in the steady state limit for closely degenerate excited levels. By including three ingredients in our analysis: nonequilibrium baths, nondegeneracy of levels, and asymmetry of pathways, we show that quantum coherences are generated and sustained in the V model in the steady-state limit if three conditions, conjoining thermal and coherent effects are simultaneously met: (i) The two baths are held at different temperatures. (ii) Bath-induced pathways do not interfere destructively. (iii) Thermal rates do not mingle with the control parameter $\alpha$ to destroy interferences through an effective local equilibrium condition. We find that coherences are maximized when the heat current is suppressed. On the other hand, the secular Redfield quantum master equation is shown to fail in a broad range of parameters. Although we mainly focus on analytical results in the steady state limit, numerical simulations reveal that the transient behavior of coherences contrasts the steady-state limit, suggesting that different mechanisms are at play in these two regimes. Enhancing either the lifetime of transient coherences or their magnitude at steady state thus requires the control and optimization of different physical parameters.

The V level scheme is one of the simplest models that can be analyzed in revealing the role of quantum coherences in dynamical and steady state properties. This model, and the related Λ model, have been analyzed in great details in the context of quantum optics [32][33][34][35]41,[45][46][47][48][49][50][51] , and have more recently received renewed attention motivated by insights they offer in the study of light-harvesting systems [52][53][54][55][56][57][58] . The V model is composed of three levels: A ground state and two quasidegenerate excited states, with coherences generated between the two excited levels due to driving by coherent or incoherent environments. Many variants of this model exist; in our work, the ground state is coupled to both excited states through two heat baths, see Fig. 1. Each transition (ground state ↔ each excited state) leads to heat flow between the hot and cold baths. However, the two transitions (paths) can be tuned to interfere constructively or destructively with each other.
Long-lived transient quantum coherences were shown to be generated in the V model when coupled to a single heat bath, see e.g., Refs. [53][54][55][56][57][58] . Though of significant magnitude and lifetime, these quantum coherences eventually decay and vanish in the steady-state regime, limiting potential applications in e.g. autonomous quantum heat engines. Quantum coherences were demonstrated to be sustained in the steady-state limit once the V model was coupled to nonequilibrium reservoirs or when detailed balance was broken by other means 35,52,[59][60][61][62][63][64] . In both transient and steady state regimes, coherences can be generated due to incoherent processes instead of benefiting from coherent sources such as laser or cavity confinement 65 , in which case transfer of coherence (e.g., from the laser to the system) more appropriately describes the phenomenon. Despite this similarity, transient and steady state coherences are controlled by different parameters, as we shall demonstrate in this work.
In electron transport junctions, steady-state quantum coher-arXiv:2207.09512v1 [quant-ph] 19 Jul 2022 ences have been explored intensively, both theoretically [66][67][68][69][70][71][72] and experimentally [73][74][75][76][77][78] ; references here are examples of a rich literature. In a typical setup, two quantum dots are coupled to two metals through separate arms, a reference arm, and an adjustable arm, constituting the realization of an electron interferometer. The phase difference between the two paths can be tuned by a magnetic flux via the Aharonov-Bohm effect 79 . The type of interference between the two arms controls the extent of quantum coherences between the quantum dots, and thus the steady state charge current. Using the nonsecular Redfield master equation, which is a microscopic technique that treats quantum coherences faithfully, here we demonstrate the connection and analogy between the V model and the double-quantum-dot setup.
In this work, we use the V model with quasi-degenerate excited states as a thermal-conducting element, mediating heat transport between two thermal baths. Since heat is transported through two "arms", coherences develop between the excited states. By controlling the balance between the arms, interference effects are modulated to enhance or suppress the heat current. The V model thus demonstrates coherence-control of thermal energy transport. We describe this effect in both the transient and steady state limits. Our key new results are: (i) We derive closed-form expressions for coherences, populations and the heat current in the steady state regime using a nonsecular quantum master equation, identifying conditions and the extent of control of heat transport via coherences. As we show below, coherences between excited states are nonzero in steady state due to the interplay of coherent and thermal effects. Their real part satisfies with the imaginary part proportional to the real part. Here, α, a dimensionless real parameter, controls the asymmetry between the interfering arms. n h,c (ν) is the Bose-Einstein occupation factor of the h, c baths and J h,c (ν) is the spectral density function of the heat baths. Both functions are calculated at the energy gap ν, between ground and excited states. Similarly, we show that the steady state heat current can be tuned by α, Note that the coherences and the heat current are nonmonotonic in α, as we show in Sec. V with full expressions and simulations.
(ii) We confirm numerically the thermodynamic consistency of the nonsecular Redfield equation in the present model through an analysis of the dynamical map for the explored range of parameters.
(iii) We demonstrate stark differences between transient dynamics of the V model and its steady state behavior, suggesting that different parameters are at play for either generating or sustaining coherences: While the energy splitting between the closely-degenerate excited stated dictates the lifetime of transient dynamics, this energy has only a secondary, small role in dictating steady state values.
This work is organized as follows. In Sec II, we describe the V model. In Sec. III, we outline the Redfield quantum master equation that we use to study dynamics and steady state transport. We present simulations of the transient dynamics in Sec. IV and analytic results in the steady state limit in Sec. V. In Sec. VI we discuss (i) the local-site picture and its physical interpretation of quantum coherences, (ii) the thermal diode effect, and (iii) the thermodynamic consistency of our approach. We conclude in Sec. VII. Figure 1. Schematic representation of the V model with quasidegenerate excited states, ∆ ν. The two heat baths (T h > T c ) allow transitions between the ground state |1 and the two excited states, |2 and |3 , thus realizing a two-arm interferometer: Heat is transferred from the hot (T h ) to the cold (T c ) bath through two paths, which are distinguishable when |α| = 1, for vanishing ∆. The sign and magnitude of α, a real-valued tunable parameter, determines whether the interference pattern in steady state is perfectly constructive (α = −1) or destructive (α = 1).

II. MODEL
The nonequilibrium V model (Fig. 1) is a minimal model that exhibits nontrivial quantum coherences in its operation as a thermal junction, at weak system-bath coupling. In the energy basis, the Hamiltonian of the model is given bŷ We adopt here natural units,h = 1, k B ≡ 1. To facilitate interference effects, we work in the limit of ∆ ν. Transitions between states of the system, |1 ↔ |2 and |1 ↔ |3 are enacted by two independent bosonic heat baths denoted by h (hot) and c (cold) with T h > T c . The Hamiltonian of bath m iŝ Here,b † k,m (b k,m ) is the creation (annihilation) bosonic operator of a mode k of frequency ν k,m in bath m. The system-bath interaction Hamiltonian is given in a bipartite form, with a system operatorŜ m coupled to a bath operatorB m , λ k,m describes the system-bath coupling energy between mode k in the mth bath and the system. To control the interference pattern between the two arms, we introduce an asymmetry in transition couplings with a real-valued dimensionless parameter α, so thatŜ Here, h.c. is a hermitian conjugate. Note that a more generalized form for this coupling would replace α by a complex number αe iθ . In the spirit of the Aharonov-Bohm interferometer,Ŝ h corresponds to the reference arm whileŜ c is the adjustable arm. In sum, the total Hamiltonian is given bŷ Numerous studies had approached this system in the energy basis ofĤ S 35,52-60,63 . However, as we discuss in Sec. VI, interference effects can also be visualized in the local basis, as was done, e.g., in Ref. 80 .

III. METHOD: REVIEW OF THE REDFIELD EQUATION
The main objective of this work is to understand how quantum coherences can be used to tune and control heat currents in a minimal model of a quantum heat conductor. This is achieved by deriving closed-form expressions for quantum coherences and heat transport in the steady state limit. These expressions are expected to serve as the groundwork for understanding more complex systems (e.g., those amenable to strong coupling effects and non-Markovian dynamics).
To allow analytical work, we contain ourselves with the second-order Markovian quantum master equation of the Redfield form, which is not secular, i.e., population and coherences are coupled 81 . The assumptions behind Markovian Redfield equations are that (i) the system and the baths are weakly coupled thus their states are separable at any instant (Born approximation), and (ii) that the baths' time correlation functions, which dictate transition rates in the system, quickly decay, thus the dynamics of the system is Markovian. The resulting Markovian Redfield equation for the reduced density matrix (RDM) σ (t) in the Schrödinger representation is given byσ Here, ω ab = E a − E b are the Bohr frequencies with E a being an eigenenergy of the system, see Eq.
(3). The index m identifies the bath at issue (cold, hot). The Redfield dynamics can be also written in a compact form aṡ , with the dissipators D m (σ (t)). Terms in the dissipators are given by half Fourier transforms of autocorrelation functions of the baths, Here, operators are given in their interaction representation with respect toĤ S + ∑ mĤB,m . The matrix elements (S m ) ab are given in Eq. (6). The real and imaginary parts of the correlation functions can be written in terms of the spectral density function of the bath, J m (ω) = π ∑ k λ 2 k,m δ (ω − ν k,m ). For the interaction model Eq. (5), the real part k m (ω) 2 reduces to As for the imaginary part, Z m (ω), we neglect it in this work assuming that it is small at high enough temperatures. The expressions that we derive below do not rely on a specific form for the spectral density function, J m (ω). In numerical simulations we assume Ohmic spectral density functions with an infinite high frequency cutoff, J m (ω) = γ m ω with γ m a dimensionless coupling parameter taken identical for both baths. The occupation function of the baths follows the Bose-Einstein distribution, n m (ω) ≡ [e β m ω − 1] −1 with β m = 1/T m as the inverse temperature. Locally, the transition rates obey the detailed balance relation, k m (ω) = e β m ω k m (−ω).
To calculate heat currents, we study the energy of the system, Ĥ S (t) ≡ Tr S [σ (t)Ĥ S ]. Using the dynamical equation for the RDM, the rate of energy change of the system is Since the dissipators are additive, we identify the heat current flowing into the system from the hot bath as j h (t) = Tr S Ĥ S D h σ (t) . A similar definition holds for j c (t). In the long time limit, j h = Tr S Ĥ S D h σ and j h = − j c ; the heat current is defined positive when flowing from the mth bath towards the system. Note that we suppress the time variable to indicate steady state quantities. To evaluate the current in the steady state limit, we first solve the algebraic equations for the RDM by settingσ = 0 with the normalization condition Tr[σ ] = 1, then substitute the results into the dissipator.
In what follows, we study the behavior of the system in both the transient regime and the steady state limit, with and without coherences: (i) We solve numerically Eq. (8) and gather both the dynamics of the system and its steady state behavior. (ii) In the steady state limit, we simplify the nonsecular Redfield equation and derive closed-form expressions for σ and j h in the limit of small, yet nonzero ∆. (iii) We further make the secular approximation and obtain analytic results in the steady state limit. Comparisons between (i) and (ii) illustrate fundamental differences between transient dynamics and steady state behavior. Comparing results from (ii) to (iii) reveals the role of coherences on steady state heat transport in the V model.

IV. TRANSIENT DYNAMICS
We adopt the Redfield equation (8) to simulate the system's dynamics; Appendix A includes the specific equations of motion for the V model. As the initial condition, we assume that only the ground state is populated σ 11 (t = 0) = 1 while all other elements of the RDM are zero. However, we confirmed with simulations that observations were generic for other initial conditions. (We are not considering here exotic initial con-ditions that lead to anomalous open-system relaxation as in the Mpemba effect 82 ). Results are presented in Fig. 2, where we study σ R 32 (t) and the heat current j h (t) from the initial time to the steady state limit.
First, in Fig. 2(a), we focus on the behavior of coherences between the excited states. We note an intriguing contrast between the transient regime and the steady state limit: When transient coherences are long-lived at α = 1, they eventually vanish in the long-time limit. In contrast, for other values of α such as for α = −1, temporal coherences are smaller in magnitude compared to for α = 1, yet coherences survive in the steady state regime. Long-lived temporal coherences were explored in e.g. Ref. [53][54][55][56][57][58] , yet under equilibrium (single bath) conditions. As such, these studies could not observe the curious contrast shown here, between temporal dynamics and steady-state behavior.
The temporal behavior of the heat current is presented in Fig. 2(b). We observe the following: (i) Long-lived coherences temporarily support higher currents, which decay to their steady state values only when the coherences are suppressed. This effect is most clearly demonstrated for α = 1, but it is still visible for other values, α = 0.8, 1.2. (ii) As for the steady state behavior: The current diminishes to a small value for α = −1. In fact, in Sec. VI we show that for α = −1 the current scales as ∆ 2 , thus when the levels approach degeneracy, j h → 0.
We continue and study in Fig. 3 the excited state populations, focusing on the population difference σ 33 (t) − σ 22 (t).
Since E 3 > E 2 , a positive difference indicates population inversion. Interestingly, we once more observe contrasting behavior between the transient and the steady state regimes: For |α| ≤ 1 population inversion appears only in the steady state limit, but not in the transient regime. In contrast, for |α| > 1, population inversion develops at short time, but disappears in the long time limit.
It was shown for the V model in equilibrium that the lifetime of transient effects scales inversely with ∆ 2 , see e.g., Refs. [53][54][55][56][57][58] . As we show next, in contrast, ∆ has a marginal role on steady state values. Regardless, one must include it for properly deriving results.

V. ANALYTIC RESULTS IN STEADY STATE
The Redfield equation for the V model (Appendix A) is cumbersome to solve analytically-even in the steady state limit. To make progress, we simplify it as follows: Since ∆ ν, it is reasonable to assume that the decay rates from either excited states to the common ground state, activated by the mth heat bath, are about equal, k m (ω 21 ) ∼ k m (ω 31 ), evaluated at ν. An analogous approximation holds for the excitation rates. Recall that since ν > 0, the decay rate is given by see equation (10). To simplify the notation, henceforth we withhold the frequency value and use the short notation k m to denote the relaxation rate induced by the mth bath. Excitation rates are written using the detailed balance relation. The resulting equations for the RDM, together with the condition of conservation of population, arė Here, σ R 32 (t) corresponds to the real part of coherences. For details, see Appendix A. In the long time limit, this set of equations can be solved analytically as an algebraic system by settingσ = 0 together with the population normalization condition.

A. Quantum Coherences
In the steady state limit, we find from Eq. (13) that the real (R) and imaginary (I) parts of the coherences are related via Therefore, taking the real part of Eq. (13) at steady state, we get Recall that σ without an explicit time dependence refers to the density operator in the steady state limit.
Equations (14) and (15) in the long time limit, together with Eq. (17) and the population normalization condition provide a closed-form expression for quantum coherences between the two excited states, Here, Equation (18) can be written using the microscopic expression for the rates, Equation (18), or in its microscopic form, Eq. (1), is the first main result of this work. It reveals that quantum coherences are non-vanishing in the steady-state limit if the following three conditions are met: (i) Quantum interferences are nondestructive, with α = 1.
(ii) The system is at a nonequilibrium steady state, T h = T c .
(iii) k h = αk c . The equality defines a special point in parameter space. At this point, the hot and cold bath-induced rates, between the ground state and level |3 , are equal, defining a local equilibrium condition. This nontrivial special point can be reached either by controlling the temperatures of the baths and their spectral properties, or by tuning the strength of the adjustable "arm" through the α parameter.
Thus, to maintain steady state coherences, conditions (i) and (ii) separately address the state of the baths (out of equilibrium) and the interference pathway (α = 1). In contrast, condition (iii) requires that the two aspects, α-control and the out-of-equilibrium setting dot not collectively compensate each other and create a local equilibrium situation.
Numerical results are presented in Fig. 4. Destructive interferences take place at α = 1 reflected by a dip whose width is controlled by ∆. When ∆ vanishes, the (α − 1) 2 factor cancels and condition (i) is obscured. Moreover, steady-state quantum coherences displayed by fully degenerate V-systems may arise through a distinct mechanism, interpreted as transient coherences that never decay, see Refs. [53][54][55][56][57][58] , even in equilibrium settings. It is therefore critical to perform the analysis at nonzero ∆ to correctly capture the interference behavior around α = 1. Coherences reach a maximum due to a constructive interference effect when α = −1. The special point k h = αk c at which the coherences vanish is arrived here at α ≈ 1.8.
Our expression Eq. (18) demonstrates clearly that a proper description of nonequilibrium steady state coherences requires the inclusion of all three aspects: (i) Nonequilibrium settings, T h = T c . (ii) Nondegeneracy of the excited states with ∆ small but nonzero. (iii) Asymmetry of the setup, incorporated courtesy of the control parameter α adjusting the arms. Although prior studies succeeded in deriving closed-form solutions of nonequilibrium steady state coherences, they did not capture the complete solution Eq. (18) as some of these aspects were neglected. For example, asymmetry was not included in Ref. 61 . As a result, their solution for steady state coherences only depended on thermal occupations of the baths and the eigenenergies of the V-system, missing the conditions for nonvanishing coherences α = 1 and α = k h /k c . Ref. 63 similarly assumed α = 1, although coherences persisted in steady state through breaking detailed balance via polarized incoherent light. Other studies assumed full degeneracy of the excited states, see e.g., Ref. 62 ; we stress that setting ∆ = 0 implies that the dip feature near α = 1, as shown in Eq. (18) and Fig. 4, is missed, and thus one may infer incorrectly that steady state nonequilibrium coherences is nonvanishing even in symmetrical couplings. The three aspects discussed here: nonequilibrium, nondegeneracy, and asymmetry were included in Ref. 60 , but a closed-form solution was not presented there.
Complementing this subsection, in Appendix B we show that the previously-analyzed equilibrium V model is a limiting case of our equations of motion. In Appendix C we further look at the special symmetrical point α = 1, where steady state coherences vanish even away from equilibrium.

B. Population
We now analyze excited state populations in steady state and demonstrate their dependence on the adjustable interference parameter, α. Solving Eqs. (14) and (15) in steady state, together with Eq. (17), we obtain and The steady-state populations, σ 22 , σ 33 , and thus σ 11 are insensitive to the sign of α. Furthermore, adding equations (20) and (21) reveals an analog to the so-called "populationlocked-states" effect discussed e.g. in Ref. 35 , A fraction of the population, given by Eq. (22), is always shared between the two excited states, |2 and |3 , independent of α.
Complementing this analytical result, which is valid when ∆ → 0, in Fig. 5 we present numerical simulations for P ex performed at finite ∆. As ∆ diminishes, we observe a sharp spike when approaching α = 1, at a point where strong interferences between two pathways takes place leading to zero coherences. Interestingly, at this point the total excited state population P ex coincides with the prediction of the secular equation. This anomaly is otherwise hidden in the O(∆ 2 ) term in Eq. (22) (the width of the spike depends on ∆). As we increase ∆, we find that the prediction of the secular method approaches the nonsecular result, while the analytic result Eq. (22), which relies on ∆ → 0, naturally fails to capture the correct behavior of P ex .
From Equations (20) and (21) we find that both excited states are equally populated when α = 1 or α = −1, up to the first order in ∆. . Therefore, we can write the state the of the nonequilibrium system as the algebraic average of the two thermal states, HereĤ S is the system Hamiltonian of the fully degenerate V model; C contains the off-diagonal terms of the RDM. Furthermore, when T h = T c , we conclude from Eqs. (20) and (21) that Thus, α can be tuned to force the system to exhibit population inversion at steady state (although the effect is small with our parameters as shown in Fig. 3). This observation is another nontrivial result of this work. When prepared in the ground state, dynamical simulations in Fig. 3 also show that there exists a time interval τ where σ 22 (τ) < σ 33 (τ) before reaching σ 22 > σ 33 , and vice-versa.   26) and (27), respectively. Parameters are the same as in Fig. 2.

C. Heat current
The heat currents j h,c are evaluated at the hot and cold contacts, respectively, using the dissipators, see Eq. (11). Specifically at the hot contact, Inserting the solutions for the coherences, Eq. (18), and the populations, Eqs. (20) and (21) into the dissipator, we get By making the substitution Eq. (12) we highlight that the steady state current is non-vanishing when T h = T c and α = −1, and obtain Eq. (2). Note that the α-dependence cancels out in Eq. (26) when αk c = k h , a point at which one of the pathways becomes equilibrium-like. Equation (26), or in its microscopic form, Eq. (2), complements our expressions for the coherences and populations, and it is another central result of this work. We expect that it will be used as the groundwork for more involved studies, e.g., when studying thermal energy transport at strong coupling. We plot the steady-state heat current in Fig. 6, showing that the above expression is accurate for quasi-degenerate excited states. Inspecting Eq. (26) and Fig. 6, we observe the following: (i) The nonsecular heat current is not an even function of α, a clear indication that quantum coherences are at play, with the sign of α playing a role in transport. In contrast, under the secular approximation, the current is symmetric in α.
(ii) The [e −β h ν − e −β c ν ] factor in the numerator exposes the trivial dependence of the heat current on the temperature difference.
(iii) Given the structure of the denominator, it is clear that exchanging the temperatures of the two baths should lead to a thermal diode effect, since the current is not identical under the k h ↔ k c operation (except for |α| = 1). We discuss this aspect in Sec. VI.
(iv) The steady state heat current is high when quantum coherences are low, and vice versa; compare Fig. 4 to Fig. 6.
(v) To the lowest order, the steady state heat current does not depend on the excited-state splitting ∆. Similarly, the levels' population and their coherences are independent of ∆ to lowest order (assuming ∆ → 0). Parameters that do dictate these functions are the spectral density function, temperature of the baths, and the overall splitting ν. In stark contrast, the rise of coherences and the development of an extended transient region over which they survive is controlled by a timescale inversely proportional to ∆ 2 , see Refs. [53][54][55][56][57][58]63,83 . Thus, the lifetime of transient coherence, and their magnitude at steady state are dictated by a different set of parameters.

D. The secular approximation
The secular approximation decouples population dynamics and coherent dynamics, amounting to crossing out σ 32 (t) in Eqs. (14) and (15). It is justified when ∆ is large such that the characteristic time of coherent oscillations is much shorter than the population decay time. The steady state solution thus involves a 2 × 2 coefficient matrix for the population of the excited states. At steady state, the heat current is given by (S stands for secular) Most notably, this expression is even in α, unlike the nonsecular result, Eq. (26), thus reflecting that only the magnitude of α controls the secular current, rather than its sign. This behavior is illustrated in Fig. 6. In the special cases |α| = 1, Eq.
This expression recovers the result obtained in Ref 80 where α = −1 was used throughout. Interestingly, Fig. 6 shows that the secular approximation works reasonably well when α 1. Another nontrivial observation from Fig. 6 is that for a broad range of parameters, quantum coherences reduce the heat current in the V model, with the secular current exceeding the nonsecular one. Coherences are either noninfluential to the current, or they suppress it. This can be viewed as a coherent heat current flowing against the temperature difference 61 .

A. Quantum interferometers and local-basis mapping
The V model as depicted in Fig. 1 is analogous to a double quantum dot interferometer. However, instead of electron transport through two parallel dots, here thermal energy transfer takes place through quasi-degenerate excited states of frequencies ν − ∆ and ν. In our work, the driving force for the transfer process are incoherent baths, rather than a coherent laser, but the mechanism is similar to Fano interferences 84 , also considered in Refs. 16,17,22,[53][54][55][56][57][58]63 .
The effects observed in this work can be also rationalized without requiring interference arguments, by working in the local basis 80 , referred to as the V-Energy Transfer System (VETS). The V model, Fig. 1, can be unitarily transformed into a local (L) picture, The three levels in this system are the ground state |1 (same as in the energy basis) and the fully degenerate excited states |a and |b with the tunneling energy ∆/2. This system couples to heat baths viâ see Appendix D for further details. The two limiting configurations, α = 1 and α = −1, are sketched in Fig. 7. The α = 1 case corresponds to a "sidecoupled" model in the context of e.g. double quantum dots 85 . In this case, the transition |1 ↔ |a is driven by both heat baths, but level |b is accessible only through a coherent coupling from level |a . In contrast, the α = −1 case is reminiscent of a serial double dot model. Here, the ground state is coupled to level |a through the hot heat bath. Since the two excited states are coherently coupled, energy transfers coherently from |a to level |b , followed by a relaxation to the ground state and energy transfer to the cold bath usingŜ L c . In this picture, it is obvious that the heat current vanishes for α = −1 once ∆ = 0.
Consistent with this description, the steady-state heat currents in Fig. 8 demonstrate distinct scaling with ∆ for α = −1 and α = 1: In the former case, the current scales as ∆ 2 . On the other hand, the α = 1 model only minimally depends on ∆. Note that when ∆ is very large (approaching the total gap ν), deviations from these scaling are observed due to the large modification of energy levels, see Ref. 80 .

B. Thermal rectification
Thermal rectification is a transport phenomenon whereby the magnitude of the heat current is different upon exchanging the direction of the applied temperature bias. In the quantum domain, early proposals for a thermal diode (rectifier) were built on the nonequilibrium spin-boson model 86,87 . Recently, a quantum thermal diode was experimentally realized using superconducting quantum circuits 88 .
The diode effect is realized in systems that include (i) anharmonicity and (ii) asymmetry [86][87][88][89][90] . Anharmonicity is captured in the V model due the three-level structure, distinct from the harmonic oscillator spectrum of the baths. Spatial asymmetry is introduced into the V model due to the imbalance between the arms, courtesy of the coherence parameter α 86 . The figure of merit used to quantify the thermal diode effect is the rectification ratio, R = j h (∆T ) j h (−∆T ) , with rectification marked by R = 1.
For quasi-degenerate levels, ∆ → 0, the nonequilibrium V system acts as a thermal diode once α = ±1: For α = 1 (see Fig. 7), the model is spatially symmetric and thus trivially it cannot act as a diode. For α = −1, the current itself vanishes when ∆ → 0. The impact of interferences on thermal rectification and heat amplification was investigated in the fully degenerate nonequilibrium V model, by setting ∆ to zero 62 . In contrast, in Fig 9 we study the diode effect in a nondegenerate system with tunable α. We observe the following: First, Fig 9(a) shows that steady-state coherences are asymmetric in α, as was also observed in Fig. 4. Specifically, there are striking differences in the behavior of coherences with ∆T = T h − T c between α = −0.5 and α = 0.5: While in the former, coherences increase monotonically with ∆T , this is not the case for the latter. Second, a larger temperature bias produces greater coherences in absolute value. Third, there is an asymmetry in the magnitude of coherences with respect to temperature biases: a positive bias mostly produces larger coherences than a negative temperature bias.
Next, in Fig 9(b1)-(b2) we focus on the heat current and the rectification ratio. For α = ±1, no rectification is observed. For intermediate values of α, rectification behavior is observed, and is maximal when α = 0. Regardless, rectifica- tion ratios achieved in this model are rather modest, reaching at most R ∼ 1.7. It remains a challenge to realize a significant nanoscale quantum thermal rectifier.

C. Thermodynamic consistency of the Redfield equation for the V model
Quantum master equations in the Redfield form are known to exhibit pathologies [91][92][93] . Most notably, the Redfield equation is in principle not a completely-positive and trace preserving map 94 , allowing violations of thermodynamic laws as observed in the transient domain 95,96 . To support our analysis in this paper, in Fig. 10 we study properties of the Redfield quantum master equation for the V model, demonstrating numerically the thermodynamic consistency of our results.
First, in Fig. 10(a) we test the conservation of energy by studying the sum of the two heat currents, j h (t) + j c (t) for different values of α. Recall our sign convention: currents are positive when flowing towards the system. In the steady state limit, we find that this sum approaches zero, validating the first law of thermodynamics of energy conservation. Next, in Fig. 10(b) we test the second law by displaying the entropy production rateΣ(t) as a function of time, Throughout the time evolution, the entropy production rate is positive, and thus the system obeys the second law. Finally, in Fig. 10 (c1)-(c3) we test whether for the V model the Redfield equation is a completely-positive trace-preserving map. We find that throughout the dynamical evolution, the eigenvalues of the reduced density matrix are strictly positive, ensuring complete positivity. The sum of the eigenvalues is unity, implying trace preservation. Pathologies associated with the Redfield equation therefore do not manifest in the V model for the studied parameters.

VII. SUMMARY
We studied a quantum thermal junction in which quantum coherences play a decisive role. While in the popular nonequilibrium spin-boson nanojunction coherences in the system do not alter transport (in the weak system-bath coupling limit), in the V model coherences are generated by the nonequilibrium heat baths, and their transient and steady state behavior go hand in hand with that of the heat current. Our main achievements are: (i) In the limit of quasi-degenerate excited states we used the nonsecular Redfield equation and derived closed-form expressions in steady state identifying conditions for nonzero coherences and heat currents with a tunable interference parameter α. We found that quantum coherences were nonvanishing in the steady-state limit due to a combination of coherent and incoherent factors: (a) the V model was in a nonequilibrium steady-state, (b) the parallel transport pathways did not destructively interfere, and (c) the parallel pathways did not cancel out due to an effective equilibrium condition on one of the arms. Quantum coherences were generated and sustained in our model due to the incoherent baths, and this method may find practical applications in emerging quantum information technologies.
(ii) The Redfield approach is more cumbersome to handle than the Lindblad equation, and in principle is not always granted the luxury of complete positivity. Here we numerically showed that for the V model and the examined parameters, the Redfield map was completely positive and trace preserving.
(iii) By studying the behavior of the heat current under the secular approximation we identified physics missed by making this common approximation.
(iv) While our focus has been on the steady-state limit, we had further performed numerical simulations in the transient domain. We found that the behavior of coherences in the transient regime was markedly different compared to the steady state limit, indicating that different physical factors were at play. Indeed, while the timescale of the transient dynamics is inversely related to ∆ 2 , see Refs [53][54][55][56][57][58]63,83 , steady state coherences, population and currents only secondarily-negligibly depended on ∆ in the quasi-degenerate limit.
In sum, we derived a closed-form expression for the heat current, Eq. (26), which depended on the quantum interference control parameter (α). We expect that this analytical result will serve as a starting point for additional investigations of quantum coherence effects in thermal transport. Open experimental challenges include the construction of phase-controlled quantum thermal devices realizing this system with generalized complex couplings 97 . From the the-ory side, open questions concern the interplay of quantum coherences and non-Markovianity, understanding the impact of strong system-bath coupling on quantum coherences, and uncovering the relationship between quantum coherences and fluctuations.

ACKNOWLEDGMENTS
DS acknowledges the NSERC discovery grant and the Canada Research Chair Program. The work of FI was funded by the University of Toronto Excellence Award. We acknowledge discussions with Marlon Brenes.
where m = h, c is the index for the heat bath and the frequencies are calculated from the eigenenergies of the system, ω i j = E i − E j . Note that by construction, there are no coherence terms between the ground state and the excited states. To calculate the dissipator, we use Eq. (9) ignoring the so-called Lamb shifts, an approximation that should be valid when T > ∆. We further make use of the detailed-balance relation and write all rates in terms of bath-induced decay rates, k(ω) with ω > 0 [see Eq.
(10)]. The coherences satisfẏ The excited state population dynamics includes gain-loss terms from-to the ground state, as well as the coupling to coherences, A similar equation of motion holds for the other excited state, We recall that these equations rely on (i) the Born-Markov approximation and (ii) the omission of the imaginary part of the rates ("Lamb shifts"). However, unlike the Lindblad form, the equations arise from microscopic principles and they are nonsecular. Simulations presented in the main text are based on solving numerically equations (A4)-(A6).
To solve these equations analytically, further simplifications are necessary. Working in the limit of nearly-degenerate excited levels, ∆ → 0, we make an additional approximation for the thermally-induced decay rates, k h (ω 31 ) ≈ k h (ω 21 ), denoted in short by k h . A similar approximation is made for k c . All rates are now evaluated at the transition energy ν, see Fig. 1. However, the unitary term i∆σ 32 (t) in Eq. (A4) is kept intact. The resulting simplified equations are presented in Eqs. (13)- (15). In steady state, after using the population normalization condition, we recast the system in an algebraic form as with the steady state solution v ≡ [σ 22 σ 33 σ R 32 ] T . The matrix is and the constant terms are The steady state solution is obtained by inverting the coefficient matrix M.
To obtain the analogous set of secular equations, coherent terms are set to vanish in Eqs. (A5) and (A6) then we proceed as with the nonsecular analog to obtain Eq. (27). For completeness, we provide here expressions for population, and where recall that the S superscript denotes the secular limit.

APPENDIX B: EQUILIBRIUM SINGLE-BATH V MODEL
Here, we investigate the V model when coupled to a single bath. In this case, the system is expected to equilibrate to a Gibbs state in the long time limit. Explicitly, we set k c = 0 in Eqs. (13)- (15) and k h = k as the decay rate, These equations of motion map to a limiting configuration of Eqs. (1a) and (1b) of Ref. 54 as follows. First, in their notation, r = r 1 = r 2 is the pumping rate while γ = γ 1 = γ 2 is the spontaneous decay rate.
APPENDIX C: THE NONEQUILIBRIUM V MODEL AT α = 1 As we show in Fig. 4 and Eq. (18), out of equilibrium baths generate and sustain coherences in the excited states in the steady state limit, but this coherence nullifies when α = 1. We now show that the equations of motion at α = 1 map to the single-bath case.
We simplify Eqs.

APPENDIX D: LOCAL PICTURE FOR THE V MODEL
In this Appendix we provide more details on the transformation between the V model in the energy basis and a particular local basis picture as discussed in Ref. 80 .
In the local basis, |1 denotes the ground state and |a and |b are exactly-degenerate upper states, which are coherently coupled with the tunneling energy ∆/2. It can be shown that the local Hamiltonian The eigenvalue matrix corresponds to the V model ("global") Hamiltonian Eq. (3). The system operators that couple to the baths go through the same transformation yielding the system-bath coupling operators in the local basis, While the hot bath allows only the transition from the ground state to |a , depending on α, the cold bath may allow transitions to both excited states. The limit α = 1 corresponds to the "side-coupled model" as level |b does not couple to the ground state. In contrast when α = −1 we reach the scenario of Ref. 80 where each excited state couples to the ground state through a different bath. In the latter case the tunneling energy ∆ is essential for allowing energy transfer thus the current vanishes when ∆ → 0 80 .