Quantum interference between independent reservoirs in open quantum systems

When a quantum system interacts with multiple reservoirs, the environmental effects are usually treated in an additive manner. We show that this assumption breaks down for non-Markovian environments that have finite memory times. Specifically, we demonstrate that quantum interferences between independent environments can qualitatively modify the dynamics of the physical system. We illustrate this effect with a two level system coupled to two structured photonic reservoirs, discuss its origin using a non-equilibrium diagrammatic technique, and show an example when the application of this interference can result in an improved dark state preparation in a $\lambda$ system.


I. INTRODUCTION
The study of environmental effects on a quantum mechanical system is at the heart of quantum information science. An environment generally results in quantum decoherence [1,2]. At the same time, it can be used to control the dynamics of a system through quantum reservoir engineering [3][4][5][6][7], or generate interesting new phases in many-body systems [8][9][10][11][12]. The theoretical paradigm to study open quantum systems often relies on the master equation (ME) approach [13]. In practice, a system often couples to multiple reservoirs as is the case of cavity quantum electrodynamics (CQED) [14][15][16], Jaynes-Cummings lattices [17], photon-ion interfaces [18], ion chain systems [19], phonon-induced spin squeezing [20], the dynamical Casimir effect [21], and so on. While it is commonly believed that different baths are additive in the ME in the Markovian limit, this assumption is not always valid for the more general non-Markovian physical situations [22][23][24][25][26][27][28][29][30][31][32][33]. There are various definitions and measures of non-Markovianity. The identification of non-Markovianity from the structure of a ME can be found in Ref. [34].
In this paper, we demonstrate and explain the emergence of quantum interference between independent baths and the corresponding impact on the reduced system. We show that independent baths can produce correlations that cannot be described by an additive ME, in both short and long time regimes, except if all the baths are Markovian. This can be understood in the framework of Keldysh nonequilibrium diagrammatic theory [35], in which an additive ME misses all diagrams that represent interferences among multiple baths.
As an example, we show that one can achieve a better dark state preparation of a driven system by taking the bath interference effects into account.
Our paper is outlined as follows. In Sec. II, we first illustrate this interference effect using a specific example of a two-level system (TLS) coupled to two independent and structured photonic environments. This is followed by a general theoretical description in Sec. III that shows how the environmental interference emerges dynamically when an open quantum system interacts with multiple reservoirs. We will discuss both the diagrammatic approach and the projection operator ME technique. In Sec. IV, we present the applicable example of dark state preparation. We summarize in Sec. V and provide some technical details in the Appendix.

A. Additive ME
It is usually assumed that the dynamics of a system coupled to independent baths is governed by an additive ME. An additive ME for the reduced density matrix ρ s (t) of the system has the formρ where the superoperatorL i describes the noise effect due to bath i. For a Markovian environment,L is a standard timeindependent Lindblad operator, while in the non-Markovian regime,L gains time dependence and can be a time integral of the memory kernel operator [13]. The issue of an additive ME can be best illustrated by a simple theoretical model: a two-level system interacting with two uncoupled photonic reservoirs as shown in Fig. 1(a). The canonical interaction is given by where g i,k is the coupling constant, b † i,k creates a photon in bath i, σ − = |g e|, and ω 0 is the TLS energy splitting. The photonic environments are independent, so that [b 1,k ,b † 2,k ] = 0. The two baths can have very different decay rates, linewidths, and central frequencies. When the baths are initially in the vacuum state, the dynamics of the TLS is determined by the spectral densities defined as J i (ω) = dt 2π k g 2 i,k e i(ω−ω k )t [36]. A constant J (ω) corresponds to the Markovian limit and a non-Markovian spectral density in general has a frequency dependence, e.g., the dynamics of a CQED system can be modeled by a Lorentzian J (ω). In the Each photonic bath is characterized by a non-Markovian spectral density J i (ω). (b, c) The population dynamics of an initially excited TLS based on the additive ME and the exact solution, respectively, is presented. The ME solution produces a much suppressed oscillation due to the missing of the bath interference. The parameters used are λ 1 /γ 1 = 0.01, λ 2 /γ 1 = 0.02. Similar discrepancies are observed for other parameters. (d) The error (t) = ρ exact s (t) − ρ add s (t) of the additive ME using λ 1 = λ 2 = λ and γ 1 = γ 2 is plotted. As λ increases, memory times are reduced and the additive ME result starts to approach the exact solution in the Markovian limit.
case of a single photonic bath (say g 1,k = 0, but g 2,k = 0), this model can be solved using an exact ME that is applicable in both Markovian and non-Markovian regions [37].
Since the two baths share a similar structure, one would expect the decoherence effects to be similar and just add in the ME. If we intuitively assume that the presence of a second independent photonic bath simply contributes an additional term to the ME, the additive ME will take the forṁ where the time-dependent decay rate i (t) due to bath i originates from the non-Markovian properties of the baths [13,37]. We note that a non-Markovian ME in general has an integrodifferential structure. The single bath case here is a specific example where the ME can take a time-convolutionless form [13]. Figure 1(b) shows the evolution of the excited state population ρ ee (t) of an initially excited TLS in the presence of two independent photonic baths evaluated by this additive ME [Eq. (3)]. We consider Lorentzian spectral densities J i (ω) = where λ i can be interpreted as the inverse memory time of the photonic bath i. This model can be mapped to the CQED, where λ i and √ γ i λ i /2 play the roles of cavity decay rate and vacuum Rabi frequency, respectively. Using this spectral density the time-dependent decay rate becomes [13,37]. In the absence of the second bath (γ 2 = 0), the TLS displays the expected damped vacuum Rabi oscillation. This correct behavior attributes to the exact and time-dependent i (t) in the ME when there is only one bath. When we add another photonic bath (γ 2 = 0), the additive ME solution leads to a severely suppressed dynamics of the TLS.
The reduction of the coherence of the additive ME result is caused by the additive assumption. This can be shown by comparing with the exact solution of the Schrödinger equation of this model system that includes the entanglement dynamics between the system and the environments. The total systemenvironment wavefunction reads where |1 k ,0 denotes the state of a photon in bath 1 and vacuum in bath 2. Figure 1(c) shows the exact result using the same physical parameters (with d = 0). Apparently, the coherence of the vacuum Rabi oscillation is maintained despite the existence of an extra bath. For baths with similar spectral widths λ i , the addition of a second bath simply modifies the vacuum Rabi frequency, rather than the decay rate. These non-Markovian coherent features are missing in the additive ME solution. Unlike previous works [38,39] that showed the failure of the additive ME based on a small parameter expansion of the superoperator, we emphasize that the non-Markovian decay rates in Eq. (3) are nonperturbative and exact in the single bath situation.
It is interesting to note that the discrepancies between the additive ME and the exact solution start to fade when both baths become Markovian. Figure 1(d) gives the difference between these two results as a function of the spectral width λ = λ 1 = λ 2 . As λ increases, the spectral densities are flattened and the memory times are shortened. The revival behavior of the TLS is replaced by an exponential decay, and in this Markovian limit, the system can be well approximated by an additive ME.

B. Nonadditive ME
Within this model, we have just demonstrated the deviation of the additive ME solution that assumes no bath interference from the exact solution that contains the full system-bath entanglement. We now show how the exact result can lead to a nonadditive ME that contains explicit interference terms between the two bath parameters. Recalling the total systemenvironment wavefunction in Eq. (4), one can show that the coefficients are governed by the following exact differential equations: where k i (t) = k g 2 i,k e i(ω 0 −ω k )t gives the memory kernel due to bath i and is related to the spectral density through k i (t) = dωe i(ω 0 −ω)t J i (ω). Following the procedures in Refs. [13,37,38], one can obtain a ME by tracing out the bath degrees of freedom from the total entangled wavefunction. The TLS oscillates between the excited and ground states with a non-Markovian decay rate 12 (t). In the additive ME, the overall decay rate is approximated by an incoherent sum exact 12 (t) ≈ additive The difference between exact 12 (t) (solid blue line) and additive Considering c e (0) to be real, an exact and nonadditive ME for this model system iṡ where the non-Markovian and time-dependent decay rate exact Using the Lorentzian spectral densities J i (ω) used in the previous section, the memory kernels become k i (t) = γ i λ i e −λ i |t| /2 and Eq. (5) can be solved analytically. The resultant decay rate becomes exact where s i are the three roots of the equation s(s + λ 1 )(s + λ 2 ) + γ 1 λ 1 (s + λ 2 )/2 + γ 2 λ 2 (s + λ 1 )/2 = 0. This exact decay rate clearly displays a mixing of physical parameters from the two baths. When the decay rate of either one of the baths vanishes, exact 12 reduces to the single bath situation as we described before: where d 1 = λ 2 1 − 2γ 1 λ 1 and a similar expression for 2 (t) can be obtained in the same way. Now, we can see that the additive ME [i.e., Eq. (3)] approximates the overall decay rate by a sum of the decay rates of the two baths and thus misses the interference effect additive 12 Figure 2(b) shows the difference between exact 12 (t) and additive 12 (t). Since the density matrix is oscillating, the non-Markovian decay rate can in general be negative and divergent. Here exact 12 (t) is periodic in time, implying a coherent oscillation between the excited state and the ground states. On the other hand, additive 12 (t) loses this feature because it originates from a sum of incoherent decays as shown in Fig. 2(a).

A. Diagrammatic perspective
To trace the origin of the issue of the additive ME, a microscopic theory that possesses the entanglement between the system and the environment is required. Here, we adopt the diagrammatic approach based on the Keldysh nonequilibrium Green's function technique [35,[40][41][42]. Different from the conventional ME approach that traces away the environmental degree of freedom, the diagrammatic technique has the advantage that maintains the full dynamics of the entangled systemenvironment wavefunction. For a general system interacting with two environments through V (t) = V 1 (t) + V 2 (t), the dynamics of a system observable O S can be formally expressed as where |ψ SE (0) is the initial system-environment state and the generalization to multiple reservoirs is straightforward. A diagrammatic expansion of Eq. (10) allows the visualization of the underlying physical processes and the justification of the additive ME.
To understand the correlated decoherence effect of the two baths, we depict the lowest-order interference diagrams in Figs. 3(a) to 3(c). These diagrams are extracted from the expression To simplify the physical picture, we have considered a separable initial state |ψ SE (0) = |ψ S (0) ⊗ |ψ E 1 (0) ⊗ |ψ E 2 (0) , and assumed ψ E i (0)|V i |ψ E i (0) = 0. Each diagram represents how the physical state |ψ SE evolves from the time 0 (left) to t (right) while undergoing system-bath interactions (dashed or dotted lines). It is clear that the two baths, though being independent, can affect the quantum system with finite time overlaps. The problem of the additive ME is ascribed to the missing of this kind of correlated influence from both baths.
For comparison purposes, we can perform a similar diagrammatic expansion of the additive ME governed by Eq. (1). This can be achieved by finding the perturbative solution of the reduced density matrix from the additive ME using V t as a small parameter. The lowest-order results are shown in the right-hand side of Figs. 3(a) to 3(c). We find that while the additive ME retains the nonoverlapping diagram [ Fig. 3(a)], it forbids the overlapping effect of the two baths [ Fig. 3(b)] and excludes all crossing terms [ Fig. 3(c)]. In general, the additive ME assumption corresponds to a ladder-like diagrammatic structure as shown in Fig. 3(d). The missing diagrams are of the same order as the noncrossing counterparts and should be included to recover the true dynamics of the reduced physical system. This analysis can be easily extended to the higher orders. The inclusion of higher-order terms in the expansion of eachL i [ρ s ] cannot compensate the missing interference diagrams.
It is known that the presence of an initial systemenvironment correlation can modify the ME [13,[43][44][45]. We note, however, that our correlation effect emerges dynamically and is present even when the system and environments are initially separable.
The disappearance of this interference bath effect in the Markovian limit can be understood with the aid of diagrams as well. For two Markovian baths with zero memory times, i.e., V 1(2) (t i )V 1(2) (t j ) ∼ δ(t i − t j ), the two ends of each dashed (dotted) line would collapse to the same point in the diagram. In this regime, all the crossing diagrams vanish and the diagrams are reduced to have a ladder structure. The baths can no longer affect the system with time overlaps and the additive ME becomes a good approximation. We stress that this situation requires all the baths being Markovian. The correlation effect still exists in the case of one non-Markovian bath plus multiple Markovian reservoirs, as we will see in Sec. IV.

B. Projection operator ME framework
Another theoretical framework for an open quantum system with non-Markovian environments can be carried out using the Nakajima-Zwanzig projection operator technique [38]. This would help us to have a parallel understanding of how the bath interference terms arise in the corresponding nonadditive ME, compared to the diagrammatic treatment above.
The idea of this approach is to start with an exact Liouville equation that describes the dynamics of the total system-environment density matrix ρ SE (t), then project the state on relevant and irrelevant parts, and trace out the environment at the end to obtain a formal ME. The derivation procedure is detailed in Ref. [13] and we shall not repeat it here. Defining Lρ = −i[V ,ρ] and the projection operator Pρ SE (t) = Tr E {ρ SE (t)} ⊗ ρ E (0), the resultant ME is given by where we again assumed a separable initial state and Tr E {V } = P LP = 0. The memory kernel superoperatorK(t,t ) can be expressed aŝ Equation (12) is formally exact, but hard to implement in practice due to the complicated memory kernel. Yet we can still observe how the bath interference terms appear from this formal ME treatment. Taking V = V 1 + V 2 and expanding Eq. (12) in the powers of V , we find an nonadditive MĖ The first terms correspond to the Born (without Markov) approximation. The second terms include interference (∼V 2 1 V 2 2 ) and noninterference (∼V 4 1 or V 4 2 ) contributions. These interference terms share the same physical origin as those in Eq. (11) obtained from the diagrammatic calculation.
When both environments are Markovian, i.e., , it is straightforward to show that the interference terms would vanish; while when one of the baths is non-Markovian, the interference remains, in agreement with the diagrammatic result. In general, it is difficult to solve Eq. (14) due to its multi-integro-differential structure. It would be more direct to calculate the perturbative solution as we will see in the example of the following section.

A. Dephasing of a system
We have just shown how the bath interference can affect the dynamics of the reduced system and how it emerges dynamically due to the system-environment entanglement. We now turn to an explicit example in which the inclusion of interferences between independent baths could yield a reduced decoherence of the physical system. Consider a driven three-level system under the influence of Markovian atomic relaxations and non-Markovian dephasing as shown in Fig. 4(a). The resonantly driven Hamiltonian and the system-bath interactions where σ ij = |i j |, σ z = σ 22 − σ 11 , and a k and b k are bosonic bath operators. The electric-dipole couplings g i,k are responsible for the corresponding population relaxations, and we adopt a bosonic dephasing model [46] with coupling constant λ k . In the dressed state picture, the optical drive together with the atomic decays i pump the system to the dark state |D = (|1 − |2 )/ √ 2; meanwhile, the dephasing term would decohere it to a classical mixture of dark and bright states |B = (|1 + |2 )/ √ 2. In the following, we shall examine the dark state fidelity as a consequence of these two competing dissipative processes with and without the presence of bath interferences.
The Markovian relaxation and non-Markovian dephasing rates are defined through where K(τ ) is the dephasing memory function. To the second order of γ ⊥ , an additive ME for this model is given bẏ where On the other hand, we can obtain the full equation of motions for the system and bath operators. Upon eliminating the bath operators, we arrive at an exact and reduced Heisenberg equation of motion where (17), the presence of the two time correlation functions in this exact quantum equation allows the interference between the relaxations and the dephasing processes. These two approaches are generally different except in the Markovian dephasing limit such that K(τ ) = δ(τ )/γ ⊥ . In the following, we are going to compare the additive ME [Eq. (17)] and the quantum solution [Eq. (18)]. We solve Eq. (18) perturbatively up to O(γ 2 ⊥ ) and compare it with the solution of Eq. (17), which is valid to the same order of γ ⊥ . This corresponds to a perturbative solution of the problem that is second order in the non-Markovian dephasing interaction (V 2 ), but infinite order in the Markovian relaxation interaction (V 1 ). Figures 4(d) and 4(e) provide the time and field dependences of the dark state population computed by solving the Heisenberg equations (quantum), the additive ME, and the case of Markovian dephasing. We consider a long memory time such that K(τ ) ≈ 1 for the non-Markovian calculations in the time regime of interest. The quantum solution demonstrates a much-suppressed dephasing than that of the additive ME, which behaves similarly to the Markovian result. For 1 ≈ 2 , the dark state fidelities are ∼O[(γ ⊥ / 1 ) 2 ] and ∼O[γ ⊥ / 1 ] for the quantum and additive ME solutions, respectively. Together with the TLS example in Fig. 1, the result of this dephasing example indicates the problem of the additive assumption in both the long and short time scales.
The diagrammatic representation offers a visualization of the difference of the underlying mechanisms between these two approaches. Figure 4(b) shows the virtual dark state evolution. To the second order of γ ⊥ , the dephasing interaction flips the initial state from |D to |B , which is then pumped back to |D , and this procedure repeats the second time. In this picture, the relaxations and dephasing can act simultaneously on the system. However, in the additive ME that neglects such an interference, the pumping is forbidden during the virtual dephasing process as presented in Fig. 4(c). This discrepancy eventually leads to a reduced dark state fidelity of the additive ME solution. A more detailed diagrammatic representation based on the reduced density matrix and another additive ME that includes the dressing effect of the drive can be found in the Appendix.

B. Applicability for general quantum systems
The dephasing problem above provides a particular scenario where the interference between environments is important and can help the reduction of decoherence. It is difficult to carry out the full ME, like Eq. (14), to analyze such an effect given this relatively simple example. For general and more complex physical circumstances, the bath interference can either assist or deteriorate the coherence of the physical system. The diagrammatic procedure can be straightforwardly applied to study this interference effect in open quantum systems that involve at least one non-Markovian bath. This would include many-body systems or problems with multilevel structures, where an exact solution is not available and it is difficult to carry out a full ME simulation. For simple few-level systems coupled to Markovian baths only, the bath interference is negligible and the conventional additive ME approach would be a simpler theoretical tool to calculate the corresponding decoherence dynamics.
For nonzero temperature reservoirs, since the interference effect relies on the existence of a non-Markovian environment, it would fade when the non-Markovianity is blurred by the high temperature effect [47].

V. CONCLUSION
We show that the additive assumption of independent environments in the ME is valid only when all the baths are Markovian. Because of the dynamical entanglement between the system and the environment, quantum interferences between multiple environments can emerge and strongly alter the decoherence mechanism for non-Markovian circumstances. A diagrammatic technique allows a direct justification of the ad-ditive assumption. We also provide a realistic example in which the inclusion of this interference effect could lead to a lower decoherence of the physical system. We believe the important role of the system-environment entanglement provides a new direction to understand and suppress decoherence for general open quantum systems.

ACKNOWLEDGMENTS
We thank Tony Lee and Shenshen Wang for very helpful discussions. This work is supported by NSF, CUA, ITAMP, AFOSR, and MURI. C.K.C. is supported by the Croucher Foundation.

APPENDIX A: ALTERNATIVE ADDITIVE ME FOR THE DEPHASING PROBLEM
In the main text, we describe one possible additive ME for the driven dephasing system. The equation is rewritten here: Not only does this equation assume the relaxations ( i ) and dephasing (γ ⊥ ) processes are additive, but also treats the driving term in the same manner. An alternative way is to include the dressing effect of the optical drive on the dissipative processes. We will study this alternative in the following. To obtain this alternative additive ME, we work in the interaction picture (by treating H as the bare Hamiltonian) such that the dissipative operators are dressed by the driving Hamiltonian, and then perform the standard Born approximation for the dephasing interaction. The resultant additive ME (we denote as "add2") iṡ ρ add2 FIG. 5. (Color online) A comparison of the dark state populations of the driven system using various approaches. (a) The quantum solution, (b) the additive ME that neglects the bath interference, and (c) the alternative additive ME that includes the dressing of the optical drive, but still omits the bath interference. Parameters used: where the driving propagator U (t) = e −iH t . Note that U (t) only enters through the non-Markovian dephasing memory kernel rather than the Markovian relaxation terms, and this equation still contains no interference between the two dissipative processes.
Equations (A1) and (A2) can be solved directly by the Laplace transformation technique. Figure 5 presents typical time evolutions obtained from these two additive ME approaches and the quantum method that includes the bath interference [i.e., Eq. (18)]. We see that this alternative solution shows an unphysical result [Fig. 5(c)] that it does not display the dephasing and it violates the positivity of the reduced density matrix. The violation of the positivity is not uncommon for non-Markovian MEs [13]. In this case, the failure is attributed to the fact that Eq. (A2) only includes the interference between the drive and the dephasing, but still omits that between the two dissipative processes. During the virtual dephasing process of the system, this equation allows the influence of the optical drive, but not the relaxation, resulting in an oscillatory dark state population. Therefore, this alternative approach still gives an incorrect dynamics of the reduced density matrix. The diagrammatic construction in the next section provides a better visualization of the underlying virtual processes. interference is neglected. They both lead to incorrect virtual states during the dephasing processes, regardless as to whether we include the dressing of the optical drive or not. On the other hand, the diagrams in the Markovian dephasing limit [i.e., K(τ ) = δ(τ )/γ ⊥ ] given in Fig. 6(d) shows some resemblance to those of the additive ME without the virtual states in Fig. 6(b). This explains why these two approaches have similar dark state populations as shown in Figs. 4(d) and 4(e).