Introduction

When a physical system is simultaneously subjected to both predictable and random energy exchanges, what dictates the likelihood of a given state to be found? From a classical thermodynamic perspective, where energy exchanges are classified as work or heat, the concept of dissipative adaptation has recently been put forward by J. England as the expected answer, at least in the context of driven self-assembly1. Qualitatively speaking, dissipative adaptation establishes that, given a certain drive (an external work source), the most adapted (most likely and lasting) states of a fluctuating system tend to be those with a history of exceptionally high work absorption followed by heat dissipation to the environment. Because heat dissipation is an irreversible process, the higher the dissipation, the less likely the reverse trajectory is. In the long run, the system may appear to us as self-organized in this drive-dependent state of highest energy-consuming history. Dissipative adaptation is the recent theoretical development of a long search for the emergence of order from disorder, as inspired by life-like behavior2,3. Examples revealing this general mechanism of energy-consuming irreversible self-organization span diverse systems, environments, lengths and timescales, as shown both theoretically4,5,6 and experimentally7,8,9,10,11,12.

The dissipative adaptation phenomenon has been originally formulated in terms of fluctuation theorems. Fluctuation theorems are equalities relating out-of-equilibrium processes with thermal-equilibrium variables, giving evidence that the fluctuations present in many realizations of a process can provide useful knowledge when summed up. The pioneering example is the so-called Jarzynski equality13, \(\left\langle \exp (-\beta {W}_{{\rm{abs}}})\right\rangle =\exp (-\beta {{\Delta }}F)\), where β = 1/(kBT) is the inverse temperature, Wabs is the work absorbed by the system as described by a time-dependent Hamiltonian, ΔF is the variation in the Helmholtz free energy and the brackets is the ensemble average over realizations of the process, initially departing from thermal equilibrium. The Jarzynski equality has been lately derived by Crooks from what he called a microscopically reversible condition14,

$$\frac{{p}_{i\to j}(t)}{{p}_{j\to i}^{* }(\tau -t)}={e}^{\beta {Q}_{{\rm{diss}}}},$$
(1)

where the forward \({p}_{i\to j}(t)\) and backward \({p}_{j\to i}^{* }(\tau -t)\) probabilities for the system to follow paths linking states i and j are related with the amount of heat stochastically dissipated to the environment, Qdiss (a functional of the phase-space trajectory and of the driving protocol performing work on the system). \({p}_{j\to i}^{* }(\tau -t)\) is computed with the reversed time-dependent protocol. βQdiss relates the statistical irreversibility with the thermodynamic entropy production. Let us call Eij = Ej − Ei the energy difference between the specific final and initial states. Energy conservation during each stochastic realization,

$${W}_{{\rm{abs}}}={Q}_{{\rm{diss}}}+{E}_{ij},$$
(2)

gives the hint on how the work source fuels the dissipative adaptation. The higher the absorbed work, the more heat can be released, hence the more irreversible the path can become. To emphasize this adaptation-work relation, England rearranges Eq. (1) as1

$$\frac{{p}_{i\to j}(t)}{{p}_{i\to k}(t)}={e}^{-\beta {E}_{kj}}\frac{{p}_{j\to i}^{* }(\tau -t)}{{p}_{k\to i}^{* }(\tau -t)}\frac{{\left\langle {e}^{-\beta {W}_{{\rm{abs}}}}\right\rangle }_{ik}}{{\left\langle {e}^{-\beta {W}_{{\rm{abs}}}}\right\rangle }_{ij}},$$
(3)

where the angle brackets denote a weighted average over all microtrajectories with fixed start (i) and end (j, k) points (the coarse graining over microtrajectories describes experimentally accessible states of out-of-equilibrium self-organizing classical many-body systems; see Perunov et al.15 for further details). Equation (3) establishes the classical theoretical framework behind dissipative adaptation by showing how a given final state can be statistically privileged by work consumption.

Here, we extend the concept of dissipative adaptation to the quantum realm. Our main goal is to test the robustness of the dissipative adaptation concept beyond its original theoretical framework discussed above. From a technical viewpoint, at vanishing temperatures (β → ∞), where quantum fluctuations and correlations usually prevail, Eqs. (1) and (3) are ill defined. We employ a system-plus-reservoir approach to derive the exact equations of motion of a three-level lambda (Λ) system driven by a single-photon pulse added to a zero-temperature environment. We find that the most adapted (self-organized) quantum state of the lambda system is indeed that with a history of highest work absorption followed by maximal heat dissipation, thus characterizing a dissipative adaptation. As a consequence of our work, we establish the starting point of a quantum thermodynamics of driven self-organization, so far unexplored to the best of our knowledge. We hope that the notion of a quantum dissipative adaptation may also provide fresh insights to quantum biology16, not only because adaptation and self-organization are concepts inspired in life-like behavior, but also because our results may find applications to discussions on quantum signatures in photosynthesis17,18,19,20.

Results

Self-organized quantum state

To achieve our main goal, we look for the most elementary scenario where the quantum state of a certain physical system irreversibly self-organizes, as induced by the absorption of energy from an external drive, the excess of which is dissipated to the environment. We choose a three-level system in Λ configuration, labeled as \(\left|a\right\rangle\), \(\left|b\right\rangle\), and \(\left|e\right\rangle\) (being \(\left|e\right\rangle\) the most excited state, respectively, with transition frequencies ωa,b) (see Fig. 1). To keep the model exactly solvable for both the system and the environment, we assume that the drive source is a single-photon pulse added to the vacuum state of an infinite environment at zero temperature, T = 0. The environment induces spontaneous emission rates, Γa and Γb. Our main results are that

$${p}_{a\to b}(\infty ) =\frac{{{{\Gamma }}}_{b}}{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}\,\frac{{\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}}{\hslash {\omega }_{a}},\,{\rm{and}}\\ {p}_{b\to a}^{* }(t) =0,$$
(4)

showing that adaptation likelihood at long times, pa→b(t → ∞), is linearly proportional to the average absorbed work \({\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}\) from a single-photon pulse of arbitrary shape, resonant with ωa. We call attention to the fact that this result is not immediately expected, since (i) the excitation probability is minimized (pa→e(t) ≪ 1, at all times) when the work absorption is maximized (pa→b(∞) → 1), so the process cannot be thought of as a simple absorption-plus-emission picture, and (ii) the final self-organized state is not restrained to be the ground state, which amounts to say that \({\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}\) does not depend on ωb, so the work is not related to ℏ(ωa−ωb), in particular. The second line in Eqs. (4) characterizes the irreversibility of the process. The absorbed work is partly dissipated to the environment in the form of heat, \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\), satisfying energy conservation

$${\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}={\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}+\left\langle {H}_{S}(\infty )\right\rangle -\left\langle {H}_{S}(0)\right\rangle .$$
(5)

HS here is the Hamiltonian of the system. Finally, we find an exact expression for the informational entropy of the environment at long times, SE(∞), as a function of the average dissipated heat, \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\). The entropy analysis here provides us with both a clearer physical picture of the process and an additional signature of the dissipative adaptation.

Fig. 1: Quantum dissipative adaptation in a driven self-organized quantum state.
figure 1

a Three-level system in lambda configuration is described by the density operator ρS(t) at time t, starting at t = 0. The energy eigenstates are \(\left|a\right\rangle\), \(\left|b\right\rangle\), and \(\left|e\right\rangle\) (blue horizontal bars). The transition frequencies are ωa,b (full gray arrows). The environment-induced spontaneous emission rates are Γa,b (dashed gray arrows). The initial state is a mixture between \(\left|a\right\rangle\) and \(\left|b\right\rangle\), \({\rho }_{S}(0)={p}_{a}^{(0)}\left|a\right\rangle \left\langle a\right|+{p}_{b}^{(0)}\left|b\right\rangle \left\langle b\right|\) (smaller black dots). b A single-photon pulse \(\left|{1}_{a}\right\rangle\) (the work source) drives the lambda system dynamics, inducing the time-dependent transition probability pa→b(t) from \(\left|a\right\rangle\) to \(\left|b\right\rangle\) (the full black horizontal arrow represents the forward dynamics from time t = 0 to t → ∞). The backward probability (with a time-reversed pulse), \({p}_{b\to a}^{* }(t)=0\), vanishes at zero temperature (the dashed gray horizontal arrow represents the prohibited time-reversed transition). c The driven lambda system undergoes an ideal irreversible self-organizing dynamics (pa→b(∞) = 1 and \({p}_{b\to a}^{* }(t)=0\), in b), so that the asymptotic state is pure, \({\rho }_{S}(\infty )=\left|b\right\rangle \left\langle b\right|\) (larger black dot), conditioned to maximizing the work absorbed and the heat dissipated in the \(\left|a\right\rangle\) to \(\left|b\right\rangle\) transition.

Full size image

Let us suppose that our generic three-level system in Λ configuration, with lowest-energy states \(\left|a\right\rangle\) and \(\left|b\right\rangle\) and excited state \(\left|e\right\rangle\), is initially in a non-driven steady state, in contact with the environment at temperature T = 0. In this case, at precisely T = 0, the steady state is not uniquely defined, even for non-degenerated cases. We can choose its initial quantum state to be a mixture of the lowest-energy states, as described by the density operator

$${\rho }_{S}(0)={p}_{a}^{(0)}\left|a\right\rangle \left\langle a\right|+{p}_{b}^{(0)}\left|b\right\rangle \left\langle b\right|,$$
(6)

where \({p}_{a,b}^{(0)}\) depend on the preparation scheme. To give a concrete example, in the preparation by means of a spontaneous emission starting at \(\left|e\right\rangle\), one has \({p}_{a,b}^{(0)}={{{\Gamma }}}_{a,b}/({{{\Gamma }}}_{a}+{{{\Gamma }}}_{b})\). Now we look for the most elementary out-of-equilibrium stochastic process that drives the system from this (generally mixed) initial state into a final (ideally pure) target state, let us say into state \(\left|b\right\rangle\),

$${\rho }_{S}(0)\to {\rho }_{S}(\infty )=\left|b\right\rangle \left\langle b\right|.$$
(7)

In order to guarantee that the process is irreversible, in the light of Crooks condition, we should also apply the time-reversed drive on the system departing from state \({\rho }_{S}(\infty )=\left|b\right\rangle \left\langle b\right|\) and find that \({\rho }_{S}^{* }(\infty -t)\,\ne\, {\rho }_{S}(0)\) for t→∞.

We employ a system-plus-reservoir approach, where we assume a global time-independent Hamiltonian of the system and its environment,

$$H={H}_{S}+{H}_{I}+{H}_{E}.$$
(8)

As we show in what follows, we find a self-organized quantum state in the well-known dipolar model of light-matter interaction in the rotating-wave approximation21,22,

$${H}_{I}=-i\hslash \sum _{\omega }({g}_{a}{a}_{\omega }{\sigma }_{a}^{\dagger }+{g}_{b}{b}_{\omega }{\sigma }_{b}^{\dagger }-{\rm{H.c.}}).$$
(9)

Here, \({\sigma }_{k}=\left|k\right\rangle \left\langle e\right|\) (for k = a, b) and H.c. is the Hermitian conjugate. We consider a continuum of frequencies, ∑ω→∫dωϱω, with density of modes ϱω. Modes \(\left\{{a}_{\omega }\right\}\) and \(\left\{{b}_{\omega }\right\}\) are the quantized field modes interacting with the transitions \(\left|a\right\rangle\) to \(\left|e\right\rangle\) and \(\left|b\right\rangle\) to \(\left|e\right\rangle\), respectively. The continuum of frequencies allows us to employ a Wigner–Weisskopf approximation to obtain dissipation rates \({{{\Gamma }}}_{a}=2\pi {g}_{a}^{2}{\varrho }_{{\omega }_{a}}\) and \({{{\Gamma }}}_{b}=2\pi {g}_{b}^{2}{\varrho }_{{\omega }_{b}}\). Finally, \({H}_{E}={\sum }_{\omega }\hslash \omega ({a}_{\omega }^{\dagger }{a}_{\omega }+{b}_{\omega }^{\dagger }{b}_{\omega })\) and

$${H}_{S}=\hslash {\omega }_{a}\left|e\right\rangle \left\langle e\right|+\hslash {\delta }_{ab}\left|b\right\rangle \left\langle b\right|,$$
(10)

where ℏδab = ℏ(ωa−ωb) is the energy difference between states \(\left|b\right\rangle\) and \(\left|a\right\rangle\). It is worth emphasizing that our main results in this paper are independent of δab.

To keep the model exactly solvable for both the system and the environment, we choose the drive as provided by a propagating pulse containing a single photon. We choose the photon to initially populate only the continuum of modes \(\left\{{a}_{\omega }\right\}\), so the vacuum state of \(\left\{{b}_{\omega }\right\}\) allows us to avoid depleting our target state \(\left|b\right\rangle\). The initial state of the field is

$$\left|{1}_{a}\right\rangle =\sum _{\omega }{\phi }_{\omega }^{a}(0)\,{a}_{\omega }^{\dagger }\left|0\right\rangle .$$
(11)

\(\left|0\right\rangle ={\prod }_{\omega }\left|{0}_{\omega }^{a}\right\rangle \otimes \left|{0}_{\omega }^{b}\right\rangle\) is the vacuum state of all the field modes. The initial state of the Λ system is the mixed ρS(0) given in Eq. (6). The global quantum state is given by

$$\rho (t)=U\left({\rho }_{S}(0)\otimes \left|{1}_{a}\right\rangle \left\langle {1}_{a}\right|\right){U}^{\dagger },$$
(12)

where \(U=\exp (-iHt/\hslash )\) for the (time-independent) global Hamiltonian H [Eq. (8)]. The quantum states of the system and the environment are obtained by the partial traces \({\rho }_{S}(t)={{\rm{tr}}}_{E}[\rho (t)]\) and \({\rho }_{E}(t)={{\rm{tr}}}_{S}[\rho (t)]\).

We obtain analytical expressions for the probabilities \({p}_{k}(t)=\left\langle k\right|{\rho }_{S}(t)\left|k\right\rangle\), for k = a, b, e. Equation (6) allows us to write the reduced state of the system as \({\rho }_{S}(t)={p}_{a}^{(0)}{{\rm{tr}}}_{E}\left[U\left|a,{1}_{a}\right\rangle \left\langle a,{1}_{a}\right|{U}^{\dagger }\right]+{p}_{b}^{(0)}{{\rm{tr}}}_{E}\left[U\left|b,{1}_{a}\right\rangle \left\langle b,{1}_{a}\right|{U}^{\dagger }\right]\). It proves useful to write pk(t) in terms of transition probabilities, pn→k(t), as

$${p}_{k}(t)=\sum _{n}{p}_{n}^{(0)}{p}_{n\to k}(t),$$
(13)

where n = a, b, and

$${p}_{n\to k}(t)\equiv \left\langle k\right|{{\rm{tr}}}_{E}\left[U\left|n,{1}_{a}\right\rangle \left\langle n,{1}_{a}\right|{U}^{\dagger }\right]\left|k\right\rangle .$$
(14)

Equations (13) and (14), which follow directly from the mixed-state structure of the initial state, unravel a close analogy with the notation used in Eq. (1) for the microscopically reversible condition. Equation (13) has no relation with the degree of Markovianity in the dynamics of the reduced state ρS(t), though. On the contrary, Valente et al.23 characterize non-Markovianity in a quite similar scenario. The photon initially at modes \(\left\{{a}_{\omega }\right\}\) does not interact with the three-level system initially at \(\left|b\right\rangle\) [as derived from Eqs. (6), (9), and (12)], so we find that pb→a(t) = pb→e(t) = 0 and pb→b(t) = 1. Because the transition probabilities vanish identically, regardless of the initial pulse shape, the backwards probability also vanishes,

$${p}_{b\to a}^{* }(t)=0.$$
(15)

We have defined the reverse drive protocol here as the mirrored shape of the initial pulse (further discussed below). Put simply, we revert only the drive, not the entire universe. The final global state of the system plus the environment is obviously reversible in our model since it is governed by the global unitary operator U in Eq. (12). Equation (15) explains the second line in Eqs. (4). We now compute

$${p}_{a\to k}(t)=\left\langle k\right|{{\rm{tr}}}_{E}\left[\left|\xi (t)\right\rangle \left\langle \xi (t)\right|\right]\left|k\right\rangle ,$$
(16)

where \(\left|\xi (t)\right\rangle \equiv U\left|a,{1}_{a}\right\rangle\). Since our H conserves the total number of excitations, we can restrict our model to the one-excitation subspace, \(\left|\xi (t)\right\rangle =\psi (t)\left|e,0\right\rangle +{\sum }_{\omega }{\phi }_{\omega }^{a}(t){a}_{\omega }^{\dagger }\left|a,0\right\rangle +{\phi }_{\omega }^{b}(t){b}_{\omega }^{\dagger }\left|b,0\right\rangle .\) We find that

$${p}_{a\to e}(t) =| \psi (t){| }^{2},\\ {p}_{a\to b}(t) =\sum _{\omega }| {\phi }_{\omega }^{b}(t){| }^{2}\\ =\frac{1}{2\pi \varrho c}\mathop{\int}\nolimits_{\!\!\!\!-\infty }^{\infty }dz\,| {\phi }_{b}(z,t){| }^{2},$$
(17)

and similarly with pa→a(t). Without loss of generality, we have defined a one-dimensional real-space representation for the amplitudes, \({\phi }_{n}(z,t)\equiv {\sum }_{\omega }{\phi }_{\omega }^{n}(t)\exp (i{k}_{\omega }z)\), which characterizes the pulse shape. We have also employed a linear dispersion relation, ω = ckω, and approximated the density of modes by a constant, ϱω ≈ ϱ. Note that, in Eq. (9), we have implicitly assumed the three-level system to be positioned at zS = 0; otherwise, a phase term such as \(\exp (i{k}_{\omega }{z}_{S})\) should have been included in the sum over modes.

To go a step further, as to obtain explicit expressions for the transition probabilities, we solve the Schrödinger equation for \(\left|\xi (t)\right\rangle\) (see Methods). Our intention with keeping modes \(\left\{{b}_{\omega }\right\}\) initially in the vacuum state, as we did in Eq. (11), was to minimize excitations promoting the unwanted backward (\(\left|b\right\rangle \to \left|a\right\rangle\)) transitions. Now that we have defined the amplitudes ϕn(z, t), we see that Eq. (11) implies ϕb(z, 0) = 0, which we combine with (34) and substitute in (17). After a change of variables, we find that

$${p}_{a\to b}(t)={{{\Gamma }}}_{b}\mathop{\int}\nolimits_{\!\!\!\!|0}^{t}dt^{\prime} | \psi (t^{\prime} ){| }^{2}.$$
(18)

Although the main results in this paper do not depend on the initial pulse shape in modes \(\left\{{a}_{\omega }\right\}\) (i.e., on the choice of ∣ϕa(z, 0)∣), it is worth working out an explicit example. To that end, we now set ϕa(z, 0) to have an exponential envelope profile of linewidth Δ and a central frequency ωL (see Eqs. (36) and (37) in the Methods), as typical in spontaneous emission. We are particularly interested in the resonant condition ωL = ωa (see the Methods for the more general solution). We find, in the monochromatic limit, Δ ≪ Γa + Γb, and at long times t ≫ Δ−1, that

$${p}_{a\to b}(\infty )={{{\Gamma }}}_{b}\,\frac{4{{{\Gamma }}}_{a}}{{({{{\Gamma }}}_{a}+{{{\Gamma }}}_{b})}^{2}}.$$
(19)

When Γa = Γb, we have that pa→b(∞) = 1. Equation (19) reveals the ideal driven self-organization we have been looking for. As stated before, the self-organization in our model does not depend on δab. As a final step, we combine Eq. (19) and pb→b(t) = 1 to see that, in the ideal (monochromatic and resonant) regime,

$${p}_{b}(\infty )={p}_{a}^{(0)}{p}_{a\to b}(\infty )+{p}_{b}^{(0)}{p}_{b\to b}(\infty )=1,$$
(20)

that is, \({\rho }_{S}(\infty )=\left|b\right\rangle \left\langle b\right|\) and \({\rho }_{S}^{* }(\infty -t)=\left|b\right\rangle \left\langle b\right|\) for all times t. Hence, in the t → ∞ limit, \({\rho }_{S}^{* }(\infty -t)\,\ne\, {\rho }_{S}(0)\), as was to be shown.

Energetics of the self-organization

We now need to verify that the self-organized quantum state we have found can indeed be classified as dissipative adaptation. We shall find that ideal self-organization costs maximal work absorption, followed by maximal dissipation. This is not an obvious relation because in the ideal self-organization (which takes place in the monochromatic limit, as we have shown above), the excitation probability is minimized (rather than maximized),

$${p}_{a\to e}(t)=| \psi (t){| }^{2}\ll 1.$$
(21)

For instance, in the case of an incoming pulse with exponential profile (as used in Eq. (19)), \(| \psi (t){| }^{2}\le 4{{\Delta }}{{{\Gamma }}}_{a}/{({{{\Gamma }}}_{a}+{{{\Gamma }}}_{b})}^{2}\ll 1\), for Δ ≪ Γa + Γb. Extremely low-excitation probabilities may leave the false impression that no energy is absorbed neither dissipated at all. To address this issue, we need to resolve energy transfer into work and heat in our model. As mentioned earlier, work and heat can be regarded as the predictable versus the random energy exchanges, respectively. We restrict ourselves once again to the resonant case, ωL = ωa, keeping in mind however an arbitrary pulse shape. The average energy of the Λ system driven by this resonant photon is given by

$$\left\langle {H}_{S}(t)\right\rangle \equiv {\rm{tr}}[{\rho }_{S}(t){H}_{S}]={p}_{e}(t)\hslash {\omega }_{a}+{p}_{b}(t)\hslash {\delta }_{ab}.$$
(22)

The resonant condition here avoids dynamic Stark shifts that would otherwise bring extra (dispersive) energetic contributions from time-dependent frequencies, as shown by Valente et al.24. To be more precise, dispersive (refractive) energetic contributions depend on the average interaction energy24, which vanishes at resonance (ωL = ωa implies that \(\left\langle {H}_{I}(t)\right\rangle \equiv {\rm{tr}}[\rho (t){H}_{I}]=0\); see details in the Methods). That justifies why we have used only HS in Eq. (22). As inspired by Eq. (3), we are interested in the work performed on the system during the dynamical transition starting from \(\left|a\right\rangle\) and arriving at \(\left|b\right\rangle\). Therefore, we now consider \({\rho }_{S}(0)=\left|a\right\rangle \left\langle a\right|\), for which pe(t) = ∣ψ(t)∣2. From the Schrödinger equation (see Eq. (32) in the Methods), we find that

$${\partial }_{t}{p}_{e}=-({{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}){p}_{e}-2{g}_{a}{\rm{Re}}\left[{\phi }_{a}(-ct,0){\psi }^{* }(t)\right],$$
(23)

where \({\rm{Re}}\left[\bullet \right]\) stands for the real part. Equation (23) clearly shows us that the excited state of the Λ system is governed by a (predictable) drive-dependent term, related to ϕa(−ct, 0), and a (random) spontaneous emission term, proportional to Γa + Γb. This motivates us to define the total average absorbed work \({\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}\) and the total average absorbed heat \({\left\langle {Q}_{{\rm{abs}}}\right\rangle }_{a}\) in the dynamics starting at \({\rho }_{S}(0)=\left|a\right\rangle \left\langle a\right|\) as

$$\left\langle {H}_{S}(\infty )\right\rangle -\left\langle {H}_{S}(0)\right\rangle ={\left\langle {Q}_{{\rm{abs}}}\right\rangle }_{a}+{\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a},$$
(24)

where

$${\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}\equiv \hslash {\omega }_{a}\mathop{\int}\nolimits_{0}^{\infty }-2{g}_{a}{\rm{Re}}[{\phi }_{a}(-ct,0){\psi }^{* }(t)]dt$$
(25)

and

$${\left\langle {Q}_{{\rm{abs}}}\right\rangle }_{a}\equiv \mathop{\int}\nolimits_{0}^{\infty }-({{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}){p}_{e}(t)\hslash {\omega }_{a}dt+{\dot{p}}_{b}\hslash {\delta }_{ab}dt.$$
(26)

Our definition of work, Eq. (25), can also be written as \({\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}=\mathop{\int}\nolimits_{0}^{\infty }\left\langle \left({\partial }_{t}d(t)\right){E}_{{\rm{in}}}(t)\right\rangle dt\) (in the Heisenberg picture and in the rotating-wave approximation), revealing a more clear link to our classical notion of work. Here, d(t) = U†dU is the dipole operator, where d = ∑kdekσk + H.c., the driving field operator at the system’s position is \({E}_{{\rm{in}}}(t)={\sum }_{\omega }i{\epsilon }_{\omega }\left({a}_{\omega }+{b}_{\omega }\right)\exp (-i\omega t)+{\rm{H.c.}}\), the coupling is \({g}_{a}={d}_{ea}{\epsilon }_{{\omega }_{a}}/\hslash\), and the average is calculated at the initial state \(\left|a,{1}_{a}\right\rangle\). The Heisenberg picture also explains why work can be finite even though the average driving field is precisely zero at the single-photon state (i.e., \(\left\langle {1}_{a}\right|{E}_{{\rm{in}}}(t)\left|{1}_{a}\right\rangle =0\)). This shows how the dissipative adaptation can be robust to a phase-incoherent work source, in contrast to the classical forces of well-defined phases used in ref. 4 (in our model, a semiclassical driving field with a well-defined phase would correspond to an initial coherent, or Glauber, state \(\left|\alpha \right\rangle ={\prod }_{\omega }\left|{\alpha }_{\omega }\right\rangle\), where \({a}_{\omega }\left|{\alpha }_{\omega }\right\rangle ={\alpha }_{\omega }\left|{\alpha }_{\omega }\right\rangle\)). In Eq. (26), \({{{\Gamma }}}_{a,b}\propto {g}_{a,b}^{2}\) reveals that the heat is related with the variance of the interaction energy, explaining why the vanishing average interaction energy does not hinder energy exchange in the form of heat and confirming its stochastic nature. The energy conservation in Eq. (5) results, of course, from defining \({\left\langle {Q}_{{\rm{abs}}}\right\rangle }_{a}\equiv -{\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\) in Eq. (24). Finally, \(\left\langle {H}_{S}(0)\right\rangle =0\) and \(\left\langle {H}_{S}(\infty )\right\rangle ={p}_{a\to b}(\infty )\hslash {\delta }_{ab}\). We remind that, during the dynamical transition from \(\left|a\right\rangle\) to \(\left|b\right\rangle\), we have \(p_e(t) = |\psi(t)|^2\), allowing us to establish an exact adaptation-energy relation between Eqs. (18) and (26). We finally find our quantum adaptation-work relation,

$${p}_{a\to b}(\infty )=\frac{{{{\Gamma }}}_{b}}{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}\,\frac{{\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}}{\hslash {\omega }_{a}},$$
(27)

valid for a photon of arbitrary pulse shape and resonant with ωa, as well as for arbitrary δab, as stated in Eqs. (4). In the case of the exponential pulse used in Eq. (19), for instance, we find that \({\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}=\hslash {\omega }_{a}\,4{{{\Gamma }}}_{a}/({{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}+{{\Delta }})\), maximized in the monochromatic regime (Δ ≪ Γa + Γb). If in addition Γa = Γb, we get \({\left\langle {W}_{{\rm{abs}}}\right\rangle }_{a}=2\hslash {\omega }_{a}\), which is twice the initial average energy contained in the single-photon pulse. This counterintuitive result reinforces the notion that work is the amount of energy transferred during a process, rather than the average energy stored in a system at a given time (after all, unitary evolutions with time-independent Hamiltonians conserve energy; here, \({\partial }_{t}\left\langle H(t)\right\rangle =0\) implies that \(\left\langle {H}_{S}(t)\right\rangle +\left\langle {H}_{I}(t)\right\rangle +\left\langle {H}_{E}(t)\right\rangle ={p}_{b}(0)\hslash {\delta }_{ab}+\hslash {\omega }_{L}\), so at resonance and for pb(0) = 0 we find \(\left\langle {H}_{S}(t)\right\rangle \le \hslash {\omega }_{L}\)). Equation (27) is our key signature of a quantum dissipative adaptation. The seeming low-excitation issue (∣ψ(t)∣2 ≪ 1, in the monochromatic limit) has finally been clarified, given that the time integral of ∣ψ(t)∣2 (in Eq. (18)) is not only finite, but also linearly proportional to the work required for the quantum dissipative adaptation to take place.

Entropy in the self-organization

In the classical formulation of dissipative adaptation, entropy production has a key role in establishing the connection between energy transfer and statistical irreversibility, requiring no detailed knowledge on the state of the environment. Here, irreversibility is readily characterized by the asymmetry between the forward, pa→b(t), and the backward, \({p}_{b\to a}^{* }(t)\), processes. Nevertheless, we have the advantage that we keep the full description of the quantum state of the system plus the environment, ρ(t) (at the expense of a greater degree of generality in our global Hamiltonian H). With \({\rho }_{E}(t)={{\rm{tr}}}_{S}[\rho (t)]\) at hands (see Methods), we now seek to describe what happens to the environment during and after the drive interacts with the three-level system. We calculate the exact expression for the von Neumann entropy,

$${S}_{E}(t)=-{\rm{tr}}[{\rho }_{E}(t)\mathrm{ln}\,{\rho }_{E}(t)].$$
(28)

See Eqs. (45) and (46) in the Methods for the analytic expression. We have found that our classical intuition, namely, that better adaptation produces more entropy in the environment, can be found by an appropriate distinction between the classical and the quantum contributions to SE.

The idea behind this distinction is the following. Let us first suppose that the system is initially at \(\left|a\right\rangle\). A highly monochromatic incoming photon will fully induce transition from state \(\left|a\right\rangle\) to \(\left|b\right\rangle\). Hence, the outgoing photon will be detected at modes \(\left\{{b}_{\omega }\right\}\), as well. Now, by considering an initially mixed state of the three-level system as given by \((1/2)\left|a\right\rangle \left\langle a\right|+(1/2)\left|b\right\rangle \left\langle b\right|\), a highly monochromatic incoming photon will have 1/2 probability of leaving at \(\left\{{a}_{\omega }\right\}\) (in the cases where it encounters the system at \(\left|b\right\rangle\)) and 1/2 of leaving at \(\left\{{b}_{\omega }\right\}\) (in the cases where it encounters the system at \(\left|a\right\rangle\)), so the final global state would be mixed and separable,

$$\rho (\infty )\approx \left|b\right\rangle \left\langle b\right|\otimes \left(\frac{1}{2}\left|{\tilde{1}}_{b}\right\rangle \left\langle {\tilde{1}}_{b}\right|+\frac{1}{2}\left|{1}_{a}^{{\rm{free}}}\right\rangle \left\langle {1}_{a}^{{\rm{free}}}\right|\right).$$
(29)

\(\left|{\tilde{1}}_{b}\right\rangle\) and \(\left|{1}_{a}^{{\rm{free}}}\right\rangle\) have been defined in Eqs. (43) and (44) below (see Methods). This is what we call the classical contribution to the final mixed state of the field: the initial mixture of the system is fully transferred to the environment.

Let us now consider again that the system is initially at \(\left|a\right\rangle\) (i.e., \({p}_{a}^{(0)}=1\)). However, let us assume that the linewidth of the pulse is of the order of the dissipation rates of the three-level system. In that case, the final state of the global system becomes entangled,

$$\left|\xi (\infty )\right\rangle =\sqrt{{N}_{a}}\left|a,{\tilde{1}}_{a}\right\rangle +\sqrt{{N}_{b}}\left|b,{\tilde{1}}_{b}\right\rangle ,$$
(30)

Nk being the average number of photons at modes k (see Eq. (42) in the Methods). Therefore, the quantum state of the environment in this case, \({\rho }_{E}(\infty )={{\rm{tr}}}_{S}[\left|\xi (\infty )\right\rangle \left\langle \xi (\infty )\right|]\), is also mixed between modes \(\left\{{a}_{\omega }\right\}\) and \(\left\{{b}_{\omega }\right\}\). However, the mixture in this case arises from a sustained system–environment quantum entanglement rather than from the statistical mixture in the initial state of the system.

To unravel these two entropy contributions, we define the classical contribution to the environment entropy as

$${S}_{E}^{c}(t)\equiv {S}_{E}(t)-{p}_{a}^{(0)}S({{\rm{tr}}}_{S}[\left|\xi (t)\right\rangle \left\langle \xi (t)\right|]),$$
(31)

where \(S(\bullet )\equiv -{\rm{tr}}[\bullet \mathrm{ln}\,\bullet ]\) is the von Neumann entropy (see Eq. (51) in the Methods). In Eq. (31), we have taken into account that only the term proportional to \({p}_{a}^{(0)}\) in ρS(0) generates the entanglement discussed in Eq. (30). Now we focus on the long time limit, t→∞. Most interestingly, we have analytically expressed \({S}_{E}^{c}(\infty )\) as a function of pa→b(∞) (see Eqs. (45)–(49) in the Methods). Figure 2 illustrates that \({S}_{E}^{c}(\infty )\) vs. pa→b(∞) (solid black line) is a monotonic function, in contrast with the non-monotonic SE(∞) vs. pa→b(∞) (dashed blue line). This monotonic behavior shows that the more organized is the three-level system, the higher is the classical contribution to the entropy of the environment at long times t→∞. Besides having derived \({S}_{E}^{c}(\infty )\) as a function of pa→b(∞), we have also expressed pa→b(∞) as a function of the average dissipated heat in the \(\left|a\right\rangle\) to \(\left|b\right\rangle\) transition, \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\), with the help of Eq. (27). This establishes the function SE(∞) vs. \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\), as mentioned earlier. In addition, we have derived \({S}_{E}^{c}(\infty )\) as a function of \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\). This is even more meaningful, because we have found that the function \({S}_{E}^{c}(\infty )\) vs. \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\) is monotonic. To provide an example illustrating this monotonicity, let us take the degenerate case (δab = 0) with equal dissipation rates (Γa = Γb). In that case, we have that \({p}_{a\to b}(\infty )={\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}/(2\hslash {\omega }_{a})\), showing that \({S}_{E}^{c}(\infty )\) vs. \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\) is monotonic, since \({S}_{E}^{c}(\infty )\) vs. pa→b(∞) is monotonic as well (Eq. (50) in the Methods provides the more general relation). The monotonicity here strengthens the signature of dissipative adaptation. Namely, maximal adaptation (irreversible self-organization) not only costs maximal work absorption (as we have shown in the energetics analysis), but also maximizes the dissipated heat which, in turn, maximizes the classical contribution to the environment entropy.

Fig. 2: Entropy of the environment.
figure 2

The classical contribution to the environment entropy, \({S}_{E}^{c}(\infty )\), increases monotonically as a function of the transition probability from \(\left|a\right\rangle\) to \(\left|b\right\rangle\), pa→b(∞) (solid black), providing an additional signature of the dissipative adaptation. The environment entropy, SE(∞), presents a non-monotonic behavior as a function of pa→b(∞) (dashed blue), due to a sustained system–environment quantum entanglement that can contribute to the entropy production of the environment. In both curves, we refer to the asymptotic state t→∞ of the dynamics. We use Γa = Γb as the environment-induced spontaneous emission rates and \({p}_{a}^{(0)}={p}_{b}^{(0)}=1/2\) as the probabilities in the initial state of the lambda system, \({\rho }_{S}(0)={p}_{a}^{(0)}\left|a\right\rangle \left\langle a\right|+{p}_{b}^{(0)}\left|b\right\rangle \left\langle b\right|\) (the curves do not qualitatively depend on these choices). The dotted gray line is \(\mathrm{ln}\,(2)\). Note that both \({S}_{E}^{c}(\infty )\) as a function of pa→b(∞) and SE(∞) as a function of pa→b(∞) do not explicitly depend on the photon pulse shape, though pa→b(t) does; e.g., attaining pa→b(∞)→1 requires an extremely monochromatic (long) photon, whereas the pa→b(∞)→0 limit is attained in the extremely broadband (short) pulse regime.

Full size image

Discussion

Our results establish the quantum dissipative adaptation underlying the driven self-organization of a quantum state, going beyond the classical formulation. We have explored an elementary fully quantized model, exactly solvable for both the system and the environment, where the irreversibly self-organized quantum state of a three-level system is owing to the work absorption from a single-photon pulse, part of which is dissipated to the environment (as shown in Eqs. (4) and (5)). The irreversibility of this self-organization became clear from the asymmetry in the transition probabilities in the forward and backward processes. Finally, with the purpose of providing an additional signature of dissipative adaptation, we have analytically investigated the environment entropy. We have found that the classical contribution to the environment entropy variation is a monotonic function of pa→b(∞) (as illustrated in Fig. 2), showing that maximizing adaptation not only requires maximal work absorption, but also leads to maximal increase in the classical contribution to the environment entropy, due to maximal dissipated heat.

We remind that the meaning of this increase in the environment entropy is that of a statistical mixing between the field modes \(\left\{{a}_{\omega }\right\}\) and \(\left\{{b}_{\omega }\right\}\), as we have analyzed in Eq. (29). Our model’s dynamics do not lead to the thermalization of the environment, which always remains out of equilibrium. If we wish to take the analogy further, so as to give a thermodynamic meaning to the entropy increase, we will have to add extra ingredients. In practice, there could be for instance a slight increase in the local temperature of the environment in the vicinity of the Λ system. To describe that kind of effect, the model should consider a finite heat capacity environment (in contrast to our infinite-size environment) and some auxiliary light-matter coupling mechanism that could effectively create an interaction between the many frequency modes of the light field. Such features provide means for the environment state to eventually approximate a Gibbs state at some new finite temperature \(T^{\prime} \, > \,T=0\) (producing even more entropy than that we have calculated from our model), in the spirit of Timpanaro et al.25.

Before concluding, we would like to point out how our model can be significant to quantum many-body systems. Our intention is to indicate how broadly applicable the concept of quantum dissipative adaptation may become. First, we notice that the Λ structure of energy levels can arise from the quantization of collective degrees of freedom describing many interacting atoms and electrons, as happens in the so-called artificial atoms (e.g., electron-hole pairs in semiconductor quantum dots26 and the quantized magnetic flux in superconducting rings27). Lodahl et al.26 and Gu et al.27 also discuss how these artificial atoms can be driven by single-photon pulses propagating in one-dimensional waveguides, building the closest possible scenario to that in our model.

We can also envision the quantum dissipative adaptation in ensembles of non-interacting atoms or spins. That is more easily seen when we realize that our idea of a driven self-organized quantum state is notably similar to the dynamics induced by stimulated Raman adiabatic passage (STIRAP)28, a versatile and robust technique that has been performed, e.g., in ultracold gases, in doped crystals and in nitrogen-vacancy centers. STIRAP consists of an efficient population transfer between two discrete quantum states of an ensemble of emitters (usually the lowest levels of Λ systems) by coherently coupling them with two radiation fields (well-controlled classical pulses) through an unpopulated intermediate state. The connection we have in mind between the driven self-organization provided by STIRAP and dissipative adaptation becomes more evident in light of the recent proposal for using STIRAP as a tool for spectral hole burning (SHB) in inhomogeneously broadened systems29. The reason is that the mechanism behind standard SHB30 turns out to be precisely that of classical dissipative adaptation, as discussed by Kedia et al.31. Namely, those dipoles that get excited by the resonant drive (whose frequency can be swept on demand) can become irreversibly trapped in dark states (at sufficiently low temperatures). The difference concerning the newly proposed STIRAP-based SHB29, or the single-photon pulse that we have studied here, is the quantum coherent nature of the process. Understanding how STIRAP-based SHB depends on temperature seems a valuable opportunity for widening the concept of quantum dissipative adaptation (similarly, Ropp et al.10 show a groundbreaking experiment on temperature-dependent dissipative self-organization in optical space).

Optomechanical nanoresonators32 also hold the promise of displaying some kind of quantum dissipative adaptation in the lines of our results. We glean this notion from Kedia et al.31, where signatures of dissipative adaptation are shown in disordered networks of classical bistable springs. In the mechanical nanoresonators by Yeo et al.32, bistability can arise from a strain-mediated coupling between the center of mass of an oscillating nanowire and the quantum state of a single semiconductor quantum dot embedded therein. This kind of coupling between an optically controllable microscopic degree of freedom (within the quantum dot) and a mesoscopic degree of freedom (the nanoresonator center of mass) opens appealing perspectives for studying dissipative adaptation at a quantum-classical boundary of a many-body system. As a last word on resonators, it seems relevant investigating whether the vibration-assisted exciton transport found in photosynthetic light-harvesting antennae17 could be related with a quantum dissipative adaptation.

To conclude, the quantum dissipative adaptation we have found can be regarded as a proof of principle, in need of generalizations towards many directions. As a first example, it would be worth investigating dissipative adaptation in larger Hilbert spaces (i.e., in a vast energy landscape, in the words of ref. 1,4). Multistability makes self-organization in classical many-body systems a fascinating problem1,4,31. Curiously, the multistability from complex classical systems is reminiscent of our zero-temperature model. Among the huge amount of available “non-organized” stable states that our Λ system can occupy, only a single (exceptional) state is populated when the work source is optimal. In other words, although we find an infinite number of (pure or mixed) combinations of states \(\left|a\right\rangle\) and \(\left|b\right\rangle\) that are equally stationary in the non-driven zero-temperature environment, the system always ends up in the rare state \(\left|a\right\rangle\) (irrespective of the infinitely many possible initial states) once it has been driven by the suitable pulse (which maximizes the absorbed work in the \(\left|a\right\rangle\) to \(\left|b\right\rangle\) transition).

It also remains to be investigated how other sources of work and finite temperatures affect dissipative adaptation in quantum systems. In our Λ system, for instance, we expect that the forward (\(\left|a\right\rangle\) to \(\left|b\right\rangle\)) transition will keep being favored with respect to the backward (\(\left|b\right\rangle\) to \(\left|a\right\rangle\)) transition at finite temperatures, as long as the work source is strong enough. Also, the time dependence of this asymmetry should be transient or stationary, depending on the work source type (pulsed or continuous). In quantum thermodynamics, the influence of the environment temperature on entropy production at the quantum regime, on the one hand, has been the subject of many recent studies25,33,34,35,36. The focus on model-independent aspects is justified not only due to the interest in generality, but also due to the closest possible analogies with classical fluctuation theorems and with the Landauer erasure principle (that establishes how information erasure is connected to thermodynamic entropy production due to heat dissipation). On the other hand, Manzano et al.37, for instance, provide an alternative formalism for quantum fluctuation theorems that go beyond thermal-equilibrium states of the environment. It may be the case that the quantum fluctuation theorems developed in the papers above (and in the references therein) present fruitful methods in the search for a general quantum theory of dissipative adaptation. In summary, we provide the starting point towards a quantum thermodynamics of driven self-organization.

Methods

Global system-plus-reservoir quantum dynamics

To obtain explicit expressions for the transition probabilities in Eq. (17), we need to solve the Schrödinger equation \(i\hslash {\partial }_{t}\left|\xi (t)\right\rangle =H\left|\xi (t)\right\rangle\), where \(\left|\xi (t)\right\rangle =\psi (t)\left|e,0\right\rangle +{\sum }_{\omega }{\phi }_{\omega }^{a}(t){a}_{\omega }^{\dagger }\left|a,0\right\rangle +{\phi }_{\omega }^{b}(t){b}_{\omega }^{\dagger }\left|b,0\right\rangle\). After a Wigner–Weisskopf approximation and using that \({\phi }_{n}(z,t)\equiv {\sum }_{\omega }{\phi }_{\omega }^{n}(t)\exp (i{k}_{\omega }z)\), with a linear dispersion relation, ω = ckω, and ∑ω → ∫dωϱω ≈ ϱ∫dω, we find that

$${\partial }_{t}\psi (t)=-\left(\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}{2}+i{\omega }_{a}\right)\psi (t)-{g}_{a}{\phi }_{a}(-ct,0),$$
(32)

with

$${\phi }_{a}(z,t)\,= \, {\phi }_{a}(z-ct,0)\\ \,+\sqrt{2\pi \varrho {{{\Gamma }}}_{a}}{{\Theta }}(z){{\Theta }}(t-z/c)\psi (t-z/c)$$
(33)

and

$${\phi }_{b}(z,t)\,= \, {\phi }_{b}(z-ct,0){e}^{-i{\delta }_{ab}t}\\ +\,\sqrt{2\pi \varrho {{{\Gamma }}}_{b}}{{\Theta }}(z){{\Theta }}(t-z/c)\psi (t-z/c){e}^{-i{\delta }_{ab}z/c},$$
(34)

where \({{{\Gamma }}}_{k}=2\pi {g}_{k}^{2}\varrho\).

Integrating Eq. (32) for ψ(0) = 0 gives

$$\psi (t)=-{g}_{a}\mathop{\int}\nolimits_{\!\!\!\!0}^{t}{\phi }_{a}(-ct^{\prime} ,0){e}^{-\left(\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}{2}+i{\omega }_{a}\right)(t-t^{\prime} )}dt^{\prime} ,$$
(35)

which depends on the initial photon wavepacket,

$${\phi }_{a}(z,0)={\phi }_{a}^{{\rm{shape}}}(z,0){e}^{i{\omega }_{L}\frac{z}{c}}.$$
(36)

ωL is the central frequency of the pulse and \({\phi }_{a}^{{\rm{shape}}}(z,0)\) is its spatial shape. For the sake of providing an explicit and physically motivated example, we consider in Eq. (19) the shape of the initial photon wavepacket to be an exponential (raising in space, equivalent to decaying in time for a right-propagating pulse), as typical from spontaneous emission,

$${\phi }_{a}^{{\rm{shape}}}(z,0)=N{{\Theta }}(-z){e}^{\frac{{{\Delta }}}{2}\frac{z}{c}},$$
(37)

where Δ is the pulse linewidth and \(N=\sqrt{2\pi \varrho {{\Delta }}}\) is a normalization factor (considering ϕb(z, 0) = 0). Finally, and substituting \({{{\Gamma }}}_{a}=2\pi {g}_{a}^{2}\varrho\), we have that

$$\psi (t)=-{f}_{{{\Delta }}}\,{e}^{-\left(\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}{2}+i{\omega }_{a}\right)t}\left[{e}^{\left(\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}-{{\Delta }}}{2}-i{\delta }_{L}\right)t}-1\right],$$
(38)

where

$${f}_{{{\Delta }}}\equiv \frac{\sqrt{{{{\Gamma }}}_{a}{{\Delta }}}}{\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}-{{\Delta }}}{2}-i{\delta }_{L}}$$
(39)

and δL ≡ ωL − ωa.

The average interaction energy is given by

$$\left\langle {H}_{I}(t)\right\rangle \equiv {\rm{tr}}[\rho (t){H}_{I}].$$
(40)

From Eqs. (6), (12), and the definition \(\left|\xi (t)\right\rangle \equiv U\left|a,{1}_{a}\right\rangle\), we have that

$$\left\langle {H}_{I}(t)\right\rangle ={p}_{a}^{(0)}\left\langle \xi (t)\right|{H}_{I}\left|\xi (t)\right\rangle +{p}_{b}^{(0)}\left\langle b,{1}_{a}\right|{H}_{I}\left|b,{1}_{a}\right\rangle .$$
(41)

Using HI from Eq. (9), we have that \(\left\langle b,{1}_{a}\right|{H}_{I}\left|b,{1}_{a}\right\rangle =0\) and \(\left\langle \xi (t)\right|{H}_{I}\left|\xi (t)\right\rangle =2(\hslash {g}_{a}{\rm{Im}}[{\psi }^{* }(t){\phi }_{a}(-ct,0)]+\hslash {g}_{b}{\rm{Im}}[{\psi }^{* }(t){\phi }_{b}(-ct,0)])\). Choosing ϕb(z, 0) = 0 obviously makes \({\rm{Im}}[{\psi }^{* }(t){\phi }_{b}(-ct,0)]=0\). By choosing the resonance condition, δL = 0, we have from Eqs. (35) and (36) that \({\rm{Im}}[{\psi }^{* }(t){\phi }_{a}(-ct,0)]=0\). This shows that \(\left\langle {H}_{I}(t)\right\rangle =0\) at resonance.

Entropy of the environment

The quantum state of the environment in our model is obtained from the global initial state \(\rho (0)=({p}_{a}^{(0)}\left|a\right\rangle \left\langle a\right|+{p}_{b}^{(0)}\left|b\right\rangle \left\langle b\right|)\otimes \left|{1}_{a}\right\rangle \left\langle {1}_{a}\right|\). We obtain that

$${\rho }_{E}(t)\,= \, {p}_{a}^{(0)}| \psi (t){| }^{2}\left|0\right\rangle \left\langle 0\right|+{p}_{a}^{(0)}{N}_{b}(t)\left|{\tilde{1}}_{b}\right\rangle \left\langle {\tilde{1}}_{b}\right|\\ \,+{p}_{a}^{(0)}{N}_{a}(t)\left|{\tilde{1}}_{a}\right\rangle \left\langle {\tilde{1}}_{a}\right|+{p}_{b}^{(0)}\left|{1}_{a}^{{\rm{free}}}\right\rangle \left\langle {1}_{a}^{{\rm{free}}}\right|.$$
(42)

We have defined \({N}_{k}(t)\equiv {\sum }_{\omega }| {\phi }_{\omega }^{k}(t){| }^{2}\), with

$$\left|{\tilde{1}}_{k}\right\rangle \equiv {N}_{k}^{-1/2}\sum _{\omega }{\phi }_{\omega }^{k}(t){k}_{\omega }^{\dagger }\left|0\right\rangle ,$$
(43)

and

$$\left|{1}_{a}^{{\rm{free}}}\right\rangle \equiv \sum _{\omega }{\phi }_{\omega }^{a}(0)\exp (-i\omega t)\,{a}_{\omega }^{\dagger }\left|0\right\rangle .$$
(44)

To explicitly calculate the von Neumann entropy,

$${S}_{E}= \, -{\rm{tr}}[{\rho }_{E}(t)\mathrm{ln}\,{\rho }_{E}(t)],\\ = \,-\sum _{j}{\lambda }_{j}{\mathrm{ln}}\,{\lambda }_{j},$$
(45)

we need to find the eigenvalues λj of ρE. The exact diagonalization of ρE gives us four non-zero eigenvalues, namely, \({\lambda }_{1}={p}_{a}^{(0)}| \psi (t){| }^{2}\), \({\lambda }_{2}={p}_{a}^{(0)}{N}_{b}(t)\), and

$${\lambda }_{3,4}= \, \frac{{p}_{a}^{(0)}{N}_{a}+{p}_{b}^{(0)}}{2}\\ \,\pm \frac{1}{2}\sqrt{{\left({p}_{a}^{(0)}{N}_{a}-{p}_{b}^{(0)}\right)}^{2}+4{p}_{a}^{(0)}{p}_{b}^{(0)}{N}_{a}\left|\left\langle {1}_{a}^{{\rm{free}}}\right|{\tilde{1}}_{a}\right\rangle {| }^{2}}.$$
(46)

Terms Na and Nb can be interpreted as the average number of photons in each continuum of modes. Their mathematical expressions, however, are invariably connected with the three-level system probabilities, namely,

$${N}_{k}(t)={p}_{a\to k}(t).$$
(47)

We calculate the overlap between the states representing the free propagation and the reemitted photon at mode a in the real-space representation,

$$\left\langle {1}_{a}^{{\rm{free}}}\right|\left.{\tilde{1}}_{a}\right\rangle =\frac{1}{\sqrt{{N}_{a}}}\frac{1}{2\pi \varrho c}\mathop{\int}\nolimits_{-\infty }^{\infty }dz\,{\phi }_{{\rm{a,free}}}^{* }(z,t){\phi }_{a}(z,t),$$
(48)

where ϕa, free(z,t) ≡ ϕa(z − ct, 0).

With the aid of (32) and (33), and using integration by parts, we find that

$$\left\langle {1}_{a}^{{\rm{free}}}\right|\left.{\tilde{1}}_{a}\right\rangle \sqrt{{N}_{a}}=1-\left(\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}{2{{{\Gamma }}}_{b}}\right){p}_{a\to b}(\infty ),$$
(49)

valid at long times t→∞. Equation (49) provides the core connection between the overlap and the transition probability at long times that we needed. Additionally, Na(∞) = 1 − Nb(∞) = 1 − pa→b(∞) and λ1(∞) = 0. Since Eqs. (27) and (24) link pa→b(∞) to \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\),

$${p}_{a\to b}(\infty )={\left(\hslash {\omega }_{a}\frac{{{{\Gamma }}}_{a}+{{{\Gamma }}}_{b}}{{{{\Gamma }}}_{b}}-\hslash {\delta }_{ab}\right)}^{-1}{\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a},$$
(50)

we have thus analytically established the function SE(∞) vs. \({\left\langle {Q}_{{\rm{diss}}}\right\rangle }_{a}\), valid at long times t → ∞, under the resonant condition ωL = ωa, for any arbitrary pulse shape.

Finally, the classical contribution to the entropy is \({S}_{E}^{c}\equiv {S}_{E}-{p}_{a}^{(0)}S({{\rm{tr}}}_{S}[\left|\xi (t)\right\rangle \left\langle \xi (t)\right|])\), where

$$S({{\rm{tr}}}_{S}[\left|\xi (t)\right\rangle \left\langle \xi (t)\right|])=-\left({N}_{a}{\mathrm{ln}}\,{N}_{a}+{N}_{b}{\mathrm{ln}}\,{N}_{b}+| \psi (t){| }^{2}{\mathrm{ln}}\,| \psi (t){| }^{2}\right),$$
(51)

for which we also have the analytic solution.