Emergent $\mathcal{PT}$ symmetry in a double-quantum-dot circuit QED set-up

Open classical and quantum systems with effective parity-time ($\mathcal{PT}$) symmetry, over the past five years, have shown tremendous promise for advances in lasers, sensing, and non-reciprocal devices. And yet, the microscopic origin of such effective, non-Hermitian models is not well understood. Here, we show that a non-Hermitian Hamiltonian emerges naturally in a double-quantum-dot-circuit-QED (DQD-circuit QED) set-up, which can be controllably tuned to the $\mathcal{PT}$-symmetric point. This effective Hamiltonian, derived from a microscopic model for the set-up, governs the dynamics of two coupled circuit-QED cavities with a voltage-biased DQD in one of them. Our analysis also reveals the effect of quantum fluctuations on the $\mathcal{PT}$ symmetric system. The $\mathcal{PT}$-transition is, then, observed both in the dynamics of cavity observables as well as via an input-output experiment. As a simple application of the $\mathcal{PT}$-transition in this set-up, we show that loss-induced enhancement of amplification and lasing can be observed in the coupled cavities. Our results pave the way for an on-chip realization of a potentially scalable non-Hermitian system with a gain medium in quantum regime, as well as its potential applications for quantum technology.


I. INTRODUCTION
For an isolated (quantum) system, the Hamiltonian is the generator of its time evolution. A fundamental postulate of the quantum theory is that this Hamiltonian is Hermitian. It ensures real energy eigenvalues, and a unitary time evolution for the system. This changed two decades ago, when Bender and co-workers discovered a large class of non-Hermitian, continuum Hamiltonians on an infinite line with purely real spectra [1]. The key feature shared by all of these Hamiltonians is an antilinear symmetry, i.e. they are invariant under the combined operations of parity (P) and time reversal (T ). This antilinear symmetry, [PT , H] = 0, ensures that the eigenvalues of the Hamiltonian H are either purely real or occur in complex conjugate pairs [2,3]. When the eigenvalue λ is real, the eigenstate |ψ λ is also a simultaneous eigenstate of the PT operator with eigenvalue one. When λ is complex, the PT operator maps the eigenstate onto the eigenstate for its complex conjugate λ * , i.e. PT |ψ λ = |ψ * λ . Therefore, a PT -symmetric Hamiltonian spectrum shows a transition, known as the PT -symmetry breaking transition, when it changes from purely real to complex. In addition to the eigenvalues, corresponding eigenstates also become degenerate at the PT -transition point, and therefore the non-Hermitian Hamiltonian is defective at this exceptional point [4,5]. Over the years, it has become clear that PT -symmetric Hamiltonians faithfully model open, classical, zero-temperature systems with balanced, spatially separated gain and loss; the PT -symmetric phase with real eigenvalues thus represents an open, quasi-equilibrium system, whereas the PT -broken phase with amplifying and decaying eigenvectors represents a system far removed from equilibrium.
The subject of non-Hermitian, PT -symmetric Hamiltonians and their exceptional point degeneracies has evolved into a rich and active field. In the classical domain, effective non-Hermitian systems have been theoretically and experimentally studied in waveguides [6,7], fiber loops [8], resonators [9,10], electrical circuits [11][12][13][14][15], mechanical oscillators [16,17], viscous fluids [18], magnonics [19][20][21], acoustics [22][23][24], optomechanics [25,26], optical lattices [27]. Occurrence of exceptional points has been used in device applications like sensing [23,24,[28][29][30], single-mode lasing [31], unidirectional invisibility [32], loss-induced transparency [33], loss-induced lasing [34] etc (see recent reviews [35][36][37][38][39] for further details and references). A wide variety of physically motivated one-dimensional (for example, [40][41][42][43][44][45][46][47][48]), and two-dimensional (for example, [49][50][51][52][53]) condensedmatter lattice models have also been studied. From a fundamental perspective, exploration of topological aspects of non-Hermitian systems has bee gaining ground both theoretically (for example, [25,40,44,[52][53][54][55][56]) and experimentally [13,27,[57][58][59]. In the quantum domain, realization of systems governed by effective non-Hermitian Hamiltonians and exploration of the effect of exceptional points has only been possible very recently. Quantum non-Hermitian systems have been experimentally realized in linear optical circuits [59,60], quantum photonics [61], ultracold atoms [62,63], a diamond NV center [64], superconducting circuits [65,66], atomlight interacting systems [67], a lossy quantum-pointcontact [68]. However, all the realizations in the quantum regime are in cases where the overall system is dissipative. In particular, none of the realizations feature a gain medium. To our knowledge, transitions across exceptional points for balanced-gain-loss systems, featuring a The equivalent circuit diagram. The set-up consists of two cavities each of frequency ω0 and losses κ1,2, connected to each other with coupling λ. This is shown in (b) as two resonators of capacitance Cc and inductance Lc, which are inductively coupled to each other via L λ and capacitively coupled to sources of loss via C1 and C2. Additionally a DQD is located in the left cavity. The DQD is modelled by two fermionic sites (charge islands LD and RD in the circuit diagram) with a detuning ε between the on-site energies (controlled by the gate voltages VLP , VRP , VL and VR), a hopping between the sites with strength tc (controlled by the gate-voltage VC ) and a strong repulsive interaction of strength V between the sites. It is coupled to source (S) and drain (D) fermionic leads with coupling Γ (controlled by VLP , VRP ) via which a source-drain voltage bias of µ1 − µ2 is applied across it. The DQD is coupled with the left cavity via dipole coupling of strength g0 (shown in the circuit diagram by capacitive coupling Cg). To model experimental conditions accurately, we also consider the coupling of the DQD to the substrate phonons of spectral density J ph . The voltage-biased DQD is configured to be population inverted, making it an effective gain medium [69][70][71][72]. We show that, ε and tc can be tuned so that the net gain in the left cavity is equal to the loss in the right cavity. Under this condition, the dynamics of the two cavities is governed by an effective PT symmetric Hamiltonian. We further explore signatures of exceptional points on tuning the coupling between cavities λ and the loss at the right cavity κ2. This set-up can be realized in state-of-the-art circuit QED platforms [73][74][75][76][77]. In-situ tuning of coupling between cQED resonators, as well as of the resonator losses have been reported previously in experiments [66,[78][79][80][81].
gain medium in the quantum regime have not been realized yet, although there has been a theoretical proposal [82].
The dynamics generated by a non-Hermitian Hamiltonian is not unitary, and such dynamics can result from an open quantum system coupled to one or more environments (baths). Non-Hermitian Hamiltonians arising out of phenomenological Lindblad equations has been theoretically explored in recent works [83][84][85][86]. However, a complete microscopic theory of an open quantum system showing the emergence of a non-Hermitian Hamiltonian with a PT -transition has not yet been explored. In microscopic theories of open quantum systems, the full setup of system+baths is taken to be an isolated system described by a Hermitian Hamiltonian. Then the evolution equation for the system degrees of freedom is obtained by integrating out the bath degrees of freedom. The resulting equation describes the non-unitary dynamics of the system. In the usual cases with a thermal bath, such an analysis leads to a dissipative (lossy) system, but, by clever bath-engineering, it can also lead to a gain [70-72, 76, 87-92]. However, due to the unitary nature of the underlying dynamics of the full set-up, each bath, in addition to providing loss or gain, generates a noise that is consistent with the fluctuation-dissipation theorem [15,[93][94][95]. In the low-temperature limit, the noise generated by the lossy bath can be ignored depending on system parameters and time scales, but due to the quantum limits on linear amplifiers [93,94], a gain-bath leads to noise down to zero temperature. Therefore, whether truly balanced-gain-loss PT -symmetric quantum systems are possible, at the field-operator level [96], at the level of bosonic-field expectation values, or at the level of expectation values of bilinears of bosonic fields, is an open question. In particular, a fundamental analysis of a quantum system with gain and loss at low temperatures must take into account quantum fluctuations.
In this work, starting from a completely microscopic description, we theoretically show that a state-of-the-art open quantum system can be tuned to observe the PTtransition, and study the signatures of quantum fluctuations due to the baths. The set-up we consider consists of a solid state double-quantum-dot (DQD) connected to two coupled circuit-QED (cQED) cavities. In recent years, a DQD in a cQED cavity has been well characterized both theoretically and experimentally [74-76, 87, 90, 91, 97-99]. It is known that when the DQD is voltage-biased via electronic leads under suitable con-ditions, it can be population inverted and can act as a widely controllable gain medium [69,71,91,92]. This has led to the realization of an on-chip laser in the microwave regime [70,72].
In our analysis, the DQD is modeled by two fermionic sites, with hopping between them and strong repulsive interactions. The hopping strength and the on-site energies of the fermionic sites can be tuned widely in state-ofthe-art realizations. Each fermionic site is also connected to an external fermionic lead (bath), which allows us to apply a voltage bias across the DQD. Further, to model the DQD accurately, one should consider its coupling to the phonons in the substrate. We show here that a set-up of two coupled cavities with a DQD in one of them can be tuned so that time-evolution of the cavity operators is governed by an effective PT -symmetric Hamiltonian. This happens when the effective gain from the DQD balances the cavity losses. However, the quantum fluctuations of the DQD leads to an additional noise term in the effective equations of motion for the two coupled cavities, which we derive from the microscopic model of the entire set-up. Thus, in this work, we propose an on-chip realization of a quantum PT -symmetric system in which the effects of quantum fluctuations of the gain medium can be controllably studied.
Our proposal is very different from a recent proposal for the realization of a PT -symmetric system in the cQED system by Quijandría et al. [82]. The latter considered two cQED cavities, each coupled to a qubit whose frequency is modulated via a coherent drive, and showed that the expectation values of the cavity field operators, in appropriate parameter regime, are governed by an effective PT -symmetric Hamiltonian. The DQD voltage bias in our model acts as an incoherent drive, and our analysis extends to equations of motion for the cavity field operators and their bilinear combinations. This allows us to study the effect of quantum fluctuations on the cavities, which was not explored in Ref. [82]. As a simple application of the PT -transition in our set-up, we further show the possibility of loss-induced enhancement of amplification and lasing. We further show loss-induced increase in average photon number in absence of any coherent drive in the cavities. This particular feature requires a gain medium with quantum fluctuations, and cannot be observed either in the classical systems with negligible fluctuations from the gain medium, or in the existing dissipative quantum system experiments without a gain medium discussed previously.
The paper is arranged as follows. In Sec. II, we describe the microscopic model of the set-up and give the parameter regime where it should be operated. In Sec. III we obtain the effective equations of motion for the field operators of the coupled cavities, and show how effective PT -symmetry can be obtained. In Sec. IV we explore the PT -symmetric dynamics of the coupled cavities and the effect of quantum fluctuations when the cavity losses are balanced by the gain from the DQD. In Sec. V we show that effects of a PT -transition can be observed even when the losses of the cavities are more than the gain from the DQD. We explore input-output experiments, effects of quantum fluctuations on them, as well as loss-induced lasing and enhancement of amplification in this section. In Sec. VI we summarize our main results and give the consequences of our work and future directions. Certain details of the analytical calculations are delegated to the Appendix.

II. MICROSCOPIC HAMILTONIAN FOR THE DQD COUPLED TO CIRCUIT-QED CAVITIES
The schematic of the set-up we consider, along with the circuit diagram, is given in Fig. 1. It consists of two, identical coupled cavities, with a DQD in the left one. The cavity with the DQD is an existing experimental set-up explored in a series of recent experiments [69-72, 76, 77, 98-101]. The parameters of the DQD are widely tunable in experiment and, under a voltage bias, the DQD can be population inverted, thereby making it a tunable gain medium for the cavity. The main idea of this paper follows from the observation that when coupled to another cavity, the gain from the DQD can be tuned, within current state-of-the-art experimental parameters, to exactly match the cavity losses, thereby realizing a PT -symmetric system.
We consider a completely microscopic model of the entire set-up. The main parts of the set-up consist of the DQD and the two cavities. These are described by the following Hamiltonianŝ where Θ(t) is the Heaviside function. Here,ĉ † 1,2 denote fermionic creation operators for sites 1 and 2 that model the DQD, ε is the energy difference between the two sites, t c is the hopping amplitude between the two sites, and V > 0 denotes the capacitive charging energy between the two sites.b 1,2 represent the bosonic creation operators for the two cavities, each with frequency ω 0 , that are coupled via a number-conserving hopping process with amplitude λ. Closely following the experimental set-ups, the DQD is dipole-coupled to the (first) cavity with strength g 0 when the cavity is switched on at time t = 0. Thus, Eq.(1) captures the microscopic model of the DQD and the two cavities.
Each of the three main components (two cavities and one DQD) is connected to multiple baths. The two bosonic baths, one for each cavity, are described by the Hamiltonian whereB † s1 (B † s2 ) is the bosonic creation operator of the sth mode of the bath attached to the first (second) cavity, and Ω s1 (Ω s2 ) is the energy that mode. The baths are coupled to the respective cavities aŝ where κ s is the coupling of the sth mode of the bath attached to the th cavity. Each fermionic site of the DQD is connected to its respective lead, are described by the Hamiltonian whereâ † s1 (â † s2 ) is the fermionic creation operator of the sth mode of the source (drain) lead attached to the first (second) dot, and E s1 (E s2 ) is the energy of that mode. The fermionic leads are bilinearly coupled to the DQD viaĤ where Γ s is its coupling of the sth mode of the fermionic lead attached to the DQD site . To model experimental set-up consistently, we also have to consider that the DQD is dipole-coupled to a phononic bath in the substrate on which it is located [98]. The phonon bath and the coupling Hamiltonians are given bŷ whereB ph s is the phonon annihilation operator for the sth mode of the phononic bath, and λ ph s is its coupling to the DQD dipole operator. While this piece of the microscopic Hamiltonian is not necessary for the physics that we will discuss, it is unavoidable, and relevant, in some of the state-of-the-art experimental set-ups [69-72, 98, 99, 101].
The two bosonic baths for the two cavities, the two fermionic leads, and the phonon bath are characterized by their frequency-dependent spectral functions, Here J (ω) is the spectral function of the bosonic bath coupled to th cavity, J f (ω) is the spectral function of the fermionic lead coupled to th site in the DQD, and J ph (ω) is the spectral function of the phononic bath. For simplicity, we consider both fermionic leads to be in the "wide-bath limit"; that gives rise to a constant spectral function denoted by Γ. The spectral functions for the bosonic baths of the two cavities give rise to individual decay rates κ . The effect of the substrate phonons on the DQD, as well as the spectral function for the phonons can vary depending on the platform. They have been well-characterized in literature both theoretically [98,102,103] and experimentally [99,104]. The spectral function J ph (ω) in Eq. (7) is known to be a good description for gate-defined DQDs on a GaAs substrate [98,102]. The frequency ω c is given by ω c = c n /d, where c n is the speed of sound in GaAs and d is the distance between the two quantum dots. The frequency ω max is the upper cut-off that provides spectral damping at frequencies much higher than the repetition rate of phonon travel between the two dots. We take ω max = 10ω c . For a DQD with GaAs substrate, c n ≈ 3 × 10 4 m/s and d ≈ 150 nm, gives ω c = 20GHz. The dimensionless parameter γ b controls the coupling between the DQD and the phonons.
We will like to mention that, although not required, microscopic models giving rise to all the spectral functions given in Eq. (7) can be written down assuming finite but large hard-cutoffs in frequency [105]. The initial state of the whole set-up is the direct product of arbitrary states of the DQD and the cavities, and equilibrium thermal states of the baths with their respective temperatures and chemical potentials. We consider all the baths at the same inverse temperature β, but the chemical potentials for the two fermionic leads are given by µ 1 = µ 2 (see Appendix A). This creates a voltage bias across the DQD, thereby driving the DQD to a non-equilibrium steady state. We consider that the connection between the DQD, the fermionic leads and the phononic substrate is switched on at a time t = t 0 0. On the other hand, the two cavities are connected to the DQD and the bosonic baths at time t = 0. Further, as we discuss below, we make assumptions on the various energy scales appearing in the full set-up.
The DQD Hamiltonian is diagonalized by rotating to the a new basiŝ where σ y is the Pauli y-matrix. In the rotated basis, the DQD Hamiltonian becomeŝ where ω q = ε 2 + 4t 2 c , andN =Â † Â ( = 1, 2). We look at the parameter regime that ensures a resonant DQD at low temperature, a weak cavity-DQD couplng, i.e.
where n photons is the average number of photons in the cavity coupled to the DQD, and In experiments, the widely controllable parameters of the DQD are ε and t c , which thereby allow us to tune ω q and θ. Thus, for a ω 0 , the resonant condition ω 0 = ω q can be satisfied by tuning ε and t c . Under this condition, since the DQD-cavity coupling is weak, we can simplify the dipole coupling Hamiltonian through rotating wave approximation tô where the effective coupling is given by It is possible to tune θ while maintaining the DQD-cavity resonance condition ω 0 = ω q by varying ε and t c along the ellipse defined by ω 2 0 = ε 2 +4t 2 c . As we will see below, this freedom allows us to tune the system such that the effective dynamics of the two cavities is described by a PT -symmetric Hamiltonian.

III. EFFECTIVE HAMILTONIAN AND EMERGENT PT SYMMETRY
We can obtain effective equations of motion for the bosonic cavity operators by integrating out the rest of the set-up (see Appendix A). When we are in the appropriate parameter regime (see Eqs. (11), (12)), this leads to where, for future reference, we have additionally added a coherent drive of strength E 0 and frequency ω d in the cavity without the DQD. The 2×2 non-Hermitian Hamiltonian is given by where δ < 1 by the choice of our parameters and ∆N ss = N 1 ss − N 2 ss , is the steady-state population inversion in the DQD in the absence of cavities. The operator ξ A (t) is embodies the noise resulting from the presence of the DQD, in accordance with fluctuation-dissipation theorem. It has a Lorenzian power spectrum, is the steady state coherence of the DQD in absence of the two cavities. Note that the properties of the noise operator are not phenomenological, but are microscopically derived from our model. Under our assumptions on energy scales (Eqs. (11), (12) In this condition, it can be checked from Eq.(16) that ξ A (t) has a negligible effect on the dynamics of the cavities as long as So, assuming this, we can set In all of the equations above, the expectation value rep- where ρ tot is the density matrix for the initial state of the whole set-up (see Appendix A). We have also numerically checked that Eq.(20) holds in our chosen parameter regime. On the other hand, as we will see, ξ † A (t)ξ A (t ) , which embodies the quantum fluctuations due to the DQD, cannot be neglected. The bosonic baths that provide dissipation for the two cavities also lead to thermal noise. However, it can be shown that at low temperatures, i.e, for βω 0 1, their contribution is negligible compared the fluctuations coming from the DQD (see Appendix A 5). So, we have ignored them in Eq. (15).
We note that H eff is non-Hermitian in two ways. First, its diagonal elements have imaginary parts that represent the coupling-to-the-bosonic bath losses κ and the gain Γδ for the first cavity which houses the DQD. Secondly, the population inversion in the DQD has suppressed the hopping from the second cavity to the first one, leading to real, asymmetric off-diagonal elements. The presence of two non-Hermiticities makes H eff different from the most commonly studied 2×2 PT -symmetric Hamiltonians. It is straightforward to obtain the quadratic characteristic equation for H eff and check that it is purely real if and only if Eq. (21) ensures that the eigenvalues of H eff are either purely real or occur in complex conjugate pairs, and thus endows H eff with some antilinear symmetry. With this constraint, the eigenvalues Λ ± of the effective, non-  Hermitian Hamiltonian, Eq.(16) are given by (for balanced gain-loss), with an exceptional point degeneracy occuring at When the gain and loss cavities are strongly coupled, λ ≥ λ EP , the system is in the PT -symmetric phase (real spectrum); a weak gain-loss coupling, λ < λ EP drives the system into the PT -broken region (complex conjugate spectrum). These properties are reminiscent of the standard PT -transition; however, H eff is not paritytime symmetric when parity exchanges the cavity labels 1 ↔ 2, and time-reversal is just the complex conjugation operation K. To obtain the corresponding antilinear PT operator for this model, we express H eff as a linear combination of identity and Pauli matrix components, It follows from Eq.(24) that the Hamiltonian can be cast into the traditional form with a rotation R x (φ) = exp(−iφσ x /2) where φ = arctan(λδ/2κ 2 ). In the transformed basis, the Hamiltonian The traditional Hamiltonian H is symmetric under interchange of cavity labels (parity) and complex conjugation operation K (time-reversal). Due to the antilinear nature of the latter, though, in the original basis, the PT operator becomes to It is straightforward to check that the antilinear operator in Eq. Eq. (21), we now show that the constraint can be satisfied. Since the steady-state population inversion ∆N ss (θ) depends on the DQD parameters, we recast Eq. (21) as where the θ-dependent quantities are only present on the left-hand side of the equation. We consider a realistic set of parameters for the system (Table I) However, the exceptional point degeneracy at is still present, and the transition across it is manifest by the two non-orthogonal eigenmodes of H eff having different decay (or amplification) rates at weak coupling λ < λ EP , and a common decay (or amplification) rate at strong coupling, λ > λ EP [12,62]. Thus, the key features of the PT symmetry breaking transition, including the coalescence of the eigenmodes, survives even if the constraint in Eq. (21) is not fulfilled [65]. Lastly, we discuss the limitations of our analysis. Our effective Hamiltonian H eff is obtained via a linearized theory that works only when the number of photons in the two cavities are not too large. Whenever at least one the eigenvalues of H ef f has a positive imaginary part, the linearized theory predicts an exponential growth for the photon number associated with that eigenmode. This can happen either because H eff is PT symmetric and λ > λ EP , or because the net gain resulting from the DQD in the first cavity exceeds the loss in the second cavity, i.e. 2Γδ > κ 1 + κ 2 . Either of these scenarios essentially point to breakdown of the linearized theory. Nevertheless, their identification is important because such behavior in the linearized theory points to onset of lasing in the actual experiment. Though not necessary, a PT -transition is thus intimately linked with onset or suppression of lasing in the two coupled cavities.
In the following, we will be looking at the dynamics of expectation values of the field operators or complex quadratures b (t) ( = 1, 2), and those of the cavityfield bilinears like cavity photon numbers n (t) = b † (t)b (t) , and the intercavity photon current The dynamics of these quantities can be obtained from Eq. (15). The formal solution of Eq.(15) is given by where I is the 2×2 identity matrix. The dynamics of the complex quadratures are obtained by taking expectation value of the above equation. The dynamics of the cavity bilinears are obtained from dynamics of the 2×2 equaltime correlation matrix, whose elements are given by The expression for C(t) is obtained by taking the transpose of Eq. (30) and multiplying on the left by its Hermitian conjugate. Assuming that the initial state of the cavities is a coherent state, which satisfies, b † (0)b (0) = b † (0) b (0) , the expression for the correlation matrix can be written as, , and M(t , t ) is a 2×2 matrix whose only non-zero element is given by the Fourier transform of the noise power spectrum, In absence of quantum fluctuations, we would have C(t) = C (cl) (t), which would be consistent with the classical predictions. Thus the connected part of the correlation functions, given by C(t) − C (cl) (t) embody the quantum fluctuations. We see from Eq. (18) that P (ω) ∝ N 1 ss , which is the steady state occupation of the higher energy mode of the DQD. If the DQD is configured to act as a gain medium, it must be population inverted, which means, N 1 ss ∼ O(1). Thus, the quantum fluctuations cannot be neglected and, as we will see below, have non-trivial effects on dynamics of cavity-field bilinears.

IV. DYNAMICS OF THE BALANCED GAIN-LOSS CAVITIES
In this section, we look at the the dynamics of cavity observables when the set-up is tuned so that the Hamiltonian H ef f is PT -symmetric, i.e, Eq.(21) is satisfied. We ensure this by fixing along with the parameters in Table I. We assume that initially the second, lossy cavity is empty and the first, DQD-gain cavity is in a coherent state, We further assume in this section that there is no additional coherent drive, i.e, E 0 = 0. Satisfying Eq.(19) under these conditions require α g Â † 2Â 1 ss /ω 0 ∼ 10 −4 for our choice of parameters. Under these conditions, we calculate the dynamics complex quadratures, cavity photon numbers n (t) = b † (t)b (t) and the intercavity photon current The most remarkable effect of having PT symmetry, Eq.(26), for the Hamiltonian H eff is the possibility that an open system shows periodic dynamics when λ > λ EP , i.e, when the cavities are strongly coupled. Figure 3(a) shows the numerically obtained results for the amplitude of the complex quadratures with dashed lines. We get oscillatory dynamics with period T Λ = 2π/ (Λ + − Λ − ). Thus, they indeed show the signatures of a PT -symmetric phase. In contrast, when it comes to cavity-field bilinears, i.e. photon numbers and the intercavity current, the oscillatory behavior is superseded by the effects of quantum fluctuations at long times.
To see this surprising result, we note that, away from the exceptional point, the non-Hermitian Hamiltonian can be diagonalized by a similarity transformation. Let H T eff = RΛR −1 , where R has right eigenvectors of H T eff and Λ is a 2×2 diagonal matrix with real eigenvalues Λ ± . In this skewed basis, the elements of the correlation ma-trixC(t) = R −1 C(t)R can be explicitly calculated. For times t Γ −1 , we get wherem = R −1 (1 2 + σ z )R/2 . We see from Eqs.(36)-(37) that whileC 12 (t) oscillates with period T Λ , the diagonal elementsC (t) grow linearly with time due to the quantum fluctuations arising from the DQD. This linear growth in the eigenmode occupation numbers will eventually make the assumption g √ n photons ≤ O(Γ) (weak cavity-DQD coupling) invalid. So, starting in the PT -symmetric region, due to quantum fluctuations the system will be eventually driven out of the region of validity of a PT -symmetric effective Hamiltonian description. For our choice of parameters in Table I, g √ n photons ≤ O(Γ) corresponds to n photons < ∼ 50. To explicitly show the effect of quantum fluctuations, we plot n (t) −| b (t) | 2 (solid lines) in Fig. 3(a). While the amplitudes of the complex quadratures show perfect periodic oscillations of period T Λ , the quantum fluctuations show small oscillations about a linear growth, as discussed above. Fig. 3(b) shows dynamics of the photonic current I(t)/λ. This also shows a linear growth but with a much smaller slope, along with large oscillations of period T . These larger oscillations and weaker linear growth in the photon current is due to the larger weight of the oscillatory off-diagonal elementsC 12 (t), Eq.(37), and smaller weight of the linearly growing diagonal ele-mentsC (t), Eq.(36).  Table I and Eq. (34). The corresponding oscillation period is TΛ = 320 ns.
Next, let us see the dynamics of system at the exceptional point λ = λ EP . Since the Hamiltonian satisfies the characteristic equation (H T eff − λ EP 1 2 ) 2 = 0, at this point, the time-evolution operator Taylor expansion terminates at first order, i.e.
Therefore, with our initial conditions, the timedependent complex quadratures become Thus, the empty, lossy cavity has a quadratic growth, i.e. | b 2 (t) | 2 ∝ t 2 . In the gain cavity, | b 1 (t) | 2 remains flat for small times t 1/κ 2 , but switches to quadratic growth at times t 1/κ 2 . Similarly, it can be shown that at long times t 1/κ 2 , photon numbers in both cavities, n 1 (t) and n 2 (t) scale as t 3 . The numerically obtained results for the quadrature magnitude (squared) (dashed lines) and quantum fluctuations (solid lines) are shown in Fig. 4. The dotted line shows the quadratic fit to the empty-cavity | b 2 (t) | 2 , Eq.(40), while the dotdashed line shows the cubic fit for the quantum fluctuation results. We remind the reader that our model remains valid only at times when n photons < ∼ 50, and that determined the time-range chosen in Fig. 4.   Table I and Eq. (34), which leads to λEP = 1 MHz. The maximum time shown here corresponds to when there are ∼ 50 photons in the cavities, above which the linearized model cannot be applied for our choice of parameters. For the initially empty, lossy cavity | b2(t) | 2 grows quadratically with time, while the quantum fluctuations in quadrature amplitude grow as t 3 . For the cavity initially in coherent state, | b1(t) | 2 remains constant at times t 1/κ2, and tends to a quadratic growth as t 1/κ2. The quantum fluctuations in quadrature amplitude of this cavity also tend to a t 3 behavior at these times.
Lastly, we consider the system's behavior in the PTsymmetry broken region, i.e for λ < λ EP . Here, due to the presence of the amplifying eigenmode that is delocalized over the two sites, the photon numbers in both gain and loss cavities grow exponentially at small times, and provide a short-time cutoff beyond which the linearized theory fails. As mentioned before, identifying this regime is important because, such behavior in the linearized system points to onset of lasing in the actual experimental set-up.

V. SIGNATURES OF EXCEPTIONAL POINTS FOR DISSIPATIVE CAVITIES
In the previous section, we have seen the dynamics of the cavities, including the effects of quantum fluctuations, when the gain from the DQD-unit is balanced by the losses of the cavities, Eq. (21). In this section, we look at dissipative cavities, where the combined loss from the two cavities exceeds the DQD gain in the first cavity, i.e Even under this condition, this system can traverse across exceptional points on tuning either the coupling between the cavities, or one of the losses of one of the cavities. We show that this can be observed in particular input-output experiments through transmission, phase response and fluctuations, and may lead to interesting applications like loss-induced lasing. Note that, in-situ tuning of coupling between cQED resonators, as well as, of the resonator losses has already been demonstrated experimentally [66,[78][79][80][81].

A. Individually lossy cavities
First, we consider the case where each cavity is individually lossy, with same dissipation rate. In particular, we will assume, starting from a PT -symmetric condition, the dissipation rate in each cavity is increased by the same amount κ. Therefore, the effective non-Hermitian Hamiltonian in Eq.(24) changes to while Eq. (21) is still satisfied for κ 1 and κ 2 . The eigenvalues of the effective Hamiltonian now have an additive imaginary part Λ ± → Λ ± − iκ/2, where Λ ± are given in Eq. (22). However, the real parts of the eigenvalues remain the same. The overall dissipation will cause the system to reach a time-independent non-equilibrium steady state. The system will show a passive PT -transition as a function of λ, with λ EP still given by Eq. (23), where the real parts bifurcate [12,62,65]. In what follows, we show that, when κ = 2κ 2 , properties of this steady state capture the passive PT -transition. This is accomplished via a standard input-output experiment where a weak, coherent drive E 0 e −iω d t with amplitude E and frequency ω d is applied to one of the cavities and the steady-state transmitted signals are observed in either cavity. We will assume that the coherent drive is in the cavity without the DQD in it. We first go to the rotating frame with respect to the drive frequency,b rot (t)e iω d t =b (t). The transmitted signals T at the two cavities are given by the expectation value of the long-time solution of Eq.(30) scaled properly by the dissipation rates and E 0 , We write the transmitted signal at th cavity as T = |T |e iφ , where |T | is the amplitude of the transmitted signal, and φ is the phase response, φ ∈ [−π/2, π/2]. Both of the transmission amplitude and the phase response can be experimentally measured [76,77,98,100]. Explicit expressions for both the amplitude and the phase response can be obtained (see Appendix. C). When λ > λ EP and ω d = ω 0 ∓ (Λ + − Λ − ), we observe that the  Table. I.
phase response of the second cavity is given by Thus, when κ = 2κ 2 , the phase response becomes ±π/2. This means, there will be a phase change of ±π when ω d is tuned across values equal to ω 0 ∓ (Λ + − Λ − ). Further, it can also be checked that det(H eff − ω d 1 2 ) has a minimum at ω d = ω 0 ∓ (Λ + − Λ − ) under this condition. So, when κ = 2κ 2 , both the transmission amplitude and the phase response of the lossy cavity will accurately capture the bifurcation of real parts of eigenvalues of H eff across a passive PT transition. This is shown in Fig. 5, where both the transmission amplitude |T 2 | and the phase repsonse φ , detected at the lossy, second cavity as a function of ω d and λ are shown with color coding. The real parts of the eigenvalues are also plotted as function of λ for comparison. The passive PT transition and the location of the exceptional point are completely clear both in the transmission amplitude and in the phase response. The transmission amplitude peaks when the drive frequency ω d is equal to the real parts of the eigenvalues of H eff , while the phase response undergoes a change of π when ω d is changed across these values. Note that even without the gain medium (i.e, in absence of the DQD), the passive PT -transition in lossy cavities can be seen, through input-output experiments, if the two losses are different. In Eq. (16), this corresponds to setting δ = 0. However, the difference of that situa- tion with the one with a gain will be seen in the quantum fluctuations of the complex quadratures. As we will see below, the steady-state photon number will encode quantum fluctuations from the DQD-gain, thus distinguishing between a purely lossy set-up and the set-up with a gain DQD.
To see the effect of quantum fluctuations, as before, we will be interested in the steady-state values, n − | b | 2 ( = 1, 2). As can be seen from Eq.(32), the quantum fluctuations are independent of the presence of the coherent drive. So it suffices to look at the case where the coherent drive is absent. In that case, because the overall cavities are dissipative, expectation values of the field operators in steady state are zero. So, we look at the steady-state photon numbers n , and the photon current I as a function of λ. We numerically calculate these quantities using Eq. (32), and obtain the result at time t 1/κ. These results are shown in Fig. 6. At zero coupling, there are no steady-state photons in the second cavity and there is no photon current. However, due to the presence of the DQD in the left cavity there is a steady-state photon occupation in the first cavity when κ 1 + κ > 2Γδ. As the coupling is increased to values λ < λ EP , steady-state values photon number in the second cavity and the photon current both increase quadratically with λ. The quadratic growth slows down and the photon number in the DQD-cavity suppressed when λ → λ EP , For g > ∼ λ > ∼ λ EP , where g is the effective coupling between the DQD and the first cav- The vertical dashdotted green line shows κ2 = 2Γδ − κ1, which is balancedgain-loss condition. For κ2 > κ th 2 , the system gets an overall gain pointing to onset of loss-induced lasing. Parameters: λ = 2MHz, κ1 = 2MHz. Other parameters are same as in Table. I. ity, the photon numbers in both the cavities approach the same value, and photon current approaches maximum (Fig. 6a). As λ is further increased, the occupation of both cavities still remains almost the same, but they decrease with λ. The current also decreases correspondingly. This is because, in this regime, the relative strength of coupling to the DQD decreases. The suppression of the photon number in the DQD-cavity with increasing λ is shown in Fig. 6 b. It, too, shows a quadratic scaling that is expected from a perturbation theory in λ. We remind the reader that the nozero steady-state values are due to the presence of quantum fluctuations from the DQD. Without them, the dissipative cavities would be empty in the steady state.

B. Loss induced lasing and amplification
In previous subsection, we investigated the passive PT -transition as a function of intercavity coupling λ. A similar transition can occur if the loss in one of the cavities is increased, giving way to the counterintuitive phenomenon of loss induced lasing [34,106]. In this subsection, we show that this phenomenon can be observed in our set-up. We further argue that it is a consequence of having a PT -symmetric phase in this set-up.
To see this, we consider the loss κ 2 in the second cavity as the tunable parameter, whereas the rest of the set-up is kept fixed. Note that, transition across exceptional point Here, the input drive is at the right cavity. The PT -transition on increasing κ2 is approximately tracked by the amplitude response of the input-output experiment, but the phase response does not tack the transition. The real parts of the eigenvalues are also shown by dashed lines for comparison. The bottom panel (c) shows the average steady state photon numbers of the two cavities and the photon current as a function of κ2, in absence of any coherent drive. Parameters are same as in Fig. 7.
via tuning the loss in coupled cQED resonators has been recently shown experimentally [66], but in absence of any gain medium. For concreteness, in our set-up, let us consider that κ 2 is adiabatically increased starting from a balanced gain-loss PT -symmetric condition. On slight increase of κ 2 the system becomes dissipative and approaches a steady state in the long time limit. However, as κ 2 is increased further, from Eq.(29), we see that, an exceptional point is reached when For (28)), and hence, the imaginary parts of the eigenvalues of H ef f are same. But for κ 2 > κ EP 2 , (Λ + − Λ − ) is imaginary. In this case, the imaginary parts of the two eigenvalues of H ef f bifurcate into two different values. On further increasing κ 2 , the system remains dissipative as long as If κ 2 is increased beyond κ th 2 , imaginary part of one of the eigenvalues of H ef f become positive. Thus, starting from a condition when the overall system is dissipative, on increasing the loss of the right cavity beyond the threshold value, the system gets an overall effective gain, which points to onset of lasing. This is the rather counter-intuitive phenomenon of loss induced lasing that occurs in this set-up. It can be shown that this is a consequence of existence of a PT -symmetric phase with balanced-gain loss. To see this, we note that the exceptional point κ EP 2 given in Eq. (45) . Under balanced gain loss, κ 2 = 2Γδ − κ 1 , and this regime of λ corresponds to the PT -symmetric phase (see Eq. (23)). This shows that, it was necessary that we started adiabatically increasing κ 2 from the PT -symmetric phase of the balanced gain-loss system. Further, going back to Eq.(28), it is possible to see that there is another exceptional point possible with κ 2 as the only tuning parameter. This occurs when Since κ 2 , λ > 0, occurrence of this exceptional point requires κ 1 − 2Γδ > 0. Under this condition, both cavities are lossy, and hence, there can neither be lasing, nor a balanced gain-loss PTsymmetric phase. This can also be seen from Eq. (46), which shows that a necessary condition to see loss induced lasing in our set-up is κ 1 − 2Γδ < 0. Thus, the phenomenon of loss-induced lasing that can be observed in this set-up is a consequence of having a PT -symmetric phase with balanced gain-loss. (Note that, in other setups loss induced lasing can occur without any exceptional point [107].) This explicitly requires a gain medium and cannot be observed in the existing experiments in the quantum regime [59][60][61][62][63][64][65][66][67][68], none of which features a gain medium. The real and imaginary parts of the eigenvalues of H ef f as a function of κ 2 are plotted in Fig. 7. The positions of κ EP 2 , κ th 2 and κ 2 corresponding to the balanced gain-loss PT -symmetric condition are shown.
Our linearized theory allows us to explore the case κ 2 < κ th 2 . In this case, since the overall system is dissipative, there is a unique time-independent non-equilibrium steady state. As before, we look at an input-output experiment with a coherent drive at the second, lossy cavity. We look at the amplitude of the transmitted signal as well as the phase response at the lossy cavity. The results for the transmitted signal at the lossy cavity are shown in Fig. 8. Fig. 8(a) shows the transmission amplitude |T 2 |. The first point to note is that |T 2 | > 1. Thus, there is amplification of the transmitted signal in the lossy cavity. The peak in the transmission amplitude as a function of ω d and κ 2 approximately tracks the transition across the exceptional point. Moreover, the amplification of signal at the lossy cavity at resonant drive ω d = ω 0 increases with increase in the loss of the cavity. Thus, in this system, we have a novel microwave amplifier [69,91,92,108] with loss-induced enhancement of amplification facilitated by transition across an exceptional point. Fig. 8(b) shows the phase response at the lossy cavity. The phase response does not capture the transition across the expectional point. Nevertheless, the phase response does show some interesting features. There occurs a phase change of π at specific values of ω d , independent of κ 2 . These values are given by On the other hand, the phase response φ 2 is zero at values which depend on κ 2 . These values are given by where κ = (κ 2 + κ 1 − 2Γδ) /2. We will like to mention that even though loss induced suppression and revival of lasing has been experimentally observed in classical set-ups previously [34,106], to our knowledge, the phase response has not been previously measured. This measurement is experimentally possible in our DQD-cQED set-up (see, for example, [76,77,98,100]). To see the effect of quantum fluctuations, as before, we plot the steady state photon numbers of the two cavities, and the photon current in absence of any coherent drive in Fig. 8(c). They show non-monotonic behavior with κ 2 . Thus, the photon number of the two cavities and the photon current between the two cavities first decrease but then increase with increase in κ 2 , even in absence of any coherent drive. We remind the reader once again that, since there is overall dissipation, without the effect of quantum fluctuations from the DQD, the steady state photon numbers would be zero in absence of any coherent drive. Thus, the loss-induced increase in average photon number of the cavities in absence of any coherent drive is a feature possible only in a non-Hermitian system with a gain medium having quantum fluctuations. This feature of the non-equilibrium steady state thereby distinguishes our proposed set-up from previous experiments involving gain in classical systems [35][36][37][38][39], as well as, the existing experiments in quantum regime which do not feature a gain medium [59][60][61][62][63][64][65][66][67][68].

VI. SUMMARY AND OUTLOOK
In this paper, we have shown that a set-up of two coupled cQED cavities with a DQD in one of them can be used to explore the physics of PT -symmetry breaking transition, including the effect of quantum fluctuations. As has already been established in several experiments [69][70][71][72]98], a DQD in a cavity can act as a highly tunable gain medium for the cavity. In this work, under a reasonable approximation on various energy scales consistent with state-of-the-art experiments, we have microscopically derived an effective non-Hermitian Hamiltonian, along with the quantum noise term, which govern the dynamics of bosonic complex quadratures that describe the cavity. This derivation requires a linearization of the non-linear interaction between the DQD and the cavity and only holds if the number of photons in the cavity is not too large. We have also kept track of the quantum fluctuations coming from the presence of the DQD. We have thus obtained the signatures of the PT transition in the expectation values of the field operators (complex quadratures), and signatures of quantum fluctuations in the expectation values of field-operator bilinears (cavity photon numbers and intercavity current). In both PT -symmetric phase and at the exceptional point, we have found that quantum fluctuations lead to a linearin-time factor over the classical predictions for photon numbers and photon current. These signatures are then generalized to dissipative cavities that reach a steady state. Here, the passive PT transition is marked by the bifurcation of the decay rates for the two dissipative eigenmodes of the system. We have shown that this can be accurately tracked in a particular input-output experiment. We have further shown that loss-induced enhancement of amplification and lasing can be observed in the set-up as a consequence of a transition across an exceptional point. Lastly, we have identified loss-induced increase of average photon number in absence of any coherent drive as a particular property of a quantum non-Hermitian system featuring a gain medium.
Our work has several important consequences. The set-up of a DQD in a cQED cavity is already state-of-theart [69-72, 76, 77, 98-101] and connecting several cavities is also a standard technique in cQED experiments [109]. We have shown that the resulting set-up provides a highly tunable system where the exotic physics of exceptional point degeneracies in the quantum regime can be explored, along with a systematic study of quantum fluctuations that accompany any gain medium [93,96]. This set-up is potentially scalable [110,111]. Moreover, since a completely controllable gain can be introduced in an arbitrary cavity simply with a DQD in that cavity, our work presents the possibility of creating non-Hermitian bosonic systems with arbitrarily distributed gain and loss in the quantum regime. This is of particular importance for extending the forefronts of non-Hermitian topological systems [13,25,27,40,44,[52][53][54][55][56][57][58][59] into the quantum domain, as well as their possible applications in quantum technology [95,112].
It is important to emphasize that a rigorous description starting from a microscopic model showing the emergence of effective parity-time symmetric Hamiltonian has been missing so far. Most models [35][36][37][38][39] start with a non-Hermitian Hamiltonian (containing complex frequencies) that is applicable only in the classical limit where the number of energy quanta is much larger than one. For dissipative open quantum systems, one can write down a Lindblad equation (that encodes losses) from which one can get a non-Hermitian Hamiltonian by ignoring quantum jumps [83][84][85]. Generalizing this approach to gainloss systems [86], however, may be problematic. This is because the emergence of a loss-gain Lindblad equation starting from a larger set-up involving system and reservoirs is far from obvious and is not guaranteed. Starting from a microscopic description involving a system and multiple reservoirs, and then carefully integrating out the reservoirs while taking into account the hierarchy of time scales circumvents the problem. Our work is therefore a major step forward from current phenomenological constructions. Such a microscopic description in terms of operator equations of motion also facilitates calculation of two-time correlations, which is often a challenge in the master equation approaches. We will investigate these directions in future works. Although, our work was based on using a voltage-biased DQD as a gain medium, our findings and methods can be adapted to other potential quantum non-Hermitian systems with a gain medium governed by population inversion. where I is the 2x2 identity matrix and σ x is the usual Pauli matrix. The matrix H S is digonalized by the following orthogonal matrix, The isolated dynamics of the two cavities is then given byb Now, we go ahead to find effective evolution equations for the cavity operators in presence of bosonic baths and the DQD, by integrating out everything other than the cavities.

Setting the relative strength of various energy scales
In order to integrate out everything other than the cavities, we need to keep track of the relative strengths of various energy scales. To keep track of relative strengths of various energy scales, we introduce a dimensionless small parameter 1 and consider Our goal is to integrate out the bosonic baths as well as the DQD unit to obtain effective equations of motion for the two coupled cavities up to O( 4 ). In the following, we do this part by part.

Integrating out bosonic baths
The cavities are assumed to be coupled to the bosonic baths as well as to the DQD unit at time t = 0. The state of the entire set-up unit at time t = 0 is taken to be where β is the inverse temperature of the baths,Ĥ The formal solution for the bosonic bath operators is given bŷ (A7) Using these formal solutions in Eq.(A6), we obtain, Because of our chosen initial state, we get, is the bose distribution. Now, using the definitions in Eq.(A2), we note that, So we obtain up to O( 4 ), Let τ B be the time in which If J (ω) is such that it is nearly flat around ω 0 on a scale much larger than λ, and is given by J (ω) κ , then to a good approximation, we can write, So, finally, we obtain, In above effective equations of motion, we have the dissipative terms coming from the bosonic baths. If there was no boson-fermion coupling, i.e, g = 0, we would have had time evolution by an effective non-Hermitian dissipative Hamiltonian. However, there would be an inhomogenous part to the equation coming from the 'noise' terms, ξ b (t) associated with the corresponding dissipations. The relavant noise correlations are given in Eq.(A10). The 'noise' and dissipation coming from the baths satisfy the fluctuation-dissipation relations and embody quantum, as well as thermal, fluctuations. The connection to the DQD unit is coming from the term i 2 gÂ † 2Â 1 in above equation. Next we will integrate out the DQD unit.

Integrating out fermions
The DQD is assumed to be connected to the fermionic leads and the phononic substrate at time t = t 0 , with t 0 0. The state of the DQD unit at time t = t 0 is taken to be s â s is the Hamiltonian of the fermionic lead attached to the th site of the DQD, N s â s is the total number operator of the lead, Z f , Z ph are corresponding partition functions, β, µ 1 and µ 2 are the corresponding inverse temperature and the chemical potentials. The initial state of the DQD is ρ DQD , which is arbitrary.
Following exactly similar steps as for the bosons, the effective equation of motion forÂ † 2Â 1 , after integrating out the fermionic baths, can be written up as where,ξ The formal solution of this equation can be written aŝ Now, we note that the first four lines of above equation give the formally exact time evolution ofÂ † 2Â 1 , when the DQD unit is not connected to the cavity unit, i.e, with g 0 = 0. Let us now definê With this definition, Eq.(A20) becomeŝ (A22) Using above equation in equation of motion ofb 1 (Eq.(A16)), we have, It can be checked thatN (t ) =N (t) + O( 2 ), so that, up to O( 4 ), we have the following simplification in the second line of the above equation, where ∆N (t) ≡N 1 (t) −N 2 (t), and again, we have used the definitions in Eq.(A2). Since the term is already O( 4 ), we can ignore the boson-fermion coupling in obtaining ∆N (t). Up to O( 0 ), ∆N (t) is a conserved quantity. This means, assuming our initial state was a product of cavity and DQD units, in above equation, we can write This can only be done if the DQD-cavity coupling is small. We are interested in time (t − t 0 ) 1/( 2 Γ). In this time regime, the DQD, when not connected to the cavity unit, will reach a non-equilibrium steady state (NESS). We thus have is the NESS expectation value of ∆N (t) with g = 0.
Thus, in Eq.(A26), we have an effective gain coming from the DQD if ∆N ss > 0. In other words, the DQD need to be population inverted. This is acheivable in the DQD set-up for the detuning ε > 0 if conditions given in Eq. (12) are satisfied, as can be separately checked. Associated with the gain is a fluctuation term, 2ξ A (t) coming from the DQD, thereby satisfying the fluctuation dissipation relation. We are interested in obtaining the dynamics of expectation values of cavity field operators and the cavity bilinears. For obtaining the expectation value of the cavity field operators, we need where ... ss refers to NESS expectation value in absence of coupling to the cavity. For obtaining the cavity bilinears, we will require the following correlation function, To evaluate the above correlation function, we note the following result which can be derived using Eq.(A20), and integrating out the fermionic leads following exactly analogous procedure as used in the previous section for integrating out the bosonic baths, with t > t . In above, we have kept the factor 2 Γ in the exponent even though it is O( 2 ). This is because ( 2 Γ) −1 gives the time scale for the DQD to reach steady state, which in turn, gives the linewidth to the power spectrum of the noise appearing in the cavity due to the DQD. As we will see below, this leads to an O( 2 ) contribution in ξ † A (t 1 )ξ A (t 2 ) . The equivalent equation for A † 1 (t)Â 2 (t) is obtained by taking Hermitian conjugate. Now, to evaluate the correlation function in Eq.(A29), we consider two cases, t 1 < t 2 and t 1 > t 2 separately. If t 1 < t 2 , we use above formula to relateÂ † 2 (t 2 )Â 1 (t 2 ) tô A † 2 (t 1 )Â 1 (t 1 ). If t 1 > t 2 , we use the Hermitian conjugate of above formula to relatedÂ † , the combined result for both cases can then be written as where we have used fermionic commutation relations for simplification and used the fact that for V µ 1 ωq 2 , the DQD cannot be doubly occupied, so N 1N2 ss = 0. In above, P (ω) is the power spectral density of the noise appearing in the cavity due to the presence of the DQD. Note that the leading contribution in the integral comes from P (ω q ) ∝ −2 . Thus, 4 ξ † A (t)ξ A (t ) ∝ 2 , as we previously said.
So, finally, after integrating out the bosonic baths, as well as the DQD unit, we have the following effective equation of motion of the cavity operators, , where we have dropped the 's for notational convenience. Above equation can be simplified to Eq.(16) of main text using Eq.(A2), the fact that ω q = ω 0 and neglectingξ b (t). We justify the reason for neglectingξ b (t) of the bosonic baths in following subsection.

Reason from fluctuations of bosonic baths being negligible
The quantitiesξ b (t) embody the thermal and quantum fluctuations due to the presence of the bosonic baths. Since ξb (t) = 0, they do not play any role in obtaining expectation values of cavity field operators. However, they do play a role in obtaining cavity bilinears via their two-time correlator (Eq.(A10)). But, at low temperatures, we can write their two-time correlator as where the upper limit of the intergral has been replaced by 1/β because the bose distribution suppresses the contribution from βω 1. Since we assume that the frequency of the cavities ω 0 satisfies βω 0 1, this means the main contribution of the fluctuations come from bath frequencies which are extremely off-resonant with the cavity frequency. On the other hand, the DQD is in resonance with the cavity frequency, since ω q = ω 0 . As a result, the contribution from fluctuations of the gain medium is always much larger that that from the loss medium, i.e, the bosonic baths. This is why, for the choice of energy scales in our problem, we can neglect ξ b (t). We have also numerically confirmed this by keeping the thermal fluctuations, and observing that they make no change in the results.
Appendix B: Obtaining NESS results for DQD-unit

The general Redfield Quantum Master Equation
To obtain the NESS results for DQD-unit without the cavity-unit, we will take the approach of the Redfield Quantum Master Equation (RQME). For an arbitrary set-up of a system connected to a bath, the full set-up can be taken as isolated and described via the full sys-tem+bath HamiltonianĤ =Ĥ S +Ĥ SB +Ĥ B . HereĤ S is the system Hamiltonian,Ĥ B andĤ SB is the system-bath coupling Hamiltonian. The system-bath coupling can be written in the general formĤ SB = Ŝ B , wherê S is some system operator andB is some bath operator and is a small parameter controlling the strength of system-bath coupling. The initial density matrix of the set-up ρ tot (0) is considered to be in product form ρ tot (0) = ρ(0) ⊗ ρ B , where ρ(0) is some initial state of the system, and ρ B is the initial state of the bath. The initial state of the bath is often taken to be the thermal state with respect to the bath Hamiltonian. The Redfield Quantum Master Equation (RQME) is obtained by writing down the equation of motion for the reduced density matrix of the system up to O( 2 ) under Born-Markov approx. If B B = 0, where ... B refers to the average taken only with respect to bath, is satisfied initially, then, the RQME is given by whereŜ m (t) = e iĤ S tŜ m e −iĤ S t ,B m (t) = e iĤ B tB m e −iĤ B t . This gives the leading order dissipative term. Note that if system-bath coupling Hamiltonian is O( ), the dissipative part of RQME is O( 2 ).
The RQME for DQD-unit, when not connected to the cavity, is obtained by using Eq.(B1) and simplifying. The RQME is given by, Hereν is the complement of ν (i.e, if ν = 1,ν = 2 and vice-versa, this convention will be followed throughout). We are interested in the regime, V µ 1 , µ 2 . In this regime, there is negligible probability for double occupancy of the DQD, so N 1N2 0. Assuming this, and taking expectation values, we have where again, we have dropped the 's for notational convenience. To obtain N 1 ss , N 2 ss , Â † 1Â 2 ss , the above equations are solved with the LHS set to zero. Â † 2Â 1 ss is complex conjugate of Â † 1Â 2 ss .