Fermi’s Golden Rule and the Second Law of Thermodynamics

Abstract

We present a Gedankenexperiment that leads to a violation of detailed balance if quantum mechanical transition probabilities are treated in the usual way by applying Fermi’s “golden rule”. This Gedankenexperiment introduces a collection of two-level systems that absorb and emit radiation randomly through non-reciprocal coupling to a waveguide, as realized in specific chiral quantum optical systems. The non-reciprocal coupling is modeled by a hermitean Hamiltonian and is compatible with the time-reversal invariance of unitary quantum dynamics. Surprisingly, the combination of non-reciprocity with probabilistic radiation processes entails negative entropy production. Although the considered system appears to fulfill all conditions for Markovian stochastic dynamics, such a dynamics violates the Clausius inequality, a formulation of the second law of thermodynamics. Several implications concerning the interpretation of the quantum mechanical formalism are discussed.

1 Introduction

The probabilistic nature of quantum physics is related to a process whose correct interpretation and mathematically sound formulation is still under debate, the so-called collapse of the wave function [1]. The quantum mechanical collapse is of fundamental importance for our common-sense concept of macroscopic reality. In most cases, collapses are invoked in the context of measurements of quantum observables. But extending the applicability range of quantum mechanics beyond the microscopic realm prompts the question whether the discontinuous change of the wave function during the measurement is a physical process or not. Schrödinger’s famous Gedankenexperiment involving a cat demonstrates drastically the consequences of the assumption that the collapse process is merely epistemic and thus concerns only the knowledge of the observer [2]. Schrödinger’s goal was to demonstrate that a paradox arises if the collapse is not considered as a microscopic process taking place independently from the presence of an observer. To explain the absence of macroscopic superpositions without the assumption of a physical collapse, the “decoherence” interpretation considers not individual systems but the density matrix of an ensemble which becomes mixed after tracing over unobserved environmental degrees of freedom [3,4,5]. Although this avoids the need for a “measurement apparatus”, the tracing operation is not physical but epistemic. Without a physical mechanism triggering a real collapse event as assumed e.g. in [6], the decoherence interpretation is equivalent to the many-world hypothesis [7].

Already in Schrödinger’s cat example, the transition from unitary and deterministic to probabilistic evolution is tied to a microscopic, unpredictable event, the decay of a radioactive nucleus. There are strong arguments from a fundamental perspective supporting the occurence of spatio-temporally localized events as objective and observer-independent equivalent of collapse processes [8].

These random events, called “quantum jumps” in the early debate between Schrödinger and Bohr [9], are generally thought to underlie the statistical character of the emission and absorption of light quanta by atoms. Einstein used arguments from the theory of classical gases to derive Planck’s formula by assuming detailed balance between the atoms and the radiation in thermal equilibrium. His derivation did not require a microscopic Hamiltonian [10]. Nevertheless, it is easy to derive the corresponding rate equations for the time-dependent occupancy $\langle n\rangle (t)$ of a light mode with frequency $\Omega$ coupled resonantly to M two-level systems (TLS) within quantum mechanics. The transition probabilities follow from the microscopic interaction Hamiltonian by employing Fermi’s golden rule [11, 12], which tacitly incorporates the collapse event by replacing the unitary time evolution by a stochastic process.

The interaction Hamiltonian is well known [13, 14]:

$$\begin{aligned} H_{\text {int}}=g\sum _{l=1}^M\left( a\sigma ^+_l+a^{\dagger }\sigma ^-_l\right) , \end{aligned}$$

(1)

where $a,a^{\dagger }$ denote the annihilation/creation operators of the radiation mode and $n=a^{\dagger }a$. The Pauli lowering/raising operators $\sigma ^-_l,\sigma ^+_l$ describe the l-th two-level system with $H_{\text {TLS}}=\hbar \Omega (\sum _l\sigma _l^+\sigma _l^--1/2)$. The probabilities for a transition between the upper and lower state of a TLS accompanied by the emission or absorption of a photon with frequency $\Omega$ are computed with the golden rule (see “Derivation of the Rate Equations (2) and (3)”) to obtain the rate equations

$$\begin{aligned} \frac{\text {d}\langle n\rangle }{\text {d}t}&= \gamma \left[ \langle m\rangle (\langle n\rangle +1)- (M-\langle m\rangle )\langle n\rangle \right] , \end{aligned}$$

(2)

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t}&= \gamma '\left[ (M-\langle m\rangle )\langle n\rangle -\langle m\rangle (\langle n\rangle +1)\right] . \end{aligned}$$

(3)

Here, $\langle m(t)\rangle$ denotes the time-dependent average number of excited TLS. The rate constants $\gamma ,\gamma '$ depend on the coupling g and the density of states of the radiation continuum around $\Omega$ (see below). These equations describe the irreversible change of average quantities and thus use the ensemble picture of statistical mechanics [15]. Nevertheless they account for the temporal behavior of a single system as well, because a typical trajectory will exhibit a fraction m(t)/M of excited TLS close to $\langle m(t)\rangle /M$ for sufficiently large M [16,17,18].

It is crucial that the rate equations (2) and (3) satisfy the detailed balance condition, which implies that they lead from arbitrary initial values $\langle n\rangle (0), \langle m(0)\rangle$ to a unique steady state characterized by

$$\begin{aligned} \langle m\rangle (\langle n\rangle +1) = (M-\langle m\rangle )\langle n\rangle . \end{aligned}$$

(4)

We have the relations

$$\begin{aligned}&P^e=\langle m\rangle /M,\quad P^g=(M-\langle m\rangle )/M, \nonumber \\&P_{e\rightarrow g}=\gamma '(\langle n\rangle +1),\quad P_{g\rightarrow e}=\gamma '\langle n\rangle , \end{aligned}$$

(5)

for the probabilities $P^g$ ($P^e$) for a TLS to be in its ground (excited) state and the probabilities for emission $P_{e\rightarrow g}$ and absorption $P_{g\rightarrow e}$. Equation (4) can therefore be written as

$$\begin{aligned} P^eP_{e\rightarrow g} = P^gP_{g\rightarrow e} \end{aligned}$$

(6)

which is the definition of detailed balance [15, 19]. Equation (4) entails Planck’s formula for thermal equilibrium between radiation and matter. If one considers (3) as an equation of motion for the probability $P^e(t)$, even a microscopic system consisting of a single TLS will thermally equilibrate with the surrounding continuum of radiation.

This surprising result rests on the fact that the TLS does not interact with the light mode exactly on resonance only but with all modes in a frequency interval of width $\Delta$ around $\Omega$ with a similar strength g [13, 14]. The energy uncertainty $\hbar \Delta$ allows for a radiation event occuring during a short time span $\tau _c\sim \Delta ^{-1}$, whereas the process itself is energy conserving (see [13], p. 419). The coupling to a continuum of modes leads therefore to real and irreversible microscopic processes, the emission or absorption of light quanta, although no macroscopic measurement apparatus is involved. Such a microscopic collapse process is tacitly assumed whenever the golden rule is employed. If, however, the TLS is embedded into a cavity and coupled only to a single radiation mode, the collapse cannot take place; the system shows the Rabi oscillations of an unitarily evolving state instead, the TLS being entangled with the bosonic mode. To this case, the golden rule cannot be applied. Therefore, the golden rule cannot be taken as an approximation to the full unitary time development given by solving the Schrödinger equation, although it corresponds formally to a perturbative computation of the unitary dynamics for short times [20]. The very concept of a transition rate implies that the deterministic evolution of the state vector is replaced with a probabilistic description of events. The central element in the computation is the overlap $\langle \psi _{\text {final}}|H_\text {int}|\psi _{\text {initial}}\rangle$, which, according to the Born rule, determines the probability for the transition from $|\psi _{\text {initial}}\rangle$ to $|\psi _{\text {final}}\rangle$. The Born rule is invoked here although the process is not a macroscopic measurement. The “macroscopic” element is provided by the continuum of radiation modes. The presence of this continuum causes the necessity of a statistical description [13].

Within the decoherence interpretation, one would require that the continuum of radiation modes acts as an environment for the TLS, the “system”. The environment is then traced out to yield the irreversible dynamics of the system alone. But if the system consists of the walls of a hohlraum, it excercises a strong influence on the “environment”, the enclosed radiation, such that they equilibrate together. Therefore, it is not possible to separate the environment from the system to explain the microscopic collapse processes driving the compound system towards thermal equilibrium.

The rate equations neglect completely the coherences of the TLS. The irrelevance of the coherences follows naturally from the interpretation of the radiation process as a collapse: each such event projects the TLS into their energy eigenbasis, just as a macroscopic measurement projects any quantum system into the basis entangled with the eigenbasis of the measurement device [21]. A macroscopic measurement device is not needed here because both subsystems, the TLS and the radiative continuum, contain a macroscopic number of degrees of freedom. This alone seems to justify a statistical description as in classical gas theory, although the interaction Hamiltonian (1) has no classical limit. Indeed, in our case the Born rule replaces the assumption of “molecular chaos”, needed in Boltzmann’s derivation of the H-theorem [15, 19]. Therefore, it seems almost natural that the rate equations (2), (3) lead to thermal equilibrium from a non-equilibrium initial state although the radiation and the collection of atoms are both treated as ideal gases. If $\langle m(0)\rangle$ and $\langle n\rangle (0)$ correspond to equilibrium ensembles with different temperatures at $t=0$,

$$\begin{aligned} T_{\text {TLS}}(0)=\frac{\hbar \Omega }{k_B\ln \left( \frac{M}{\langle m(0)\rangle }-1\right) }, \quad T_{\text {rad}}(0)= \frac{\hbar \Omega }{k_B\ln \left( \frac{1}{\langle n\rangle (0)}+1\right) }, \end{aligned}$$

(7)

the rate equations derived from the quantum mechanical interaction Hamiltonian (1) together with the golden rule entail thermal equilibration according to Clausius’ formulation of the second law of thermodynamics: the two gases exchange heat which flows from the hotter to the colder subsystem until a uniform temperature and maximum entropy of the compound system is reached [15, 19].

Certainly, the rate equations do not correspond to the exact quantum dynamics of the system, which is non-integrable in the quantum sense if the interaction (1) is generalized to a continuum of bosonic modes, rendering it equivalent to the spin-boson model [22]. To obtain the exact evolution equation for the density matrix of both the TLS and the radiative modes, one would have to solve the full many-body problem. The golden rule is then seen as a method to approximate the time-dependent expectation values $\langle n\rangle (t)$, $\langle m(t)\rangle$, which is justified by a large body of experimental evidence, but not through an analytical proof of equivalence between both approaches. It may even be that the golden rule provides a phenomenological description of microscopic collapse processes whose actual dynamics is not yet known. In this case, the full quantum mechanical calculation would not yield (2), (3), although they account correctly for the experimentally observed dynamics. The rate equations are derived from the interaction Hamiltonian (1) and the golden rule (see  “Derivation of the Rate Equations (2) and (3)” and “Derivation of the Rate Equations (17)–(19) for the Open System”), thereby demonstrating its applicability in the situation considered. This corroboration of the golden rule by experimental results is independent of any theoretical justification of its validity.

The interaction Hamiltonian (1) satisfies the detailed balance condition which is crucial for the physically expected behavior. We shall demonstrate in the next section that this condition does not follow from elementary principles like time reversal invariance or hermiticity as (1) seems to suggest. On the contrary, there are feasible experimental setups violating the detailed balance condition while satisfying all other prerequisites for the application of the golden rule. The consequence of this violation is a macroscopic disagreement with the second law of thermodynamics.

2 The Gedanken Experiment

In the Gedankenexperiment, we consider two identical cavities A and B supporting a quasi-continuous mode spectrum described by bosonic annihilation operators $a_{Aj}$, $a_{Bj}$ and frequencies $\omega _j$ (Fig. 1). They are coupled bilinearly to the right- and left-moving modes $a_{1k}$, $a_{2k}$ of an open-ended waveguide which form a quasi-continuum like the cavity modes [23]. The loss processes through the open ends of the waveguide are caused formally by a heat bath at $T=0$, leading to a reduction of the system entropy through heat transport. This coupling to the outside world is one argument for the applicability of the golden rule, the second is the already discussed quasi-continuum of modes. Even if the full continuum of radiation modes in the cavities is treated as part of the “system”, which would then be subject to purely unitary time evolution, the waveguide would still couple to a decohering “environment”, justifying the statistical description even if one denies that real collapse events take place in the system itself. Here we study a model which is commonly used in quantum optics to describe unidirectional loss processes [24]. In this way, the generally accepted arguments substantiating irreversible evolution equations can be transferred to the present situation.

Fig. 1

Layout of the Gedankenexperiment. Two reservoirs A and B containing black-body radiation are coupled via a non-reciprocal, open waveguide to a collection of two-level systems (TLS) and to each other. The right-, respectively left-moving modes in channels 1 and 2 couple with different parameters $g_1$ and $g_2$ to the two-level systems

The modes $a_{1k}$ and $a_{2k}$ of the waveguide are coupled to a collection of M two-level systems located at the center of the waveguide (see Fig. 1). The total Hamiltonian is given by

$$\begin{aligned} H=H_A + H_B + H_{wg} + H_{\text {TLS}} + H^1_{\text {int}} + H^2_{\text {int}}. \end{aligned}$$

(8)

$H_q$ denotes the Hamiltonian in cavity q for $q=A,B$,

$$\begin{aligned} H_q=\hbar \sum _{j}\omega _j a_{qj}^\dagger a_{qj}. \end{aligned}$$

(9)

The Hamiltonian of the waveguide reads

$$\begin{aligned} H_{wg}=\hbar \sum _{k}\omega _k \left( a_{1k}^\dagger a_{1k} + a_{2k}^\dagger a_{2k}\right) , \end{aligned}$$

(10)

where the modes 1 and 2 belong to waves traveling to the right and to the left, respectively. For the TLS we have $H_{\text {TLS}}=(\hbar \Omega /2)\sum _{l=1}^M\sigma ^z_l$, with the Pauli matrix $\sigma ^z$. The coupling between the reservoirs and the modes 1 and 2 of the waveguide is bilinear,

$$\begin{aligned} H^1_{\text {int}}=\sum _{q=A,B}\sum _{j,k} h_{jk}\left( a_{qj}^\dagger [a_{1k} + a_{2k}] +\text {h.c.}\right) . \end{aligned}$$

(11)

Using the rotating wave approximation, the interaction with the TLS has the standard form [13, 14] which is equivalent to (1),

$$\begin{aligned} H^2_{\text {int}}=\sum _{l=1}^M \left( \sum _{k} g_{1k}a_{1k} +g_{2k}a_{2k}\right) \sigma ^+_l + \text {h.c.}, \end{aligned}$$

(12)

where $\sigma ^+_l$ denotes the raising operator of the lth TLS. We consider in the following the (time-dependent) average occupancy per mode j, $\langle n_q\rangle (t)$ for reservoir $q=A,B$ in an energy interval around the TLS energy, $\Omega -\Delta /2<\omega _j<\Omega +\Delta /2$, where $\Delta$ is much larger than the natural linewidth of spontaneous emission from an excited TLS into the waveguide. The occupancy does not depend on j if the couplings $h_{jk}$, $g_{(1,2)k}$, the density of states $\rho _q(\hbar \omega _j)$ of the reservoirs and $\rho _{1,2}(\hbar \omega _k)$ of the waveguide are constant in the frequency interval of width $\Delta$ around $\Omega$.

It is crucial that $g_{1k}\ne g_{2k}$, which specifies that the TLS couple with unequal strengths to the right- and left-moving photons in the waveguide. Such an unequal coupling is a hallmark of chiral quantum optics [25]. Although the interaction Hamiltonian (12) appears to break time-reversal invariance, as the time-reversal operator maps left-moving to right-moving modes, this is actually not the case because the effective interaction term (12) does not contain the polarization degree of freedom. The angular momentum selection rules for light-matter interaction lead naturally to a dependence of the coupling strength on the propagation direction in engineered geometries [26], especially if the spin-momentum locking of propagating modes in nanofibers is employed [27, 28]. The unwanted coupling to non-guided modes can be effectively eliminated, leading to large coupling differences $|g_{1k}-g_{2k}|$ [29, 30].

We shall now study the temporal behavior of the two cavities, assuming at time $t=0$ separate thermal equilibria in A, B and the TLS system, all at the same temperature T. The probability $\langle m(0)\rangle /M$ for a TLS to be excited obeys the Boltzmann distribution

$$\begin{aligned} \frac{\langle m(0)\rangle }{M}=\frac{1}{e^{\hbar \Omega /k_BT}+1}. \end{aligned}$$

(13)

Considering the M TLS as independent classical objects, their Gibbs entropy is given by

$$\begin{aligned} S_M=-k_BM\big (p_e\ln p_e +(1-p_e)\ln (1-p_e)\big ), \end{aligned}$$

(14)

with $p_e(t)=\langle m(t)\rangle /M$. This approach is justified by the quick relaxation of the two-level systems by non-radiative processes, which decohere them on time scales much shorter than the time scale of spontaneous emission and quickly quench finite coherences of the TLS [31]. This argument for a classical description of the TLS is independent from the general justification of the golden rule via the mode continuum discussed above. Equation (14) provides an upper bound for the entropy of the TLS subsystem [32]. Analogously, the entropy of the compound system reads

$$\begin{aligned} S^{\text {sys}}=S_A +S_B + S_M, \end{aligned}$$

(15)

where $S_q$ denotes the v. Neumann entropy of the radiation in cavity q for $q=A,B$. The average occupation number $\langle n_q(\omega )\rangle$ per mode at frequency $\omega$ follows from the Bose distribution

$$\begin{aligned} \langle n_q(\omega )\rangle (0)=\frac{1}{e^{\hbar \omega /k_BT}-1}. \end{aligned}$$

(16)

Because the temperature depends on $\langle n_q\rangle$ as described by (16), we can define effective temperatures $T_q(t)$ by $\langle n_q\rangle (t)$ for each reservoir q under the assumption that the photons in each reservoir thermalize in the usual way quickly as a non-interacting Bose gas. Furthermore, we consider the case that $\langle n_{(1,2)k}\rangle (t)=0$ for the occupancy of the modes in the waveguide, i.e., the waveguide is populated through the sufficiently weak coupling to the reservoirs and its modes appear only as intermediate states (see “Derivation of the Rate Equations (17)–(19) for the Open System”). Neglecting the waveguide, the system is thus composed of three subsystems, each in thermal equilibrium at any time $t>0$, with locally assigned time-dependent temperatures. The subsystems interact through random emission and absorption processes which do not lead to entanglement, because in each such process the wave function undergoes a collapse towards a product state. This description is in obvious accord with the derivation of Planck’s law given by Einstein [10]. Note that the compound system is not coupled to several thermal baths which have different temperatures. In such a case, a description with separate master equations for the subsystems is inconsistent if the interaction between subsystems is still treated quantum mechanically. Using such a description, a violation of the second law has been deduced [33, 34], which is only apparent and caused by the inconsistent computation [34]. Our system differs from those models because the system dynamics is not described by a unitary evolution as in [33, 34], but by a random process, with all subsystems coupled to the same bath (the open waveguide).

The golden rule applied to the Hamiltonian (8) yields the rate equations for $\langle n_q\rangle (\Omega ,t)$ and $\langle m(t)\rangle$,

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t}&= -2\gamma _\text {dec}\langle n_{A}\rangle + \gamma _0(-\langle n_{A}\rangle +\langle n_{B}\rangle )- \gamma _1(M-\langle m\rangle )\langle n_{A}\rangle \nonumber \\&\quad + \gamma _2\langle m\rangle (\langle n_{A}\rangle +1), \end{aligned}$$

(17)

$$\begin{aligned} \frac{\text {d}\langle n_{B}\rangle }{\text {d}t}&= -2\gamma _\text {dec}\langle n_{B}\rangle + \gamma _0(-\langle n_{B}\rangle +\langle n_{A}\rangle )- \gamma _2(M-\langle m\rangle )\langle n_{B}\rangle \nonumber \\&\quad + \gamma _1\langle m\rangle (\langle n_{B}\rangle +1) , \end{aligned}$$

(18)

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t}&= -[\tilde{\gamma }_{11}(\Omega ) +\tilde{\gamma }_{12}(\Omega )]\langle m\rangle + \tilde{\gamma }_1(\Omega )\left[ \langle n_{A}\rangle (M-\langle m\rangle )-(\langle n_{B}\rangle +1)\langle m\rangle \right] \nonumber \\&\quad +\tilde{\gamma }_2(\Omega )\left[ \langle n_{B}\rangle (M-\langle m\rangle )-(\langle n_{A}\rangle +1)\langle m\rangle \right] . \end{aligned}$$

(19)

The first terms on the right hand side of (17)–(19) describe the loss of photons through the open ends of the waveguide. These terms are of first order in $|h_{jk}|^2$, resp. $|g_{1k}|^2$, $|g_{2k}|^2$. The following terms correspond to coherent processes of second order in the couplings. The effective rates $\gamma _{\text {dec},0,1,2}$ and $\tilde{\gamma }_{1r},\tilde{\gamma }_{1,2}$ used for the numerical solution of (17)–(19) shown in Figs. 2 and 3 belong to the strong coupling regime of the TLS and the waveguide with values accessible within a cavity QED framework [35]. The chiral nature of the coupling, $g_{1k}\ne g_{2k}$, entails $\gamma _1\ne \gamma _2$. In our example we have assumed $\tilde{\gamma }_2=\gamma _2=\tilde{\gamma }_{12}=0$, i.e., channel 2 is not coupled to the TLS. One sees from (17)–(19) that the chiral coupling leads to a breakdown of the detailed balance condition in second-order processes because absorption is no longer balanced by stimulated and spontaneous emission. The radiation processes generate an effective transfer of photons from reservoir A to B on a time scale $\tau _{\text {char}}$ given by the strong coupling between channel 1 and the TLS, $\tau _{\text {char}}\sim \tilde{\gamma }^{-1}=0.1\;\mu$s. This corresponds to a difference in the local temperatures calculated via (16). Figure 2 shows the temporal behavior of the temperatures of reservoirs A, B and the TLS for intermediate times.

Fig. 2

Solutions of the rate equations (17)–(19) as function of time, starting from initial thermal equilibrium. The effective temperatures of the reservoirs A and B deviate. The temperature drop of the TLS parallels that of A. Parameters used are $\gamma _\text {dec}=\gamma _0=10$ kHz, $\tilde{\gamma }_{11}=\tilde{\gamma }_1=10$ MHz, $\gamma _1=100$ kHz; $\tilde{\gamma }_{12}=\tilde{\gamma }_2=\gamma _2=0$ and $\hbar \Omega /kT(0)=1$. $\Omega$ corresponds to a wavelength of 10 $\mu$m

Although both reservoirs A and B loose photons through the open waveguide, the ratio between $\langle n_{A}\rangle$ and $\langle n_{B}\rangle$ attains a constant value for $t\rightarrow \infty$, which is depicted in Fig. 3. This figure reveals that the losses of reservoir B are characterized by a much larger time scale than the spontaneous population of B through reservoir A and the TLS.

Fig. 3

The occupations of reservoirs A and B for long times, plotted on a logarithmic scale. The ratio $\langle n_{A}\rangle /\langle n_{B}\rangle$ is asymptotically time-independent. The losses through the open waveguide are much slower than the breakdown of local equilibrium caused by the non-reciprocal interaction with the TLS

3 Conflict with the Second Law of Thermodynamics

Due to the steady loss of photons, our system is always out of equilibrium and the only steady state solution of the Eqs. (17)–(19) corresponds to empty cavities and all TLS in their ground state. It is clear that the system entropy¹ diminishes at any time due to the outgoing heat flow. The question arises how to apply the second law of thermodynamics to this situation. The second law has been formulated in several versions (see, e.g., [36,37,38]). It is not the purpose of this paper to discuss these in detail. Three representative formulations provide exemplary definitions of the second law [15, 19]:

• (1) The entropy of the universe always increases.

• (2) The entropy of a completely isolated system stays either constant or increases.

• (3) The entropy of a system thermally coupled to the environment satisfies Clausius’ inequality: $\Delta S \ge \Delta Q/T$.

Variant (1) is not subject to our Gedankenexperiment, because the entropy production including the environment is formally infinite, as the external bath has zero temperature. Variant (2) follows from variant (3) because the heat transfer $\Delta Q$ vanishes. Variant (3) applies to the present system: the heat transfer $\Delta Q$ is negative and therefore also $\Delta S$ may be negative. However, the second law in the form of Clausius’ inequality forbids the case $\Delta S \le \Delta Q/T$: The local entropy production $\sigma$ in the system must be non-negative [19].

The change of system entropy $S^{\text {sys}}(t)$ can be written as [39]

$$\begin{aligned} \frac{\text {d}S^{\text {sys}}}{\text {d}t} = -\int \text {d}{\varvec{o}}\cdot {\varvec{J}}_{\text {sys}} +\sigma , \end{aligned}$$

(20)

where the surface integral of the entropy current ${\varvec{J}}_{\text {sys}}$ accounts for the heat transfer to the environment, characterized by the rates $\gamma _\text {dec}$ and $\tilde{\gamma }_{11},\tilde{\gamma }_{12}$. We find for $\sigma (t)$ (see “Derivation of the Entropy Production (22)–(24)”)

$$\begin{aligned} \sigma (t)= \sigma _{A,B}(t) + \sigma _{A,\text {TLS}}(t) + \sigma _{B,\text {TLS}}(t), \end{aligned}$$

(21)

with

$$\begin{aligned} \sigma _{A,B}&=k_B\mathcal{N}\gamma _0(\langle n_{B}\rangle - \langle n_{A}\rangle )\big (\ln (\langle n_{B}\rangle [\langle n_{A}\rangle +1])-\ln (\langle n_{A}\rangle [\langle n_{B}\rangle +1])\big ), \end{aligned}$$

(22)

$$\begin{aligned} \sigma _{A,\text {TLS}}&=k_B\big (\tilde{\gamma }_2\langle m\rangle (\langle n_{A}\rangle +1)-\tilde{\gamma }_1\langle n_{A}\rangle (M-\langle m\rangle )\big )\big (\ln (\langle m\rangle [\langle n_{A}\rangle +1]) \nonumber \\&\quad -\ln ([M-\langle m\rangle ]\langle n_{A}\rangle )\big ), \end{aligned}$$

(23)

$$\begin{aligned} \sigma _{B,\text {TLS}}&=k_B\big (\tilde{\gamma }_1\langle m\rangle (\langle n_{B}\rangle +1)-\tilde{\gamma }_2\langle n_{B}\rangle (M-\langle m\rangle )\big )\big (\ln (\langle m\rangle [\langle n_{B}\rangle +1]) \nonumber \\&\quad -\ln ([M-\langle m\rangle ]\langle n_{B}\rangle )\big ), \end{aligned}$$

(24)

with $\mathcal{N}=\tilde{\gamma }_1/\gamma _1=\tilde{\gamma }_2/\gamma _2$. The three contributions result from the heat exchange between the three subsystems. Only $\sigma _{A,B}$ is always non-negative, because it has the form $(x-y)(\ln (x)-\ln (y))$ characteristic for systems satisfying the detailed balance condition. Because $\tilde{\gamma }_1\ne \tilde{\gamma }_2$, the two other contributions are not necessarily non-negative. Figure 4 shows $\sigma (t)$ calculated using the parameters of Fig. 2. The entropy production is negative for $t<0.75\;\mu$s. The time span during which the entropy is reduced beyond the loss to the environment is almost an order of magnitude longer than the characteristic time scale $\tau _{\text {char}}= 0.1\;\mu$s. The especially notable point is the entropy reduction which occurs although the dynamical evolution started with thermal equilibrium between the subsystems.

During the exchange of heat with the other subsystems, each subsystem remains in local thermal equilibrium. In principle, the exchange of heat between subsystems gives only a lower bound to the entropy production for irreversible processes and the actual $\sigma (t)$ could be larger than the value given in (21), which contains the contributions from the mutual heat exchange. However, the entropy change in each subsystem can be computed directly via the formulae (14) and (96), (97) as well. Doing so, one finds that no additional entropy production besides the mutual heat exchange occurs in the process described by the rate equations (17)–(19). This, obviously, is due to the fact that during this process each subsystem remains in local thermal equilibrium. It follows that the entropy change of the full system would even be negative if the initial entropy of the TLS subsystem was lower than the upper bound given in (14), because the entropy change is caused solely by mutual heat transfer between subsystems if the condition of local equilibrium is satisfied.

We conclude that the initial thermal equilibrium between reservoirs and the TLS is unstable and variant (3) of the second law is violated. This is true although the total entropy of the system plus the environment always increases. The second law, applied to the system alone, demands a non-negative local entropy production for any process driven by the coupling to the bath at $T=0$ [19, 39]. Such processes may generate local temperature gradients between the reservoirs, but $\sigma (t)$ must always be larger or equal zero to satisfy Clausius’ inequality. The violation of this inequality in our Gedankenexperiment shows clearly that the chiral coupling between TLS and the waveguide generates radiation processes which are in conflict with thermodynamics if they are treated statistically in the same way as black body radiation.

Fig. 4

The total entropy production $\sigma (t)$ within the system calculated for $\mathcal{N}\gamma _0=1$ MHz, $\tilde{\gamma }_1=10$ MHz, $\tilde{\gamma }_2=0$. For these parameters, the entropy production is negative until $t'\sim 0.75\;\mu$s, when it turns positive and stays so. This behavior entails that for $t<t'$ the second law of thermodynamics is not satisfied. Note that $t'$ is appreciably larger than $\tau _{\text {char}}$

4 Discussion and Conclusions

In discussing the experiment we first note that it is based on phenomena taking place in the realm where the quantum world interfaces classical physics. The three subsystems are clearly macroscopic, but their interaction Hamiltonian (11, 12) is purely quantum mechanical and cannot be treated in a (semi-)classical approximation. A unique quantum feature of the interaction is given by the fact that the emission rate of the TLS depends on the occupation of the final states. In “Detailed Balance for a System with an Embedded Cavity”, we demonstrate that this counter-intuitive effect leads to the restoration of detailed balance in a cavity system without non-reciprocal elements. The non-unitary, probabilistic state development of the device can neither be achieved in classical Hamiltonian dynamics nor in the unitary pure quantum regime described by the Schrödinger equation. The collapse processes that link the classical world and the quantum regime [40] are the cause of the thermal imbalance between the otherwise equivalent reservoirs A and B. Our Gedankenexperiment therefore reveals a clear conflict between thermodynamics and the probabilistic description of quantum phenomena on a macroscopic scale. The corresponding rate equations (17)–(19) do not satisfy the condition for detailed balance. Instead they predict a time-dependent state that violates the Clausius inequality, i.e. the original formulation of the second law of thermodynamics.

According to the quantum description of the statistical absorption and emission processes, the chirally coupled cavities are expected to develop unequal occupation numbers. This imbalance creates a temperature gradient between them, although no work is done on the system, which is coupled to the environment through the open waveguide only. This coupling to a heat bath at temperature zero leads to heat flow out of the system which is usually accompanied by a positive local entropy production. But in our case the entropy production within the system is negative during a well defined time interval. This interval is larger than the time characterizing the coupling between the subsystems. The violation of the second law is temporary. Because this violation is described by the rate equations, it is not caused by a statistical fluctuation. Such a fluctuation may occur in the stochastic evolution of the state vector of a single system, even when the initial state is typical [41], but cannot appear in the fully deterministic equations for averages.

The rate equations (17)–(19) have been derived under the assumption that the two channels of the waveguide are fed by A and B through the emission of wavepackets which in turn interact with the TLS in a causal fashion. The emission and absorption of single photons by the TLS are considered thus as probabilistic processes taking place within a finite time span, due to the quasi-continuum of modes available in the waveguide and the reservoirs. They therefore satisfy causality: It is not possible for a right-moving photon in channel 1 to be emitted by the TLS and be subsequently absorbed by reservoir A. The photons entering A are either generated by a fluctuation in channel 1 of the open waveguide or arrive through channel 2. In the latter case, they may come from the outside, from a TLS or from reservoir B. Pure scattering events at the TLS are neglected in this approximation because they are of higher order in the coupling constants. Their inclusion cannot restore the detailed balance broken by the chiral coupling.

Our reasoning is based on the assumption that the interaction of the TLS with the radiation continuum leads to real events [8] which must be described statistically, and are therefore caused by a collapse process. The physical mechanism of this “real” collapse plays no role in these considerations because no hypothesis beyond the golden rule enters the derivation of the rate equations. Of course, we have also assumed that the macroscopic nature of the radiation and the collection of TLS removes any detectable entanglement between the subsystems. It entails the statistical independence of the radiation processes and therefore Markovian dynamics². This assumption is corroborated by all available experimental evidence up to now. In case the second law of thermodynamics would be correct and therefore the presented statistical analysis wrong, the validity of the second law would be tantamount to the actual realization of a macroscopic superposition of states in the cavity system, although it is coupled to an unobserved environment, the open ends of the waveguide.

Interestingly, the use of the golden rule can also be justified even in a completely isolated system. If the system is isolated, it could in principle be described by the full unitary dynamics, leading to a trivial reconciliation with the second law in its restricted form: the fine-grained Gibbs/v. Neumann entropy does not change at all. However, also in this case an “environment” is present which consists of the infinitely many degrees of freedom of the photon gas, decohering the dynamics of the TLS, at least according to the opinion of the majority of physicists working in quantum optics [13]. The effective coarse-grained description of the closed system proceeds again via the golden rule (see “Derivation of the Rate Equations for the Closed System”). A similar temperature difference between A and B appears and a new steady state develops from initial thermal equilibrium for $t\rightarrow \infty$, having a lower entropy than the initial state, thus violating variant (2) of the second law. This is shown in Fig. 5. The temperature difference corresponds to a “sorting” between the reservoirs A and B in the closed system and resembles the action of a Maxwell demon [43]. No information is processed, stored or erased, neither in the TLS nor in the waveguide [44], therefore the usual arguments for positive entropy production based on information theory [45, 46] cannot be applied.

Fig. 5

Temporal development of the total entropy for the closed variant of our system shown in the inset. The initial state at $t=0$ (thermal equilibrium between all subsystems, including the waveguide) maximizes the entropy. The stable steady state for long times has a lower entropy

It has been argued that it is not possible to discern experimentally an interpretation of quantum mechanics based on probabilistic dynamics and real collapse from the decoherence interpretation which replaces the physical collapse by an epistemic operation: the tracing over environmental degrees of freedom in the full density matrix at the final observation time $t_\text {fin}$ [47]. As mentioned in Sect. 1, it is not known whether the photon densities in the reservoirs at $t_\text {fin}$, calculated with the tracing procedure, would differ from the results of “Derivation of the Rate Equations (17)–(19) for the Open System” based on the golden rule. If so, our proposed experiment, if performed with an isolated system, would allow to decide between interpretations based on real collapses and epistemic interpretations. Only the latter do not contradict the second law of thermodynamics, provided an actual solution of the full many-body problem would effectively restore the detailed balance condition in the statistical description. Such a solution would also be necessary to identify possible reasons for the breakdown of well-established tools like the golden rule or the Markov approximation in case the system shows the equilibrium state predicted by thermodynamics at all times. In any case, a difference between the full solution and the approximation by the golden rule would entail another mystery: why is the approximation valid for black-body radiation in arbitrary cavities (see “Detailed Balance for a System with an Embedded Cavity”) but not for chiral waveguides?

As a computation of the full quantum dynamics appears out of reach at present, the question can only be decided experimentally. A direct implementation of our model appears feasable with current technology [23, 25, 35].

In case that the experiment reported unequal distributions in A and B for the closed system, one would be forced to conclude that statistical processes such as spontaneous emission and absorption are able to reduce the total entropy for arbitrary large isolated systems and on average, not only for short times and small systems as expected from fluctuations [41]. Then the state of classical thermal equilibrium with maximal entropy is unstable and the system moves to steady states with a lower entropy. In this case, entropy-reducing processes would be expected to occur actually in nature, in structures differing greatly from our device³.

The contradictions presented are rooted in the still unresolved status of the measurement problem of quantum physics: quantum mechanical probabilities can only be computed under the assumption that a collapse (either real or epistemic) takes place. These stochastic probabilities are formally encoded in the Born rule. To our knowledge, neither the Born rule nor the golden rule have ever been used for a derivation of the second law of thermodynamics [16]. There have been attempts to deduce the second law from quantum mechanics by employing several “coarse-graining” prescriptions. The proof of the H-theorem given by v. Neumann employs assumptions about the density matrix of pure states and macroscopic distinguishability but excludes explicitly collapse or measurement processes from the quantum dynamics. In v. Neumann’s approach the quantum dynamics stays always unitary [52, 53]. This may seem surprising, because the collapse processes underlying the golden rule are inherently irreversible, as noted first also by v. Neumann [21]. To our opinion, the intrinsic probabilistic features of quantum mechanics encoded in the Born rule add an elementary irreversible process to the dynamical laws of nature. This process, the non-unitary collapse of the wavefunction, appears during macroscopic measurements but also as microscopic event, thereby leading to the correct statistics of a photon gas interacting with matter. Therefore, we disagree with the position put forward in [54], that the irreversibility of the measurement is just due to the macroscopic nature of the apparatus and has essentially the same origin as the irreversible behavior of macroscopic variables in classical mechanics which obeys the second law of thermodynamics. We have shown that the presence of microscopic collapse processes may lead under certain circumstances to a conflict with this law.

In conclusion, transition rates of quantum systems are commonly calculated with great success by using Fermi’s golden rule. This approach is widely accepted, as the golden rule directly results from the Born rule. Here, we have introduced practically realizable, open and closed quantum systems of coupled cavities and determined their behavior by applying the golden rule. The predicted behavior of both systems violates the second law of thermodynamics. We therefore conclude that

• (1) the statistical description of quantum mechanical transitions given by the golden rule is incorrect or

• (2) the second law of thermodynamics is not universally valid.

Appendix

1.1 Derivation of the Rate Equations (2) and (3)

The simple, albeit nonlinearly coupled rate equations (2), (3) underlie Einstein’s famous derivation of Planck’s law for black-body radiation [10]. In this section, we describe the approximations and statistical assumptions leading to the closed system of Eqs. (2), (3).

We begin by describing the quantum state of a single system within the statistical ensemble. A collection of M two-level systems with energy splitting $\hbar \Omega$ interacts with radiation modes $a_j$, characterized by their frequencies in the range $[\Omega -\Delta /2,\Omega +\Delta /2]$ and additional parameters such as momentum and polarization, summarized in the index j. The full many-body quantum state $|\Psi (t)\rangle$ at time t is then described as

$$\begin{aligned} |\Psi (t)\rangle = |s_1,\ldots s_M;\{n_j\}\rangle , \end{aligned}$$

(25)

where the $s_l=0,1$ correspond to the ground and excited state $|g_l\rangle$ and $|e_l\rangle$ of the l-th TLS and $n_j$ is the number of photons in mode j. While the states $\{\Psi \}$ form a complete basis in the full Hilbert space, we assume that at any time the system can be described by exactly one of these states, i.e. we neglect the coherence of each TLS and coherent superpositions between several TLS. We also consider the radiation field as thermal, namely each state being diagonal in the eigenbasis of the non-interacting Hamiltonian

$$\begin{aligned} H_0=\frac{\hbar \Omega }{2}\sum _{l=1}^M\sigma ^z_l +\hbar \sum _j\omega _ja^{\dagger }_j a_j. \end{aligned}$$

(26)

The interaction between the TLS and the radiation via the Hamiltonian (1) is incorporated as a stochastic process using the golden rule. It gives the probabilities to transition from state $|\Psi (t)\rangle$ to states $|\Psi '(t+\delta t)\rangle$ (emission) or $|\Psi ''(t+\delta t)\rangle$ (absorption),

$$\begin{aligned}&|s_1,\ldots ,1_l,\ldots s_M;\{n_j\}\rangle \rightarrow |s_1,\ldots ,0_l,\ldots s_M;\{n_{j\ne k}\},n_k+1\rangle , \end{aligned}$$

(27)

$$\begin{aligned}&|s_1,\ldots ,0_{l'},\ldots s_M;\{n_j\}\rangle \rightarrow |s_1,\ldots ,1_{l'},\ldots s_M;\{n_{j\ne k}\},n_k-1\rangle . \end{aligned}$$

(28)

The rate $\tau _{\text {em}}^{-1}$ for the emission process of a single TLS depends on the occupation $n_j$ of the mode j actually involved. However, as the process couples the TLS to a continuum of radiation degrees of freedom, the total probablility is given by the average $\bar{n}$ of the occupations of all modes in the interval $[\Omega -\Delta /2,\Omega +\Delta /2]$ belonging to the state $|\Psi (t)\rangle$ as $\tau ^{-1}_\text {em}=\gamma '(\bar{n}(t)+1)$ (see [13] and  “Derivation of the Rate Equations (17)–(19) for the Open System”). As we don’t know into which mode the emission occurs, we can simplify the description of $|\Psi (t)\rangle$ given in (25) as

$$\begin{aligned} \Psi (t)= |s_1,\ldots s_M;\bar{n}\rangle , \end{aligned}$$

(29)

which already is a form of coarse-graining procedure applied to the microstate (25). The next and crucial approximation is the assumption that the radiation processes of the several TLS are independent, although the coupling to the radiation field introduces statistical correlations between them. The approximation is justified by the fact that we have a continuum of radiation modes which do not interact among themselves. This implies that each individual radiation process can be considered as independent from all others; the chance of two TLS to interact with exactly the same microscopic field mode j within observation time is vanishingly small. Thus we consider the joint probability $P^e_l(t,\bar{n})$ for the l-th TLS to be in the excited state and the average mode occupation to be $\bar{n}$. This probability satisfies the master equation [19]

$$\begin{aligned} \frac{\text {d}P_l^e(t,\bar{n})}{\text {d}t}= -\gamma 'P_l^e(t,\bar{n})(\bar{n} +1) +\gamma 'P_l^g(t,\bar{n})\bar{n}, \end{aligned}$$

(30)

where $P_l^g(t,\bar{n})=P(\bar{n},t)-P_l^e(t,\bar{n})$ is the probablility that the l-th TLS is in its ground state. This relation is a consequence of the factorization property

$$\begin{aligned} P_l^e(t,\bar{n})=P^e_l(t)P(\bar{n},t), \end{aligned}$$

(31)

which follows from the same argument as the statistical independence of the TLS, i.e. the presence of macroscopically many different modes of the radiation field in the frequency interval $[\Omega -\Delta ,\Omega +\Delta ]$, rendering the average occupation $\bar{n}$ statistically independent from the state of a single TLS. Summing (30) over $\bar{n}$, we obtain

$$\begin{aligned} \frac{\text {d}P_l^e(t)}{\text {d}t}= -\gamma 'P_l^e(t)(\langle n\rangle (t) +1) +\gamma '(1-P_l^e(t))\langle n\rangle (t). \end{aligned}$$

(32)

Noting that $P^e_l(t)=\langle m(t)\rangle /M$, we obtain Eq. (3) for the ensemble average of the number m(t) of excited TLS. The corresponding rate equation (2) for the average photon number follows directly from the conservation of the total excitation number

$$\begin{aligned} N_{\text {exc}}=\sum ^M_{l=1} s_l + \sum _j n_j, \end{aligned}$$

(33)

but contains a different rate constant $\gamma$ to account for the continuous density of states around the resonance frequency, $\gamma =\gamma '/\mathcal{N}$ (see Eq. (70) in  “Derivation of the Rate Equations (17)–(19) for the Open System”). The state vector (25) satisfies a much more complicated master equation then (32) and its stochastic evolution shows deviations of the actual value m(t) from its average $\langle m(t)\rangle$ computed via (2), (3). But one may argue that a typical state will have a time evolution which deviates only slightly from its mean value for all times, $|m_{\text {typical}}(t)-\langle m(t)\rangle |/M\rightarrow 0$ for $M\rightarrow \infty$ [16, 53]. This statement can be shown analytically for simple stochastic models like the Ehrenfest model [17] and also in purely deterministic toy models, e.g. the Kac ring model [55]. Moreover there is ample numerical evidence for such a behavior of typical states in more realistic models [18].

In our case, the ultimate justification of (2), (3) rests on the fact that they lead to Planck’s law for black-body radiation, one of the cornerstones of modern physics. It would be hard to argue that these equations are incorrect but in perfect agreement with all experiments due to some unknown cancellation of errors.

1.2 Derivation of the Rate Equations (17)–(19) for the Open System

We begin with the computation of the loss rate of photons with frequency $\omega _j$ from reservoir A via channel 1 of the open waveguide (see Fig. 1). In the coupling Hamiltonian, (11), the process is of first order. The initial state reads

$$\begin{aligned} |\psi _\text {in}\rangle = |\{n_A\}n_{Aj},\{n_B\},\{0_1\}0_{1k},\{0_2\},\{s\}\rangle , \end{aligned}$$

(34)

for a certain configuration of occupations $\{n_A\}$ for modes $j'$ in A with $j'\ne j$. Here, $n_{Aj}$ is the occupation number of mode j. Similarly, $\{n_B\}$ denotes a configuration in reservoir B. Further, $\{s\}$ is the configuration of excited ($s_l=e$) or ground states ($s_l=g$) of the M TLS, $l=1,\ldots M$. The waveguide channels 1 and 2 are not occupied in $\psi _\text {in}$.

The initial state is connected via the term $h_{jk}^*a_{1k}^\dagger a_{Aj}$ to the final state

$$\begin{aligned} |\psi _\text {fin}\rangle = |\{n_A\}n_{Aj}-1,\{n_B\},\{0_1\}1_{1k},\{0_2\},\{s\}\rangle . \end{aligned}$$

(35)

The S-matrix element reads then

$$\begin{aligned} S^{t_0}_{fi}=-2\pi i\delta ^{t_0}(E_\text {fin}- E_\text {in})\langle \psi _\text {fin}|H_{\text {int}}^1|\psi _\text {in}\rangle \end{aligned}$$

(36)

with $E_\text {fin}-E_\text {in}= \hbar (\omega _k-\omega _j)$. The interaction is assumed to take place over a time interval $t_0$. The corresponding regularized $\delta$-function is (see [13])

$$\begin{aligned} \delta ^{t_0}(E) = \frac{1}{\pi }\frac{\sin (Et_0/2\hbar )}{E}. \end{aligned}$$

(37)

We have

$$\begin{aligned} \langle \psi _\text {fin}|H_{\text {int}}^1|\psi _\text {in}\rangle = h^*_{jk}\sqrt{n_{Aj}}, \end{aligned}$$

(38)

and, according to standard reasoning [13], the transition rate is given by

$$\begin{aligned} \tau _\text {dec}^{-1}(\omega _j)=\frac{1}{t_0}\sum _k 4\pi ^2|h_{kj}|^2\frac{t_0}{2\pi \hbar }\delta ^{t_0}(E_\text {fin}-E_\text {in}) =\frac{2\pi }{\hbar }|\bar{h}|^2n_{Aj}\rho _1(\hbar \omega _j), \end{aligned}$$

(39)

where $\bar{h}$ is the value of $h_{jk}$ for $\omega _j=\omega _k$, which is assumed to be constant in the interval $[\Omega -\Delta /2,\Omega +\Delta /2]$. Averaging over all modes j with frequency $\omega$, we obtain $\tau _\text {dec}^{-1}(\omega )=\gamma _\text {dec}(\omega )\langle n_{A}\rangle (\omega )$ with

$$\begin{aligned} \gamma _\text {dec}(\omega )=\frac{2\pi }{\hbar }|\bar{h}|^2\rho _{wg}(\hbar \omega ), \end{aligned}$$

(40)

assuming $\rho _1=\rho _2=\rho _{wg}$. Adding the contribution of the decay into channel 2, the master equation for $\langle n_{A}\rangle (\omega )$ describing the loss process reads

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle (\omega )}{\text {d}t}=-2\gamma _\text {dec}(\omega )\langle n_{A}\rangle (\omega )+\ldots \end{aligned}$$

(41)

A similar expression holds for reservoir B with the same rate constant $\gamma _\text {dec}$.

The TLS couple to the channels via spontaneous emission. The corresponding first-order process leads to

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t}=-(\tilde{\gamma }_{11}(\Omega )+\tilde{\gamma }_{12}(\Omega ))\langle m\rangle +\ldots , \end{aligned}$$

(42)

with $\tilde{\gamma }_{1r}(\omega )=2\pi \bar{g}_r^2(\omega )\rho _{wg}(\hbar \omega )/\hbar$ for $r=1,2$ (see below).

We compute now the coherent transfer of a photon from reservoir A to B through the chiral waveguide, which is of second order in the coupling $|\bar{h}|^2$. As this process concerns a wavepacket of finite width and takes place in a finite time interval, only the right-moving channel 1 is relevant, i.e., the coupling Hamiltonian is

$$\begin{aligned} H^1_{\text {int}}=\sum _{q=A,B}\sum _{j,k} h_{jk}\left( a_{qj}^\dagger a_{1k} +\text {h.c.}\right) . \end{aligned}$$

(43)

We denote the initial state as

$$\begin{aligned} |\psi _\text {in}\rangle = |\{n_A\}n_{Aj},\{n_B\}n_{Bl},\{0_1\}0_{1k},\{0_2\},\{s\}\rangle . \end{aligned}$$

(44)

The transition from A to B occurs via the intermediate states

$$\begin{aligned} |\psi _\text {im}\rangle = |\{n_A\}n_{Aj}-1,\{n_B\}n_{Bl},\{0_1\}1_{1k},\{0_2\},\{s\}\rangle \end{aligned}$$

(45)

towards the final state

$$\begin{aligned} |\psi _\text {fin}\rangle = |\{n_A\}n_{Aj}-1,\{n_B\}n_{Bl}+1,\{0_1\}0_{1k},\{0_2\},\{s\}\rangle , \end{aligned}$$

(46)

involving the operators

$$\begin{aligned} h^*_{jk}a^{\dagger }_{1k}a_{Aj}, \quad h_{lk}a^{\dagger }_{Bl}a_{1k}. \end{aligned}$$

(47)

The S-matrix $S_{fi}$ connecting initial and final states reads

$$\begin{aligned} S^{t_0}_{fi}= -2\pi i \delta ^{t_0}(E_\text {fin}-E_\text {in})\lim _{\eta \rightarrow 0^+} \sum _{k}\frac{V_{fk}V_{ki}}{E_{\text {in}}-E_{\text {im}} +i\eta }, \end{aligned}$$

(48)

and $E_\text {in}-E_\text {im}= \hbar (\omega _j-\omega _k)$. The matrix elements $V_{fk}$ and $V_{ki}$ are

$$\begin{aligned} V_{fk}=\langle \psi _{\text {fin}}|h_{lk}a^{\dagger }_{Bl}a_{1k}|\psi _{\text {im}}\rangle , \quad V_{ki}=\langle \psi _{\text {im}}|h^*_{jk}a^{\dagger }_{1k}a_{Aj}|\psi _{\text {in}}\rangle . \end{aligned}$$

(49)

We find

$$\begin{aligned} S^{t_0}_{fi} \approx -2\pi ^2\delta ^{t_0}(\hbar (\omega _l-\omega _j)) |\bar{h}|^2\rho _{wg}(\hbar \omega _j) \sqrt{n_{Aj}(n_{Bl}+1)}, \end{aligned}$$

(50)

and obtain for the transition rate out of state $|\psi _\text {in}\rangle$,

$$\begin{aligned} \tau ^{-1}_{jA\rightarrow B}=\frac{1}{t_0}\sum _l|S^{t_0}_{fi}|^2= \sum _l \frac{2\pi ^3}{\hbar }\delta ^{t_0}(\hbar (\omega _l-\omega _j))|\bar{h}|^4\rho ^2_{wg}(\hbar \omega _j)n_{Aj}(n_{Bl}+1). \end{aligned}$$

(51)

We denote with $\langle n_{B}\rangle (\omega )$ the average occupation number per mode in reservoir B at frequency $\omega$. It follows that

$$\begin{aligned} \tau ^{-1}_{jA\rightarrow B}= \frac{2\pi ^3}{\hbar }|\bar{h}|^4\rho _{wg}^2(\hbar \omega _j) \rho _B(\hbar \omega _j)n_{Aj}(\langle n_{B}\rangle (\omega _j)+1). \end{aligned}$$

(52)

If we average over all modes j in A with frequency $\omega$, the transition rate from A to B at $\omega$ reads

$$\begin{aligned} \frac{1}{\tau _{A\rightarrow B}(\omega )}= \frac{2\pi ^3}{\hbar }|\bar{h}|^4\rho _{wg}^2(\hbar \omega ) \rho _B(\hbar \omega )\langle n_{A}\rangle (\omega )(\langle n_{B}\rangle (\omega )+1). \end{aligned}$$

(53)

Similarly, the transition from B to A, proceeding via channel 2, is

$$\begin{aligned} \frac{1}{\tau _{B\rightarrow A}(\omega )}= \frac{2\pi ^3}{\hbar }|\bar{h}|^4\rho _{wg}^2(\hbar \omega ) \rho _A(\hbar \omega )\langle n_{B}\rangle (\omega )(\langle n_{A}\rangle (\omega )+1). \end{aligned}$$

(54)

Finally, the master equation for reservoir A characterizing direct transitions between A and B through the waveguide is given by

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle (\omega )}{\text {d}t} = -\frac{1}{\tau _{A\rightarrow B}(\omega )} +\frac{1}{\tau _{B\rightarrow A}(\omega )} = \gamma _0(\omega )(\langle n_{B}\rangle (\omega ) - \langle n_{A}\rangle (\omega )), \end{aligned}$$

(55)

where we have assumed identical densities of states in A and B, $\rho _A=\rho _B=\rho$. The rate $\gamma _0(\omega )$ is defined as

$$\begin{aligned} \gamma _0(\omega )= \frac{2\pi ^3}{\hbar }|\bar{h}|^4\rho _{wg}^2(\hbar \omega ) \rho (\hbar \omega ). \end{aligned}$$

(56)

The ratio of second and first order contributions follows as

$$\begin{aligned} \frac{\gamma _0(\omega )}{\gamma _\text {dec}(\omega )}=\pi ^2|\bar{h}|^2\rho _{wg}(\hbar \omega )\rho (\hbar \omega ). \end{aligned}$$

(57)

Next, we consider the absorption of radiation from reservoir A by the l-th TLS. The initial state is

$$\begin{aligned} |\psi _\text {in}\rangle = |\{n_A\}n_{Aj},\{n_B\},\{0_1\}0_{1k},\{0_2\},\{s\}g_l\rangle . \end{aligned}$$

(58)

The absorption of a wavepacket of finite spatial extension can only proceed via channel 1. The intermediate states read then

$$\begin{aligned} |\psi _\text {im}\rangle = |\{n_A\}n_{Aj}-1,\{n_B\},\{0_1\}1_{1k},\{0_2\},\{s\}g_l\rangle \end{aligned}$$

(59)

and the final state is obtained by absorption of the photon in channel 1 by the TLS,

$$\begin{aligned} |\psi _\text {fin}\rangle = |\{n_A\}n_{Aj}-1,\{n_B\},\{0_1\}0_{1k},\{0_2\},\{s\}e_l\rangle . \end{aligned}$$

(60)

For the S-matrix element we have

$$\begin{aligned} S^{t_0}_{fi}=-2\pi ^2\bar{h}\bar{g}_1\delta ^{t_0}(\hbar \Omega -\hbar \omega _j)\rho _{wg}(\hbar \omega _j) \sqrt{n_{Aj}}, \end{aligned}$$

(61)

where $\bar{g}_1$ is the value of $g_{1k}$ for $\omega _k=\omega _j$. The corresponding term in the master equation for $\langle m\rangle$, the average number of excited TLS, is obtained by summing over all initial states, which leads to the expression

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t} = -\frac{2\pi ^3}{\hbar }(M-\langle m\rangle )\bar{h}^2\bar{g}^2_1\rho _{wg}^2(\hbar \Omega )\rho (\hbar \Omega )\langle n_{A}\rangle (\Omega ) +\ldots \end{aligned}$$

(62)

In an anlogous manner, the radiation from reservoir B is absorbed via channel 2 by the TLS,

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t} = -\frac{2\pi ^3}{\hbar }(M-\langle m\rangle )\bar{h}^2\bar{g}^2_2\rho _{wg}^2(\hbar \Omega )\rho (\hbar \Omega )\langle n_{B}\rangle (\Omega ) +\ldots \end{aligned}$$

(63)

On the other hand, the emission from the TLS towards reservoir A must proceed via the left-moving channel 2. The initial state is now

$$\begin{aligned} |\psi _\text {in}\rangle = |\{n_A\}n_{Aj},\{n_B\},\{0_1\},\{0_2\}0_{2k},\{s\}e_l\rangle , \end{aligned}$$

(64)

the intermediate state

$$\begin{aligned} |\psi _\text {im}\rangle = |\{n_A\}n_{Aj},\{n_B\},\{0_1\},\{0_2\}1_{2k},\{s\}g_l\rangle , \end{aligned}$$

(65)

and the final state

$$\begin{aligned} |\psi _\text {fin}\rangle = |\{n_A\}n_{Aj}+1,\{n_B\},\{0_1\},\{0_2\}0_{2k},\{s\}g_l\rangle . \end{aligned}$$

(66)

A calculation completely analogous to the one for absorption above leads to the following term in the master equation for $\langle m\rangle$, this time summing over final modes in reservoir A,

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t} = \frac{2\pi ^3}{\hbar }\langle m\rangle \bar{h}^2\bar{g}^2_2\rho _{wg}^2(\hbar \Omega )\rho (\hbar \Omega )(\langle n_{A}\rangle (\Omega )+1) +\ldots \end{aligned}$$

(67)

In this expression, the term proportional to $\langle n_{A}\rangle$ is associated with stimulated emission from the TLS. It is noteworthy that the stimulated emission of the TLS is not caused by photons that have been emitted by A and then impinge on the TLS to create photons in the same mode. Here, in contrast, the stimulated emission of the TLS is induced by photons that are received by A after having been emitted by the TLS. This counterintuitive effect, which is solely due to the Bose statistics and therefore only possible in quantum physics, restores detailed balance for the case of a cavity embedded in another one, in accord with Kirchhoff’s law on black body radiation (see “Detailed Balance for a System with an Embedded Cavity”).

With the definition

$$\begin{aligned} \tilde{\gamma }_r(\Omega )=\frac{2\pi ^3}{\hbar }\bar{h}^2\bar{g}^2_r\rho _{wg}^2(\hbar \Omega )\rho (\hbar \Omega ) \end{aligned}$$

(68)

for $r=1,2$ we obtain the master equation for $\langle m\rangle$,

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t}&= -[\tilde{\gamma }_{11}(\Omega ) +\tilde{\gamma }_{12}(\Omega )]\langle m\rangle + \tilde{\gamma }_1(\Omega )\left[ \langle n_{A}\rangle (M-\langle m\rangle )-(\langle n_{B}\rangle +1)\langle m\rangle \right] \nonumber \\&\quad +\tilde{\gamma }_2(\Omega )\left[ \langle n_{B}\rangle (M-\langle m\rangle )-(\langle n_{A}\rangle +1)\langle m\rangle \right] , \end{aligned}$$

(69)

which is Eq. (19) in the main text. To compute the terms corresponding to the absorption and emission processes in the rate equations for $\langle n_{A}\rangle$ and $\langle n_{B}\rangle$, we note that the rates $\tilde{\gamma }_r$ contain a summation over initial, respectively final modes in the reservoirs. The coefficients $\gamma _r$ describing the temporal change in the average occupation number per mode, $\langle n_q\rangle$, are therefore $\tilde{\gamma }_r$ divided by the number of modes in the relevant frequency interval. Thus,

$$\begin{aligned} \mathcal{N}= \hbar \int _{\Omega -\Delta /2}^{\Omega +\Delta /2}\text {d}\omega \rho (\hbar \omega ), \quad \gamma _r=\tilde{\gamma }_r/\mathcal{N}. \end{aligned}$$

(70)

Together with (41) and (55), it follows for $\omega =\Omega$,

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t}&= -2\gamma _\text {dec}\langle n_{A}\rangle + \gamma _0(-\langle n_{A}\rangle +\langle n_{B}\rangle )- \gamma _1(M-\langle m\rangle )\langle n_{A}\rangle \nonumber \\&\quad + \gamma _2\langle m\rangle (\langle n_{A}\rangle +1), \end{aligned}$$

(71)

$$\begin{aligned} \frac{\text {d}\langle n_{B}\rangle }{\text {d}t}&= -2\gamma _\text {dec}\langle n_{B}\rangle + \gamma _0(-\langle n_{B}\rangle +\langle n_{A}\rangle )- \gamma _2(M-\langle m\rangle )\langle n_{B}\rangle \nonumber \\&\quad + \gamma _1\langle m\rangle (\langle n_{B}\rangle +1), \end{aligned}$$

(72)

which are Eqs. (17) and (18).

1.3 Derivation of the Entropy Production (22)–(24)

As no work is done during the irreversible process, we can write for the heat flow out of the system within our approximation

$$\begin{aligned} Q^{\text {ext}}=-\frac{\text {d}}{\text {d}t}\big (\langle H_A\rangle +\langle H_B\rangle +\langle H_{\text {TLS}}\rangle \big ) >0, \end{aligned}$$

(73)

because the coherences of the TLS are not important and the waveguide contains no photons on average. We consider in the following only quantities (entropy and energy) corresponding to the interval $[\hbar (\Omega -\Delta /2),\hbar (\Omega +\Delta /2)]$ because other photon modes do not interact with the TLS. The heat transferred to the environment per unit time reads for each subsystem (see (17)–(19))

$$\begin{aligned} Q^{\text {ext}}_q=2\hbar \Omega {\mathcal{N}}\gamma _{\text {dec}}\langle n_q\rangle , \quad Q^{\text {ext}}_{\text {TLS}}=\hbar \Omega (\tilde{\gamma }_{11}+\tilde{\gamma }_{12})\langle m\rangle . \end{aligned}$$

(74)

The corresponding entropy change of the system is given by

$$\begin{aligned} -\int \text {d}{\varvec{o}}\cdot {\varvec{J}}_{\text {sys}} = -\left( \sum _{q=A,B}\frac{Q^{\text {ext}}_q}{T_q} + \frac{Q^{\text {ext}}_{\text {TLS}}}{T_{\text {TLS}}}\right) , \end{aligned}$$

(75)

with

$$\begin{aligned} T_q^{-1}=\frac{k_B}{\hbar \Omega }\ln \left( \frac{\langle n_q\rangle +1}{\langle n_q\rangle }\right) , \quad T_{\text {TLS}}^{-1}=\frac{k_B}{\hbar \Omega }\ln \left( \frac{M-\langle m\rangle }{\langle m\rangle }\right) . \end{aligned}$$

(76)

The entropy of the system is always reduced due to the vanishing temperature of the bath. Although the entropy production of system plus bath is thus formally positive and infinite, the actual entropy loss of the system stays finite and is of course compatible with the second law which mandates that the local entropy production within the system must be non-negative. In our case, the local entropy production is induced by the heat exchange between the three subsystems. The heat change of reservoir A due to the interaction with reservoir B reads

$$\begin{aligned} Q^A_{A,B}=-Q^B_{A,B}=\hbar \Omega \mathcal{N}\gamma _0(\langle n_{B}\rangle -\langle n_{A}\rangle ). \end{aligned}$$

(77)

It follows for the entropy production due to this process

$$\begin{aligned} \sigma _{A,B}=\hbar \Omega \mathcal{N}\gamma _0(\langle n_{B}\rangle -\langle n_{A}\rangle )(T_A^{-1}-T_B^{-1}), \end{aligned}$$

(78)

which is Eq. (22). The Eqs. (23), (24) for the heat exchange between the reservoirs and the TLS are deduced correspondingly. Alternatively, one may compute the entropy change in each subsystem directly with the expressions (14) for the TLS entropy and (96) for the entropy of the radiation modes. This shows that the entropy change of each subsystem is not accompanied by additional entropy production but is solely due to heat transfer between subsystems, a consequence of the local thermal equilibrium in A, B and the TLS.

1.4 Derivation of the Rate Equations for the Closed System

This chapter discusses a variant of the open system characterized in  “Derivation of the Rate Equations (17)–(19) for the Open System”. In this variant, which is a closed system, the waveguide satisfies periodic boundary conditions, corresponding to a loop of length L, where L is much larger than the distance between the reservoirs A and B, see Fig. 6. We show that a description assuming real absorption and emission processes fulfills the detailed balance condition to first order in the coupling. Second order processes analogous to those described in Eqs. (58)–(67) break the detailed balance condition for the case that coherent absorption and emission is only possible along the short path between the reservoirs and the TLS, while dephasing occurs for wave packets emitted by a TLS and traveling the long way around the loop before reaching one of the reservoirs. The latter case is already accounted for by the first-order terms describing the equilibration between the reservoirs/TLS and the waveguide. Of course, if the dynamics of the closed system is considered to be unitary, corresponding to completely coherent evolution, the entropy does not change. This situation could be approximated by treating all second order processes as coherent, including those on the long path between the TLS and the cavities. Then the detailed balance condition would be satisfied, leading to stabilization of the state with maximal entropy. However, if the coherence is restricted to processes occurring along the short path, the ensuing steady state does not have maximum entropy.

Fig. 6

The system obtained by closing the open waveguide of the system shown in Fig. 1 of the main text. The occupation numbers of the chiral channels 1 and 2 do not vanish due to equilibration with the reservoirs A and B

In the closed system, the occupancy of the waveguide can no longer be assumed to be zero, as the photons cannot escape towards infinity. We describe the occupancy in channel q by $\langle n_q\rangle$ for $q=1,2$. The Hamiltonian is given in Eqs. (8)–(12). The exchange of photons between cavity A and channel 1 of the waveguide proceeds via a first-order process analogous to that given in Eqs. (34) and (35), but the matrix element reads now (suppressing the frequency arguments),

$$\begin{aligned} \langle \psi _\text {fin}|H_{\text {int}}^1|\psi _\text {in}\rangle = h^*_{jk}\sqrt{n_{Aj}(n_1+1)}. \end{aligned}$$

(79)

This gives for the transition rate from mode j in A to channel 1 of the waveguide

$$\begin{aligned} \tau _{A\rightarrow 1}^{-1}=\gamma _{\text {dec}}n_{Aj}(\langle n_{1}\rangle +1), \end{aligned}$$

(80)

and for the reverse process,

$$\begin{aligned} \tau _{1\rightarrow A}^{-1}=\gamma _{\text {dec}}\langle n_{1}\rangle (n_{Aj} +1), \end{aligned}$$

(81)

where in this case a sum over initial states has to be performed to obtain the rate of emission into the fixed mode j of A. The terms in the master equation for $\langle n_{A}\rangle$ are thus

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t} = \gamma _{\text {dec}}(\langle n_{1}\rangle -\langle n_{A}\rangle ) + \ldots , \end{aligned}$$

(82)

and for $\langle n_{1}\rangle$

$$\begin{aligned} \frac{\text {d}\langle n_{1}\rangle }{\text {d}t} = \gamma _{\text {dec}}'(\langle n_{A}\rangle -\langle n_{1}\rangle ) + \ldots , \end{aligned}$$

(83)

with the rate constant $\gamma _\text {dec}'=(\rho /\rho _{wg})\gamma _\text {dec}$. Analogous expressions are obtained for B and channel 2. Another first-order process couples the TLS and the waveguide modes. We find for this contribution to the rate equation for $\langle m\rangle$

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t} = \sum _{r=1,2} \tilde{\gamma }_{1r}[(M-\langle m\rangle )\langle n_{r}\rangle -\langle m\rangle (\langle n_{r}\rangle +1)] + \ldots , \end{aligned}$$

(84)

with $\tilde{\gamma }_{1r}=2\pi g_r^2\rho _{wg}/\hbar$ (compare Eq. (42)). The corresponding terms in the rate equations for $\langle n_{r}\rangle$ are

$$\begin{aligned} \frac{\text {d}\langle n_{r}\rangle }{\text {d}t} = \gamma _{1r}[\langle m\rangle (\langle n_{r}\rangle +1) -(M-\langle m\rangle )\langle n_{r}\rangle ] + \ldots , \end{aligned}$$

(85)

and

$$\begin{aligned} \gamma _{1r}=\tilde{\gamma }_{1r}/\mathcal{N}', \quad \mathcal{N}' = \hbar \int _{\Omega -\Delta /2}^{\Omega +\Delta /2}\text {d}\omega \rho _{wg}(\hbar \omega ). \end{aligned}$$

(86)

All these first-order terms satisfy the detailed balance condition. They lead naturally to thermal equilibration between the reservoirs and the waveguide. For the second-order terms, we have first the process described by Eq. (69),

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t}&= \tilde{\gamma }_1(\langle n_{1}\rangle +1)\left[ \langle n_{A}\rangle (M-\langle m\rangle )-(\langle n_{B}\rangle +1)\langle m\rangle \right] \nonumber \\&\quad +\tilde{\gamma }_2(\langle n_{2}\rangle +1)\left[ \langle n_{B}\rangle (M-\langle m\rangle )-(\langle n_{A}\rangle +1)\langle m\rangle \right] + \ldots , \end{aligned}$$

(87)

which depends also on the occupation numbers $\langle n_{r}\rangle$ of the waveguide. This term is accompanied by corresponding terms in the rate equations for the reservoirs. It does not satisfy detailed balance because we have only considered the short path between the reservoirs and the TLS (the only available one in the open system). Including also the long path around the circular waveguide would again reinstate the detailed balance condition. We assume that this second process is not coherent due to dephasing of the photon while traveling along the loop. Such a dephasing may, for example, be caused by scattering processes induced in the long section of the loop. The second order term given in Eq. (55) is modified in the closed system as follows,

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t}=\gamma _0\big ( (\langle n_{2}\rangle +1)(\langle n_{A}\rangle +1)\langle n_{B}\rangle -(\langle n_{1}\rangle +1)(\langle n_{B}\rangle +1)\langle n_{A}\rangle \big ) +\ldots , \end{aligned}$$

(88)

together with an equivalent term for $\langle n_{B}\rangle$. Another term of second order in $|\bar{h}|^2$ couples the two channels of the waveguide via a reservoir. For channel 1 it reads

$$\begin{aligned} \frac{\text {d}\langle n_{1}\rangle }{\text {d}t} = \gamma _3(\langle n_{A}\rangle +\langle n_{B}\rangle +2)(\langle n_{2}\rangle -\langle n_{1}\rangle ) +\ldots , \end{aligned}$$

(89)

with $\gamma _3=2\pi ^3|\bar{h}|^4\rho ^2\rho _{wg}/\hbar$. Finally, there is a term connecting the channels via an intermediate excitation of the TLS, proportional to $g_1^2g_2^2$. This term can be neglected if one of the chiral couplings $g_r$ is close to zero, as we have assumed in the numerical evaluation. Collecting all the terms, we obtain the rate equations for the photon occupation numbers,

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t}&=\gamma _{\text {dec}}(\langle n_{1}\rangle +\langle n_{2}\rangle -2\langle n_{A}\rangle ) -\gamma _1(M-\langle m\rangle )\langle n_{A}\rangle (\langle n_{1}\rangle +1) \nonumber \\&\quad +\gamma _2\langle m\rangle (\langle n_{A}\rangle +1)(\langle n_{2}\rangle +1) +\gamma _0\big ( (\langle n_{2}\rangle +1)(\langle n_{A}\rangle +1)\langle n_{B}\rangle \nonumber \\&\quad -(\langle n_{1}\rangle +1)(\langle n_{B}\rangle +1)\langle n_{A}\rangle \big ), \end{aligned}$$

(90)

$$\begin{aligned} \frac{\text {d}\langle n_{B}\rangle }{\text {d}t}&=\gamma _{\text {dec}}(\langle n_{1}\rangle +\langle n_{2}\rangle -2\langle n_{B}\rangle ) -\gamma _2(M-\langle m\rangle )\langle n_{B}\rangle (\langle n_{2}\rangle +1) \nonumber \\&\quad \quad +\gamma _1\langle m\rangle (\langle n_{B}\rangle +1)(\langle n_{1}\rangle +1) +\gamma _0\big ( (\langle n_{1}\rangle +1)(\langle n_{B}\rangle +1)\langle n_{A}\rangle \nonumber \\&\quad -(\langle n_{2}\rangle +1)(\langle n_{A}\rangle +1)\langle n_{B}\rangle \big ), \end{aligned}$$

(91)

$$\begin{aligned} \frac{\text {d}\langle n_{1}\rangle }{\text {d}t}&=\gamma _\text {dec}'(\langle n_{A}\rangle +\langle n_{B}\rangle -2\langle n_{1}\rangle )-\gamma _{11}\big ((M-\langle m\rangle )\langle n_{1}\rangle - \langle m\rangle (\langle n_{1}\rangle +1)\big ) \nonumber \\&\quad +\gamma _3(\langle n_{A}\rangle +\langle n_{B}\rangle +2)(\langle n_{2}\rangle - \langle n_{1}\rangle ), \end{aligned}$$

(92)

$$\begin{aligned} \frac{\text {d}\langle n_{2}\rangle }{\text {d}t}&=\gamma _\text {dec}'(\langle n_{A}\rangle +\langle n_{B}\rangle -2\langle n_{2}\rangle )-\gamma _{12}\big ((M-\langle m\rangle )\langle n_{2}\rangle - \langle m\rangle (\langle n_{2}\rangle +1)\big ) \nonumber \\&\quad +\gamma _3(\langle n_{A}\rangle +\langle n_{B}\rangle +2)(\langle n_{1}\rangle - \langle n_{2}\rangle ). \end{aligned}$$

(93)

The average number of excited TLS is determined by

$$\begin{aligned} \frac{\text {d}\langle m\rangle }{\text {d}t}&= \sum _{r=1,2} \tilde{\gamma }_{1r}[(M-\langle m\rangle )\langle n_{r}\rangle -\langle m\rangle (\langle n_{r}\rangle +1)] \nonumber \\&\quad -\langle m\rangle \big (\tilde{\gamma }_1(\langle n_{B}\rangle +1)(\langle n_{1}\rangle +1)+\tilde{\gamma }_2(\langle n_{A}\rangle +1)(\langle n_{2}\rangle +1)\big ) \nonumber \\&\quad +(M-\langle m\rangle )\big (\tilde{\gamma }_1\langle n_{A}\rangle (\langle n_{1}\rangle +1)+\tilde{\gamma }_2\langle n_{B}\rangle (\langle n_{2}\rangle +1)\big ). \end{aligned}$$

(94)

Fig. 7

Solutions of the rate equations (90)–(94) as function of time, starting from initial thermal equilibrium. Panel a displays $\langle n_{A}\rangle (t)$, $\langle n_{B}\rangle (t)$ and panel b $\langle n_{1}\rangle (t)$, $\langle n_{2}\rangle (t)$. The average occupations $\langle n_q\rangle (t)$ per mode reach a novel steady state with $\langle n_{A}\rangle \ne \langle n_{B}\rangle$. The photon densities in reservoir A and channel 1 fall to zero, whereas channel 2 stays occupied and reservoir B is populated. The displayed time interval corresponds to the time scale set by the coupling to the TLS. Parameters used are $\gamma _\text {dec}=\gamma _\text {dec}'=\gamma _0=\gamma _3=10$ kHz, $\tilde{\gamma }_1=\tilde{\gamma }_{11}=10$ MHz, $\gamma _1=\gamma _{11}=100$ kHz, $\tilde{\gamma }_2=\tilde{\gamma }_{12}=\gamma _2=\gamma _{12}=0$ and $\mathcal{N}=\mathcal{N}'=100$, $\hbar \Omega /k_BT(0)=1$

To compute the total entropy of the system, we note first that

$$\begin{aligned} S_M=-k_BM\big (p_e\ln p_e +(1-p_e)\ln (1-p_e)\big ), \end{aligned}$$

(95)

is the entropy of the TLS, as given by Eq. (14). The entropy $S_\text {rad}(\omega )$ of the radiation per mode in the reservoirs and the waveguide depends only on $\langle n_q(\omega )\rangle$ and reads

$$\begin{aligned} S^l_{\text {rad}}(\omega ) =\frac{\hbar \omega }{T}\langle n_l(\omega )\rangle +k_B\ln (1+\langle n_l(\omega )\rangle ), \end{aligned}$$

(96)

for $l=A,B,1,2$. Because the temperature depends on $\langle n_l\rangle$ as described by Eq. (16) in the main text, the entropy is only a function of $\langle n_l\rangle$. The effective entropies and temperatures in the reservoirs and the waveguide are computed under the assumption that the photons in each reservoir/channel thermalize in the usual way quickly as an ideal Bose gas. We approximate the total entropy by the sum of the entropies of the subsystems, as is justified for large $\langle n_q\rangle$ and M. The total entropy is then given by

$$\begin{aligned} S(t)=S_M(t) + \int _{\Omega -\Delta /2}^{\Omega +\Delta /2}\text {d}\omega [\rho (\omega )S^R(\omega ,t)+\rho _{wg}(\omega )S^{wg}(\omega ,t)], \end{aligned}$$

(97)

with

$$\begin{aligned} S^R=S^A_\text {rad}+S^B_\text {rad}, \quad S^{wg}=S^1_\text {rad}+S^2_\text {rad}. \end{aligned}$$

(98)

For constant densities of states $\rho$ and $\rho _{wg}$, we may write

$$\begin{aligned} S(t) =S_M(t) + \mathcal{N}S^R(\Omega ,t) + \mathcal{N}' S^{wg}(\Omega ,t). \end{aligned}$$

(99)

The temporal evolution of the closed system, starting with thermal equilibrium, is depicted in Fig. 7 for the same parameters as in the open system. The non-reciprocal interaction with the TLS empties reservoir A and the active (coupled) channel 1, while the occupation of reservoir B rises. The inert channel 2 is unaffected on this short timescales. It interacts via the weak couplings $\gamma _\text {dec}'$ and $\gamma _3$ with B, which manifests only on much longer timescales, as depicted in Fig. 8. This separation of timescales has the same origin in the closed and the open system.

Fig. 8

The asymptotic temporal behavior of reservoir B (panel a) and the inert channel 2 (panel b). The occupation in A and channel 1 is almost zero. The novel steady state has unequal occupations in all photonic subsystems and therefore lower entropy than the initial state. The parameters used are the ones of Fig. 7

1.5 Detailed Balance for a System with an Embedded Cavity

We shall now demonstrate that the term in (67) corresponding to radiation stimulated by the receiving reservoir leads to the detailed balance condition in case a cavity is embedded into another one. Here, detailed balance is also obtained in the second-order terms, in contrast to the chiral system treated in the previous section. We consider a closed cavity C with adiabatic walls. Inside of C there is a smaller cavity A which is coupled to C through a small opening. Besides A, a collection of M two-level systems is located in C (Fig. 9).

Fig. 9

The cavity system. The small cavity A exchanges radiation with the surrounding closed cavity C. The collection M of two-level systems interacts with the radiation modes of C

The Hamiltonian of this system is given as

$$\begin{aligned} H=H_A + H_C + H_{\text {TLS}} + H^1_{\text {int}} + H^2_{\text {int}}, \end{aligned}$$

(100)

where

$$\begin{aligned} H_q=\hbar \sum _{j}\omega _{qj} a_{qj}^\dagger a_{qj}, \qquad H_{\text {TLS}}=\frac{\hbar \Omega }{2}\sum _{l=1}^M\sigma ^z_l \end{aligned}$$

(101)

for $q=A,C$. The interaction between A and C is given by

$$\begin{aligned} H^1_{\text {int}}=\sum _{j,k}h_{jk}a^\dagger _{Aj}a_{Ck} + \text {h.c.}, \end{aligned}$$

(102)

and C interacts with M as

$$\begin{aligned} H^2_{\text {int}}=\sum _{l=1}^M \sum _{k} g_{k}a_{Ck}\sigma ^+_l + \text {h.c.} \end{aligned}$$

(103)

The rate equations are computed as above, but now the exchange between A and C is given by terms of first order in the coupling $h_{jk}$,

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t}&= \gamma _4\left[ -\langle n_{A}\rangle (\langle n_{C}\rangle +1) + \langle n_{C}\rangle (\langle n_{A}\rangle +1)\right] +\ldots , \end{aligned}$$

(104)

$$\begin{aligned} \frac{\text {d}\langle n_{C}\rangle }{\text {d}t}&= \gamma _4'\left[ -\langle n_{C}\rangle (\langle n_{A}\rangle +1) + \langle n_{A}\rangle (\langle n_{C}\rangle +1)\right] +\ldots \end{aligned}$$

(105)

with $\gamma _4'=(\rho _A/\rho _C)\gamma _4$, (compare Eq. (83)). The interaction between the TLS and C leads to the terms

$$\begin{aligned} \frac{\text {d}\langle n_{C}\rangle }{\text {d}t} = -\gamma _5\langle n_{C}\rangle (M-\langle m\rangle ) +\gamma _5(\langle n_{C}\rangle +1)\langle m\rangle +\ldots , \end{aligned}$$

(106)

which are of first order in $g_{k}$ and correspond to standard black-body radiation. Besides these first-order terms, there are also terms of second order in the couplings $h_{jk}$ and $g_k$, describing the interaction of the small cavity A with the TLS via intermediate states belonging to C. However, in contrast to the second-order terms discussed in “Derivation of the Rate Equations for the Closed System”, these second-order terms are compatible with detailed balance. The corresponding terms in the rate equation for A read

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t} = \gamma _6(\langle n_{C}\rangle +1)^2\left[ -\langle n_{A}\rangle (M-\langle m\rangle ) +(\langle n_{A}\rangle +1)\langle m\rangle \right] +\ldots \end{aligned}$$

(107)

As above, the term for stimulated emission, $\gamma _6(\langle n_{C}\rangle +1)^2\langle n_{A}\rangle \langle m\rangle$, is not related to radiation emerging from cavity A into C which would have the wrong direction (see Fig. 9), but comes from the occupation of A-modes in the final state. The rate equations for A and C are therefore

$$\begin{aligned} \frac{\text {d}\langle n_{A}\rangle }{\text {d}t}&= \gamma _4\left[ -\langle n_{A}\rangle + \langle n_{C}\rangle \right] + \gamma _6(\langle n_{C}\rangle +1)^2\left[ -\langle n_{A}\rangle (M-\langle m\rangle ) +(\langle n_{A}\rangle +1)\langle m\rangle \right] , \end{aligned}$$

(108)

$$\begin{aligned} \frac{\text {d}\langle n_{C}\rangle }{\text {d}t}&= \gamma _4'\left[ -\langle n_{C}\rangle + \langle n_{A}\rangle \right] -\gamma _5\langle n_{C}\rangle (M-\langle m\rangle ) +\gamma _5(\langle n_{C}\rangle +1)\langle m\rangle . \end{aligned}$$

(109)

These equations fulfill the detailed balance condition and lead to thermal equilibrium between A, C and M. The result is consistent with Kirchhoff’s law, which states that the interior of a thermally isolated hohlraum has no influence on the final steady state of the contained radiation. This radiation exhibits the black-body spectrum found by M. Planck.