Thermal machines beyond the weak coupling regime

The theory of thermodynamics originates from the study of thermal machines in the early industrial age, when it was of utmost importance to find out what rates of work extraction could ultimately be achieved. Early on, it became clear that the theory of thermal machines would be intimately related to topics of fundamental physics such as statistical mechanics and notions of classical information theory [1]. Here, the interplay and relations between the widely studied notions of work, entropy and statistical ensembles are the focus of attention. Concomitant with technological developments, the theory also became more intricate and addressed more elaborate situations. Famous thought experiments such as Maxwellʼs demon, Landauerʼs erasure and Slizardʼs engine have not only puzzled researchers for a long time, but today also serve as a source of inspiration for quantitative studies of achievable rates when employing thermal machines [2–5]. Indeed, with nano-machines operating at or close to the quantum level coming into reach, there has recently been an explosion of interest in the question of what role quantum effects may possibly play. The potential and limits of work extraction with physically plausible operations which respect quantum correlations are established in this work.

The role of correlations is already a challenging problem with a long history in classical thermodynamics. Thermal machines comprise a system that is brought into contact with a thermal bath. This process introduces correlations that are typically disregarded by assuming that the interaction between the system and bath is sufficiently weak. Due to the limited applicability of this assumption in practical situations, there have been great efforts toward characterizing thermal machines beyond the weak-coupling regime in specific models both in the classical and quantum settings [6–9]. However, a general framework for work-extraction beyond the weak-coupling assumption is still missing. This is mainly due to a lack of understanding of the process of evolution towards equilibrium under the effect of generic strong couplings, which has only started to be tackled to its full extent in recent years.

More specifically, let us introduce the weak-coupling assumption precisely as: a system S with Hamiltonian H_S, when put into weak thermal contact with a thermal bath B, equilibrates towards the state

$$\begin{eqnarray}&&{{\rho }_{S}}=\omega ({{H}_{S}})\end{eqnarray} \tag{ 1 }$$

with $\omega ({{H}_{S}})\;:=\;-{{{\rm e}}^{\beta {{H}_{S}}}}/Z$, $Z\;:=\;{\rm tr}(-\beta {{H}_{S}})$ and $\beta \gt 0$ being the inverse temperature. That is, S equilibrates to the usual Gibbs ensemble. Note that this notion of weak coupling can in general differ from the one sometimes used in the study of open quantum systems leading to Markovian dynamics of the sub-system S [3, 10]. The precise conditions on the coupling so that (1) is fulfilled have been recently tackled in the quantum setting: the strength of the coupling Hamiltonian V—measured in an adequate norm—has to be negligible in comparison with the intensive thermal energy scale ${{\beta }^{-1}}$ [11]. This formalizes the usual derivation of the canonical ensemble from the micro-canonical one in classical statistical mechanics where the coupling energy is neglected. Note that the interaction strength typically scales as the boundary of the sub-system S. Hence, in spatial dimensions higher than one, the weak-coupling assumption cannot hold true if one increases the system size. This will be the case regardless of the strength of the coupling per particle or the relative size between S or B [11, 12]. Therefore, the weak-coupling assumption is arguably inapplicable not only to realistic situations, but also to idealized systems whose constituents interact weakly.

Recently, relaxation towards equilibrium in the strong-coupling case has been addressed from the perspective of canonical typicality. The idea is that closed non-integrable many-body systems, described, however, by a unitarily evolving pure state, are generically expected to equilibrate [12–17]. Such systems behave—the overwhelming majority of times—as if they were described by a thermal state when considering the expectation values of local observables [11, 15, 18]. The eigenstate thermalization hypothesis [18–20] gives further substance to this expectation. This means that when a sub-system S is put in contact with a bath B with Hamiltonian H_B, the equilibrium state is not (1), but the reduction of the global Gibbs state of S and B

$$\begin{eqnarray}&&{{\rho }_{S}}={\rm t}{{{\rm r}}_{B}}(\omega ({{H}_{SB}})),\end{eqnarray} \tag{ 2 }$$

where ${{H}_{SB}}={{H}_{S}}+{{H}_{B}}+V$ [15].

In this work, we incorporate these recent insights to describe equilibration in the analysis of thermodynamics beyond the weak-coupling regime. We provide rigorous bounds on the optimal work extraction in the presence of thermal baths whose effect is to drive systems to an equilibrium state of the form (2). Our approach considers protocols of work extraction by performing quantum quenches on sub-systems in strong coupling with thermal baths.

We show that the strong coupling between the system and bath may induce an unavoidably irreversible component in the process and we discuss to what extent this results in a limitation on the optimal work extraction. We are able to quantify this deficit in terms of standard thermodynamic functions as the free energy and we show that it prevents one from saturating the second law of thermodynamics. Our results are completely general in the sense that they do not make use of any specific model for the description of the system or bath.

The work extraction problem requires at least the following elements.

• A system S. This is the part of the machine upon which one has control, i.e., it is possible to engineer its Hamiltonian H_S. By no means is the conservation of energy violated in this prescription.
• A battery W. This models energy storage and accounts for the energy supplied and extracted from the system S. It can be seen as a lifted weight. Any Hamiltonian with a suitably dense spectrum will be suitable.
• A thermal bath B. When the system S is put into contact with the thermal bath, S is assumed to thermalize in the sense of (2), with ${{H}_{SB}}={{H}_{S}}+{{H}_{B}}+V$, where V is the interaction that couples the system and bath. The interaction V is assumed to be fixed and not tuneable by the operator of the machine. No assumptions are made on the state of SB.

A scheme of the setting is shown in figure 1. The problem of work extraction consists of maximising the energy of the final state of the battery given an initial state S, an initial Hamiltonian ${{H}^{(0)}}$ and a set of operations. In our case, the set of operations are Hamiltonian transformations and thermalizations.

Zoom In Zoom Out Reset image size

What we refer to as a Hamiltonian transformation is a change of the Hamiltonian of the system and/or the switch on/off of the interaction V between the system and bath. Hence, at the end of each transformation, the Hamiltonian of SBW takes the form

$$\begin{eqnarray}&&{{H}^{(i)}}=H_{S}^{(i)}+{{V}^{(i)}}+{{H}_{B}}+{{H}_{W}},\end{eqnarray} \tag{ 3 }$$

and ${{H}^{(i)}}$ is taken to ${{H}^{(i+1)}}$, while ${{V}^{(i)}}$ takes values from $\{0,V\}$, for $i=0,...,n-1$. In order for this Hamiltonian transformation to be meaningful and to allow for a fair account of the work extracted, however, we require the following two natural conditions to be fulfilled:

• (i)  
  Quenches. The reduced state on SB does not change,

  $$\begin{eqnarray}&&\rho _{SB}^{(i)}=\rho _{SB}^{(i+1)},\end{eqnarray} \tag{ 4 }$$

  modelling the behaviour of the system when the Hamiltonian acting on that sub-system is changed abruptly.
• (ii)  
  Energy conservation. The mean total energy is preserved, i.e., for each transformation the energy change of the system (due to the change of its Hamiltonian) has been supplied or stored by the battery

  $$\begin{eqnarray}&&{\rm tr}({{\rho }^{(i)}}{{H}^{(i)}})={\rm tr}({{\rho }^{(i+1)}}{{H}^{(i+1)}}).\end{eqnarray} \tag{ 5 }$$

The first condition merely states that on the time scale of the dynamics taking place, the sudden approximation holds true in system SB, or in other words, that one performs a quench. The latter is not an assumption, but rather a necessary condition for a fair account of all the energy supplied or extracted from the thermal machine SB. Note that we do not impose that the machine SB is energetically isolated, which is obviously not the case since we consider time-dependent Hamiltonians that obviously do not preserve energy. Condition (5) merely states that the energy gained/lost by SB is supplied from or stored in the battery, which plays the role of the usual lifted weight in thermodynamics. The average work ${{\langle W\rangle }_{i\to i+1}}$ extracted in the quench $i\to i+1$ is the average energy change in the battery when the quench is performed. From (4) and (5) one obtains

$$\begin{eqnarray}&&{{\langle W\rangle }_{i\to i+1}}={\rm tr}(\rho _{SB}^{(i)}(H_{SB}^{(i)}-H_{SB}^{(i+1)})).\end{eqnarray} \tag{ 6 }$$

This is the standard way of accounting work as the energy difference of the combined system and bath due to the time-dependent Hamiltonian giving rise to the evolution; in this specific case, a quench [21]. Note also that no assumption is made on the global state of SBW and on the possible correlations between the battery W and the system S after implementing a quench. For the case of a unitary implementation, this issue is discussed in detail in the appendix.

A thermalization map is a map that models the effect of putting the system into actual contact with the heat bath B by thermalizing it as described in (2). This transformation can be applied only when the system and bath are interacting, that is ${{V}^{(i)}}=V$. This family of maps is physically motivated by the realistic behaviour of evolution under generic Hamiltonians¹ . However, within the abstract level of the set of operations it can be regarded simply as any completely positive map $\mathcal{T}$ acting on quantum states of SB with the defining property

$$\begin{eqnarray}&&{\rm t}{{{\rm r}}_{B}}(\mathcal{T}(\rho ))={\rm t}{{{\rm r}}_{B}}(\omega ({{H}_{SB}})).\end{eqnarray} \tag{ 7 }$$

Applying the map $\mathcal{T}$ has no effect on the battery and hence it does not have any work cost. The evolution towards equilibrium is reached by the dynamics of SB alone once they are interacting, without having to supply or extract energy from the battery or implement any change in the Hamiltonians of SB.

A sequence of such operations is called a protocol, which we denote by $\mathcal{P}$, and is specified by: (i) a list of Hamiltonians $\{{{H}^{(i)}}\}_{i=1}^{n-1}$ of the form (3) and (ii) a set of instructions specifying when the thermalization maps are realized. In order to avoid the energy in the battery originating from a change of the system Hamiltonian, we consider protocols with the final Hamiltonian being equal to the initial one, $H_{S}^{(n)}=H_{S}^{(0)}$. The work extracted in expectation by a protocol $\mathcal{P}$ for an initial state $\rho _{S}^{(0)}$ and an initial and final Hamiltonian $H_{S}^{(0)}$ is defined as the energy increase of the battery

$$\begin{eqnarray}&&\langle W\rangle (\mathcal{P},{{H}^{(0)}},{{\rho }^{(0)}})\;:=\;{\rm tr}((\rho _{W}^{(n)}-\rho _{W}^{(0)}){{H}_{W}}).\end{eqnarray} \tag{ 8 }$$

Altogether, the set of operations we consider is a generalization of the one considered in [23–25] in that we allow the changing of the eigenbasis of the Hamiltonian H_S. Importantly, the thermalization process model is not restricted to the weak-coupling regime, but also actually includes quantum correlations, which alters the situation considerably. Nonetheless, more general transformations than the ones restricted by condition (4) could be considered [26], in particular energy preserving unitaries that change the state of S, in the spirit of [27–30]. However, Hamiltonian quenches fairly capture operational capabilities in realistic situations, rather than arbitrary unitaries, and are also sufficient to cover the standard weak-coupling limit [23, 25, 28]. In the supplementary material we discuss possible ways of generalizing our approach using expectation values in even further general settings and issues related to the role of coherence in the battery [31].

Given the previous set of operations, the following theorem introduces a bound on the amount of work that can be extracted. (Bounds on work extraction).

Theorem 1 Given an initial state ${{\rho }^{(0)}}=\rho _{SB}^{(0)}\otimes \rho _{W}^{(0)}$ and an equal initial and final Hamiltonian ${{H}^{(0)}}$, the work that can be extracted by means of any protocol within the set of allowed operations is bounded by

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)}) & \leqslant & F({{\tilde{\rho }}_{SB}},H_{SB}^{(0)})-F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)}) \\ {} & {} & -\mathop{{\rm min} }\limits_{{{\tilde{H}}_{S}}}\left[ F({{\tilde{\rho }}_{SB}},{{\tilde{H}}_{SB}})-F(\omega ({{\tilde{H}}_{SB}}),{{\tilde{H}}_{SB}}) \right] \\ \end{array}\end{eqnarray} \tag{ 9 }$$

where ${{\tilde{H}}_{SB}}\;:=\;{{\tilde{H}}_{S}}+V+{{H}_{B}}$, ${{\tilde{\rho }}_{SB}}$ is any state such that ${\rm t}{{{\rm r}}_{B}}(\tilde{\rho }_{SB}^{(0)})={\rm t}{{{\rm r}}_{B}}(\rho _{SB}^{(0)})=\rho _{S}^{(0)}$ and $F(\rho ,H)\;:=\;{\rm tr}(\rho H)-{{\beta }^{-1}}S(\rho )$ is the free energy of the state ρ, with respect to the Hamiltonian H and inverse temperature β. Furthermore, for any initial state $\rho _{S}^{(0)}$, there exists a protocol ${{\mathcal{P}}^{*}}$ which saturates the bound. This optimal protocol ${{\mathcal{P}}^{*}}$ consists of a quench to the Hamiltonian ${{\tilde{H}}_{S}}$ that minimizes the difference $F({{\tilde{\rho }}_{SB}},{{\tilde{H}}_{SB}})-F(\omega ({{\tilde{H}}_{SB}}),{{\tilde{H}}_{SB}})$ followed by sequence of thermalizations and small quenches that emulates an isothermal reversible process to come back to the initial Hamiltonian.

Proof. Note that any protocol can be expressed as a concatenation of Hamiltonian transformations and thermalizations, and that the energy of the battery only changes in the Hamiltonian transformations. In the first quench, the energy stored or supplied by the battery reads

$$\begin{eqnarray}&&{{\langle W\rangle }_{0\to 1}}={\rm tr}(\rho _{S}^{(0)}(H_{S}^{(0)}-H_{S}^{(1)}))={\rm tr}(\tilde{\rho }_{SB}^{(0)}(H_{SB}^{(0)}-H_{SB}^{(1)})).\end{eqnarray} \tag{ 10 }$$

where $\tilde{\rho }_{SB}^{(0)}$ is any state such that ${\rm tr}(\tilde{\rho }_{SB}^{(0)})=\rho _{S}^{(0)}$. The rest of quenches are performed after a thermalization, hence the work extracted from them can be written as

$$\begin{eqnarray}\begin{array}{rcl} {{\langle W\rangle }_{i\to i+1}} & = & {\rm tr}(\rho _{SB}^{(i)}(H_{SB}^{(i)}-H_{SB}^{(i+1)}))={\rm tr}(\rho _{SB}^{(i)}(H_{S}^{(i)}-H_{S}^{(i+1)})\otimes {{\mathbb{I}}_{B}}) \\ {} & = & {\rm tr}(\omega (H_{SB}^{(i)})(H_{S}^{(i)}-H_{S}^{(i+1)})\otimes {{\mathbb{I}}_{B}}). \\ \end{array}\end{eqnarray} \tag{ 11 }$$

The total work extracted by a protocol is the sum of the work extracted in every Hamiltonian transformation, that is,

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)}) & = & {\rm tr}\left( \tilde{\rho }_{SB}^{(0)}(H_{SB}^{(0)}-H_{SB}^{(1)}) \right) \\ {} & {} & +\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}\left( \omega (H_{SB}^{(i)})(H_{SB}^{(i)}-H_{SB}^{(i+1)}) \right), \\ \end{array}\end{eqnarray} \tag{ 12 }$$

where n is the number of steps of the protocol. By using the identity $F(\omega (H),H)={\rm tr}(H\rho )+{{\beta }^{-1}}{\rm tr}({\rm ln} (\omega (H)\rho ))$ (note that this is valid for any ρ), the extracted work is rewritten as

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)}) & = & {\rm tr}(\tilde{\rho }_{SB}^{(0)}(H_{SB}^{(0)}-H_{SB}^{(1)})) \\ {} & {} & +\mathop{\sum }\limits_{j=1}^{n-1}\left( F\left( \omega (H_{SB}^{(i)}),H_{SB}^{(i)} \right)-F\left( \omega (H_{SB}^{(j+1)}),H_{SB}^{(j+1)} \right) \right), \\ {} & {} & -{{\beta }^{-1}}\mathop{\sum }\limits_{j=1}^{n-1}{\rm tr}\left( {\rm ln} (\omega (H_{SB}^{(i)}))\omega (H_{SB}^{(i)})-{\rm ln} (\omega (H_{SB}^{(i+1)}))\omega (H_{SB}^{(i)}) \right). \\ \end{array}\end{eqnarray} \tag{ 13 }$$

After identifying the relative entropy $S(\rho \parallel \sigma )={\rm tr}(\rho ({\rm log} \rho -{\rm log} \sigma ))$ in the last sum of the previous equation, the total work becomes

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)}) & = & F(\tilde{\rho }_{SB}^{(0)},H_{SB}^{(0)})-F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)}) \\ {} & {} & -\left( F(\tilde{\rho }_{SB}^{(0)},H_{SB}^{(1)})-F(\omega (H_{SB}^{(1)}),H_{SB}^{(1)}) \right) \\ {} & {} & -\frac{{\rm ln} 2}{\beta }\mathop{\sum }\limits_{i=1}^{n-1}S(\omega (H_{SB}^{(i+1)})\parallel \omega (H_{SB}^{(i)})). \\ \end{array}\end{eqnarray} \tag{ 14 }$$

Finally, the positivity of the relative entropy and the inequality $F(\omega (H),H)\leqslant F(\rho ,H)$ complete the proof. The existence of a protocol ${{\mathcal{P}}^{*}}$ saturating the bound is shown in appendix A.

Using (14), we can identify which protocol maximizes the work extracted and saturates the bound arbitrarily well. We need to minimize its two negative terms, that is, (i) the second difference of free energies and (ii) the sum of relative entropies. The minimum of (i) is attained by choosing the first quench to the appropriate Hamiltonian ${{\tilde{H}}_{S}}$. The term (ii) can be made arbitrarily small by performing quenches that represent a minimal change of the Hamiltonian between individual applications of thermalization maps, at the expense of performing many of them. This sequence of quenches and thermalizations precisely emulates an isothermal reversible process. Thus, theorem 1 not only introduces a fundamental bound for the maximum extracted work but also tells us which protocol attains that maximum arbitrarily well. These attainable bounds complement the findings presented in [6], in which the impact of correlations to the attainable work extraction is considered for harmonic potentials as well as for weakly anharmonic potentials within a Fokker–Planck approach.□

Note that (9) contains as a particular case the well-known bounds on expected work extraction in the weak-coupling regime [23, 29, 30] (see supplementary material). When V is weak in comparison with the energy gaps of ${{H}_{B}}+{{H}_{S}}$ to an extent that is negligible in an idealized treatment and the thermalization process is such that ${{\rho }_{S}}=\omega ({{H}_{S}})$, then the maximum work extracted is given by the difference of free energies

$$\begin{eqnarray}&&\langle {{W}^{{\rm wc}}}\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)})\leqslant F(\rho _{S}^{(0)},H_{S}^{(0)})-F(\omega (H_{S}^{(0)}),H_{S}^{(0)}).\end{eqnarray} \tag{ 15 }$$

Furthermore, expression (9) has an insightful physical interpretation. We will show that the second line in (9) vanishes if and only if the optimal protocol is reversible. Otherwise, the strong coupling between the system and bath induces an unavoidable dissipation in the thermalization process that makes the protocol irreversible and limits the work that can be extracted.

We call a protocol $\mathcal{P}$ of work extraction reversible if $\langle W\rangle (\mathcal{P},H_{SB}^{(0)},\rho _{S}^{(0)})$ $=-\langle W\rangle ({{\mathcal{P}}^{-1}},H_{SB}^{(0)},\rho _{S}^{(n)})$ where ${{\mathcal{P}}^{-1}}$ just inverts the order of the list and Hamiltonians $\{{{H}^{(i)}}\}_{i=1}^{n-1}$ and thermalizations of $\mathcal{P}$ and $\rho _{S}^{(n)}$ are the final state after applying $\mathcal{P}$ to $\rho _{S}^{(0)}$. In other words, if $\mathcal{P}$ is a protocol that brings the system to equilibrium while extracting work, ${{\mathcal{P}}^{-1}}$ supplies work in order to bring an equilibrium state out of equilibrium. It is well-known that in the weak-coupling case, when the processes are optimal, ${{\mathcal{P}}^{*}}$ and ${{\mathcal{P}}^{*-1}}$ extract/supply the same amount of work. Here we show that this is not the case in the strong-coupling case. One can show, by a similar argument as used in the proof of theorem 1, that

$$\begin{eqnarray}&&\langle W\rangle ({{\mathcal{P}}^{*-1}},H_{S}^{(0)},{{\omega }_{S}}(H_{SB}^{(0)}))=-F({{\tilde{\rho }}_{SB}},H_{SB}^{(0)})+F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)})\ =:\ \Delta {{F}_{{\rm rev}}}.\end{eqnarray} \tag{ 16 }$$

Note that the optimal forward and reversed protocol differ exactly in the second line in (9), which for this reason we refer to as the irreversible free-energy difference

$$\begin{eqnarray*}&&\Delta {{F}_{{\rm irrev}}}\;:=\;-\mathop{{\rm min} }\limits_{{{\tilde{H}}_{S}}}\left[ F({{\tilde{\rho }}_{SB}},{{\tilde{H}}_{SB}})-F(\omega ({{\tilde{H}}_{SB}}),{{\tilde{H}}_{SB}}) \right].\end{eqnarray*}$$

Hence, even a close to optimal protocol is surprisingly in general far from being reversible. The reason for the irreversibility is that when $\rho _{S}^{(0)}$ cannot be expressed as the reduced state of thermal state $\omega ({{\tilde{H}}_{SB}})$, then it is impossible that a protocol ${{\mathcal{P}}^{-1}}$ brings ${{\omega }_{S}}(H_{SB}^{(0)})$ into $\rho _{S}^{(0)}$. This is precisely the case when $\Delta {{F}_{{\rm irrev}}}$ is not zero. Note that in the weak-coupling regime this is never the case, as any state $\rho _{S}^{(0)}$ can be expressed as a thermal state at any temperature, given that one can choose the Hamiltonian. Therefore, in contrast to our case, in the weak-coupling case the optimal protocol is reversible.

The existence of a reversible protocol saturating the work extraction is well known to be related to the saturation of the second law of thermodynamics. Let us recall Clausius' theorem, that in a commonly expressed variant states that

$$\begin{eqnarray}&&\frac{\Delta Q}{T}=\int _{i}^{f}\frac{\delta Q}{T}\leqslant \Delta S={{S}_{f}}-{{S}_{i}},\end{eqnarray} \tag{ 17 }$$

where Q is the heat defined as the energy lost by the bath and S is the thermodynamic entropy. Most importantly, equality (saturation of the second law) holds only when the process is reversible. If one relates the thermodynamic entropy with the von Neumann entropy, (17) it can easily be shown to imply the bounds of work extraction in the weak-coupling regime (15), where the bound is indeed saturated for reversible processes.

In the strong-coupling regime, one can see that Clausius' version of the second law (17) implies that

$$\begin{eqnarray}&&\langle W\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)})\leqslant -\Delta {{F}_{{\rm rev}}},\end{eqnarray} \tag{ 18 }$$

differing from the bound of theorem 1 precisely in $\Delta {{F}_{{\rm irrev}}}$. This clarifies the role of the strong coupling in thermodynamics: the entanglement between the system–bath induces unavoidable irreversibility that is an obstacle against saturating the second law of thermodynamics. It is only for particular initial states (those that look like reduced states of the thermal states of a larger system) that reversible protocols can be implemented and the second law can be saturated. This striking limiting effect of entanglement contrasts with previous works in alternative scenarios [32, 33], where entanglement was regarded as an enhancer of work extraction or power.

The bounds on work extraction of our formalism coincide, in the special case of a weak-coupling regime (15), with previous results that employed a different set of operations based on unitary transformations [28–31]. There, optimal protocols employed system–bath interactions mediated by fine-tuned unitaries that differ substantially from what one would expect nature to implement generically. On the contrary, in our formalism the system–bath coupling is only required to thermalize the system following (2), being arguably the case for most interactions. This explains the ubiquity of work extraction machines which are far from needing microscopically engineered unitaries. Here, in order to connect our work with this other approach, we formulate an embedding of our set of operations into a unitary formalism.

5.1. Quenches in a unitary formulation

The standard way of describing a quench is by the sudden change of a parameter of the Hamiltonian. By solving the time-dependent Schrödinger equation it can be proven that if the change of such a parameter is performed fast enough, the state of the system remains the same immediately before and after the quench. Nevertheless, as we show next, this description of the quench has strong implications on the properties of the systemʼs environment.

Let us consider a two-level system R with energy eigenstates $|0\rangle$ and $|1\rangle$ and energy levels E₀ and E₁, respectively. In order to perform a quench, it is also necessary to consider an environment or battery with Hamiltonian H_W that supplies (stores) the lack (excess) of energy required by the quench. The Hamiltonian of the whole set-up is then $H={{H}_{R}}+{{H}_{W}}$. Let us now consider a unitary process U that performs a level transformation of the excited state $|1\rangle$ from E₁ to ${{E}_{1}}+\Delta$. More explicitly,

$$\begin{eqnarray*}\begin{array}{rcl} |0{{\rangle }_{R}}|0{{\rangle }_{W}}\ \mapsto \ U|0{{\rangle }_{R}}|0{{\rangle }_{W}} & = & |0{{\rangle }_{R}}|0{{\rangle }_{W}} \\ |1{{\rangle }_{R}}|0{{\rangle }_{W}}\ \mapsto \ U|1{{\rangle }_{R}}|0{{\rangle }_{W}} & = & |1{{\rangle }_{R}}|0-\Delta {{\rangle }_{W}}. \\ \end{array}\end{eqnarray*}$$

Because of the linearity of the unitary that implements the quench, we can also transform an initial superposition state. For instance, the $|+{{\rangle }_{R}}$ becomes

$$\begin{eqnarray}&&\left( |0{{\rangle }_{R}}+|1{{\rangle }_{R}} \right)|0{{\rangle }_{W}}\ \mapsto \ |0{{\rangle }_{R}}|0{{\rangle }_{W}}+|1{{\rangle }_{R}}|0-\Delta {{\rangle }_{W}}.\end{eqnarray} \tag{ 19 }$$

Hence, while the initial state of the set RW was a product state, the state after the quench is entangled in a superposition for the battery of having and not having supplied energy. This implies that it is impossible to perform quenches that leave the system unchanged if the battery is initially in an energy eigenstate. Although this conclusion seems a bit odd, it can be circumvented by having a battery with non-distinguishable states.

In order to clarify the above situation, let us think of a partition of the experiment into a system R, a battery W and a control qubit Q. The total Hamiltonian of the system is

$$\begin{eqnarray}&&H=H_{R}^{(0)}\otimes {{\mathbb{I}}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}}+H_{R}^{(1)}\otimes {{\mathbb{I}}_{W}}\otimes |1\rangle \langle 1{{|}_{Q}}+{{\mathbb{I}}_{S}}\otimes {{H}_{W}}\otimes {{\mathbb{I}}_{Q}}\end{eqnarray} \tag{ 20 }$$

where $H_{R}^{(0,1)}$ are arbitrary Hamiltonians with eigenvectors $\{|{{i}^{(0,1)}}\rangle {{\}}_{i}}$ and eigenvalues ${{\{E_{i}^{(0,1)}\}}_{i}}$ with $i=1,...,{{d}_{R}}$, and d_R is the dimension of the Hilbert space ${{\mathcal{H}}_{R}}$ of system R. Note that the control qubit Q dictates which is the Hamiltonian of the system.

We consider now the action of a global unitary U supported on RWQ on an initially uncorrelated state $\rho _{RWQ}^{(0)}=\rho _{R}^{(0)}\otimes \rho _{W}^{(0)}\otimes |0\rangle \langle 0{{|}_{Q}}$, in a way such that the final state can be written as

$$\begin{eqnarray}&&\rho _{RWQ}^{(1)}=U\rho _{RWQ}^{(1)}{{U}^{\dagger }}=\rho _{RW}^{(1)}\otimes |1\rangle \langle 1{{|}_{Q}}.\end{eqnarray} \tag{ 21 }$$

In this way, according to (20), the effective Hamiltonian acting on R has changed from $H_{R}^{(0)}$ to $H_{R}^{(1)}$. In addition, we impose the following three natural constraints on the unitary transformation:

• (i)  
  Energy conservation. The unitary commutes with the Hamiltonian $[H,U]=0$.
• (ii)  
  Battery translational invariance. We assume the battery to have a sufficiently dense equally spaced spectrum $\mathcal{W}$ (or a continuous one) with the property that the unitary U commutes with ${{\mathbb{1}}_{R}}\otimes {{\Gamma }_{W}}(e)$ for all $w,w+e\in \mathcal{W}$, where ${{\Gamma }_{W}}(e)|w\rangle =|w+e\rangle$ is the translation operator on the battery. This merely reflects the invariance of the transformation under changes of the energy-origin of H_W [30].
• (iii)  
  Quenches. The unitary transformation is such that there exists an initial state of the battery $\rho _{W}^{(0)}$ is such that $\rho _{R}^{(1)}=\rho _{R}^{(0)}$ for every initial state $\rho _{R}^{(0)}$ and any Hamiltonian transformation $H_{R}^{(0)}\mapsto H_{R}^{(1)}$.

Conditions (i) and (ii) are not present in the abstract formalism of work extraction of section 2. We highlight that this is a desired feature of our approach: the general formalism that provides the above mentioned bounds avoids as many assumptions as possible for the sake of general applicability. However, the particular protocol that attains the maximum fulfils further conditions of physical relevance. In particular, assumption (i) allows one to extend this analysis to a single-shot work extraction, as considered in [23, 25, 28]. We leave these analyses open for further work. The following theorem shows that under these additional assumptions, the unitary performing the quench exists and is unique.

(Uniqueness of unitary realizations).

Theorem 2 Consider unitary transformations such that $U{{\rho }^{(0)}}{{U}^{\dagger }}=\rho _{RW}^{(1)}\otimes |1\rangle \langle 1|$ for any initial state ${{\rho }^{(0)}}=\rho _{RW}^{(0)}\otimes |0\rangle \langle 0|$. The unitary that fulfils conditions (i)–(iii) is unique and can be written as

$$\begin{eqnarray}\begin{array}{rcl} U & = & \mathop{\sum }\limits_{i,j}|{{j}^{(1)}}\rangle \langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}{{|}_{R}}\otimes {{\Gamma }_{W}}(E_{i}^{(0)}-E_{j}^{(1)})\otimes |1\rangle \langle 0{{|}_{Q}} \\ {} & {} & +\mathop{\sum }\limits_{i,j}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|{{j}^{(1)}}\rangle \langle {{j}^{(1)}}{{|}_{R}}\otimes {{\Gamma }_{W}}(E_{j}^{(1)}-E_{i}^{(0)})\otimes |0\rangle \langle 1{{|}_{Q}}. \\ \end{array}\end{eqnarray} \tag{ 22 }$$

Proof. The fact that the unitary flips the state of Q implies that

$$\begin{eqnarray}&&U=U_{RW}^{{\rm on}}\otimes |1\rangle \langle 0{{|}_{Q}}+U_{RW}^{{\rm off}}\otimes |0\rangle \langle 1{{|}_{Q}}.\end{eqnarray} \tag{ 23 }$$

Let us first consider the case in which the initial state $\rho _{RWQ}^{(0)}=|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\otimes |w\rangle \langle w|\otimes |0\rangle \langle 0|$ is an eigenstate of H with energy $E_{i}^{(0)}+w$. By using condition (i), the final state after performing the unitary is also an eigenstate of H with the same energy,

$$\begin{eqnarray}&&(H_{R}^{(1)}+{{H}_{W}})U_{RW}^{{\rm on}}|{{i}^{(0)}}\rangle \otimes |w\rangle =(E_{i}^{(0)}+w)U_{RW}^{{\rm on}}|{{i}^{(0)}}\rangle \otimes |w\rangle .\end{eqnarray} \tag{ 24 }$$

Hence, the state $U_{RW}^{{\rm on}}|{{i}^{(0)}}\rangle \otimes |w\rangle$ is contained in the sub-space spanned by $\{|{{j}^{(1)}}\rangle \otimes |E_{i}^{(0)}+w-E_{j}^{(1)}\rangle {{\}}_{j}}$, that is

$$\begin{eqnarray}&&U_{RW}^{{\rm on}}|{{i}^{(0)}}\rangle \otimes |w\rangle =\mathop{\sum }\limits_{j}{{R}_{j,i,w}}|{{j}^{(1)}}\rangle \otimes |{{\Delta }_{i,j}}+w\rangle ,\end{eqnarray} \tag{ 25 }$$

where ${{\Delta }_{i,j}}\;:=\;E_{i}^{(0)}-E_{j}^{(1)}$ and ${{R}_{j,i,w}}$ are its coefficients.

By imposing condition (ii), $[U_{RW}^{{\rm on}},{{\mathbb{I}}_{R}}\otimes {{\Gamma }_{W}}(E)]=0$, one gets

$$\begin{eqnarray}&&\mathop{\sum }\limits_{j}({{R}_{j,i,w}}-{{R}_{j,i,w+E}})|{{j}^{(1)}}\rangle \otimes |{{\Delta }_{i,j}}+w+E\rangle =0,\end{eqnarray} \tag{ 26 }$$

which implies that the ${{R}_{j,i,w}}$ are independent of w, hence, ${{R}_{j,i,w}}={{R}_{j,i}}$.

In order to exploit condition (iii), let us consider that the state of the battery that allows for quenches, i.e. $\rho _{R}^{(1)}=\rho _{R}^{(0)}$ for any $\rho _{R}^{(0)}$, is pure and denoted by $|{{\Psi }^{(0)}}{{\rangle }_{W}}={{\sum }_{w}}{{B}_{w}}|w\rangle$. We choose the global initial state to be $\rho _{RWQ}^{(0)}=|\phi \rangle \langle \phi {{|}_{R}}\otimes |{{\Psi }^{(0)}}\rangle \langle {{\Psi }^{(0)}}{{|}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}}$, with $|\phi \rangle ={{\sum }_{i}}{{c}_{i}}|{{i}^{(0)}}\rangle$. Then,

$$\begin{eqnarray}&&U_{RW}^{{\rm on}}|\phi \rangle \otimes |{{\Psi }^{(0)}}\rangle =\mathop{\sum }\limits_{j,i}{{c}_{i}}{{R}_{j,i}}|{{j}^{(1)}}\rangle \otimes \Gamma ({{\Delta }_{i,j}})|{{\Psi }^{(0)}}\rangle .\end{eqnarray} \tag{ 27 }$$

From (27), we can compute the reduced state of R

$$\begin{eqnarray}&&\rho _{R}^{(1)}=\mathop{\sum }\limits_{i,i^{\prime} ,j,j^{\prime} }{{c}_{i}}c_{i^{\prime} }^{*}{{R}_{j,i}}R_{j^{\prime} ,i^{\prime} }^{*}K_{j,j^{\prime} }^{i,i^{\prime} }|{{j}^{(1)}}\rangle \langle j{{^{\prime} }^{(1)}}|,\end{eqnarray} \tag{ 28 }$$

where $K_{j,j^{\prime} }^{i,i^{\prime} }\;:=\;\langle {{\Psi }^{(0)}}|\Gamma ({{\Delta }_{i,j}}-{{\Delta }_{i^{\prime} ,j^{\prime} }})|{{\Psi }^{(0)}}\rangle$. Imposing $\rho _{R}^{(1)}=|\phi \rangle \langle \phi {{|}_{R}}$, we obtain

$$\begin{eqnarray}&&{{R}_{j,i}}R_{j^{\prime} ,i^{\prime} }^{*}K_{j,j^{\prime} }^{i,i^{\prime} }=\langle j|i\rangle \langle i^{\prime} |j^{\prime} \rangle ,\quad \forall j,\ j^{\prime} ,\ i,\ i^{\prime} .\end{eqnarray} \tag{ 29 }$$

Multiplying the previous equation by its conjugate and summing over j and $j^{\prime}$, one obtains

$$\begin{eqnarray}&&\mathop{\sum }\limits_{j,j^{\prime} }|{{R}_{j,i}}{{|}^{2}}|{{R}_{j^{\prime} ,i^{\prime} }}{{|}^{2}}|K_{j,j^{\prime} }^{i,i^{\prime} }{{|}^{2}}=1=\mathop{\sum }\limits_{j,j^{\prime} }|{{R}_{j,i}}{{|}^{2}}|{{R}_{j^{\prime} ,i^{\prime} }}{{|}^{2}},\quad \forall i,i^{\prime} .\end{eqnarray} \tag{ 30 }$$

Notice that because of the $\Gamma (x)$ being unitary, $|K_{j,j^{\prime} }^{i,i^{\prime} }|\leqslant 1$. Hence, the only way that condition (iii) can hold true requires that

$$\begin{eqnarray}&&|K_{j,j^{\prime} }^{i,i^{\prime} }|=|\langle {{\Psi }^{(0)}}|\Gamma ({{\Delta }_{i,j}}-{{\Delta }_{i^{\prime} ,j^{\prime} }})|{{\Psi }^{(0)}}\rangle |=1,\quad \forall j,j^{\prime} ,i,i^{\prime} .\end{eqnarray} \tag{ 31 }$$

This can be only satisfied for every choice of $E_{j}^{(1)},E_{i}^{(0)}$, (i.e., for every choice of energy level of the initial and final Hamiltonian) if ${{B}_{w}}={{B}_{w+\Delta E}}$ for every possible value of

$$\begin{eqnarray}&&\Delta E\leqslant \mathop{{\rm max} }\limits_{j,j^{\prime} }|{{\Delta }_{i,j}}-{{\Delta }_{i^{\prime} ,j^{\prime} }}|,\end{eqnarray} \tag{ 32 }$$

which in turn implies that $K_{j,j^{\prime} }^{i,i^{\prime} }=1$ $\forall j,j^{\prime} ,i,i^{\prime}$. This, together with (28), implies that ${{R}_{j,i}}=\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle$. This leads to

$$\begin{eqnarray}&&U_{RW}^{{\rm on}}=\mathop{\sum }\limits_{i,j}|{{j}^{(1)}}\rangle \langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}{{|}_{R}}\otimes {{\Gamma }_{W}}(E_{i}^{(0)}-E_{j}^{(1)}).\end{eqnarray} \tag{ 33 }$$

This argument can be straightforwardly extended for the case of a mixed state of the battery,

$$\begin{eqnarray}&&\rho _{W}^{(0)}=\mathop{\sum }\limits_{w}{{p}_{w}}|\Psi _{w}^{(0)}\rangle \langle \Psi _{w}^{(0)}|.\end{eqnarray} \tag{ 34 }$$

Also, a symmetric argument can be applied to $U_{RW}^{{\rm off}}$ by considering an inverse quench $H_{R}^{(1)}\mapsto H_{R}^{(0)}$ that must leave invariant the initial state of R. Altogether, we arrive at (22). □

One observation of the previous proof is that in order for the unitary to keep the state of the system invariant, the battery must be in an energy coherent state with a much larger uncertainty than the operator norm of the Hamiltonians H_R. This is encapsulated in (31). Hence, coherence is a resource required to implement quenches. Contrary to the conclusions of [31], coherence is destroyed due to time evolutions of the battery with H_W (see the supplementary material). This suggests that the catalytic role of coherence in [31] may be a consequence of disregarding time evolution as the mechanism for thermalization.

5.2. Thermalization after the quench driven by the unitary time evolution

In this work we have introduced a framework to study work extraction in thermal machines. Our formalism considers quantum Hamiltonian quenches as the fundamental operations and analyses the effect of strong couplings between the system and the thermal bath. Strikingly, system–bath entanglement seriously limits the amount of work extractable and induces irreversibility in the process, which in turn prevents one from saturating the second law of thermodynamics. This is relevant since any finite-time approach to quantum thermodynamics necessarily has to take correlations and non-zero interactions into account. Also, we introduce a formalism to embed Hamiltonian quenches into a unitary formalism. Under a set of reasonable assumptions, we show that the unitary embedding is unique and coherence is required as a resource to implement the quenches. It should be clear that the mindset presented here can also be applied to a variety of related problems in quantum thermodynamics such as Landauerʼs principle [34, 35], whenever correlations are expected to be non-negligible. Hence, this work opens new venues to understand the role of quantum effects such as entanglement and coherence in thermodynamics.

We would like to thank the EU (Q-Essence, SIQS, RAQUEL, COST), the ERC (TAQ), the BMBF, the FQXi, and the AvH for support.

Here, we show that the bound (9) can be arbitrarily well approximated. This can be most easily seen in a 'continuum limit' of protocols, where an arbitrarily large number n of operations are performed. The first step of the protocol that arbitrarily well saturates the bound (9) is to perform a quench on SBW from ${{H}^{(0)}}$ to ${{H}^{(1)}}=H_{S}^{(1)}+{{H}_{B}}+V+{{H}_{W}}$, where $H_{S}^{(1)}$ is the Hamiltonian that attains the minimum in the second term of (9). Applying (4) straightforwardly one finds

$$\begin{eqnarray}&&{{\langle W\rangle }^{0,1}}={\rm t}{{{\rm r}}_{SB}}((H_{SB}^{(0)}-H_{SB}^{(1)})\rho _{SB}^{(0)})={\rm t}{{{\rm r}}_{S}}((H_{S}^{(0)}-H_{S}^{(1)})\rho _{S}^{(0)}).\end{eqnarray} \tag{ A.1 }$$

Consider now a differentiable parametrized curve ${{H}_{S}}:[0,1]\to \mathcal{B}({{\mathcal{H}}_{S}})$, where $\mathcal{B}({{\mathcal{H}}_{S}})$ denotes the bounded operators on the Hilbert space associated with SB. This function fulfils ${{H}_{S}}(0)=H_{S}^{(1)}$ and ${{H}_{S}}(1)=H_{S}^{(0)}$. Given an integer n, one defines a sequence of $n-1$ Hamiltonians as

$$\begin{eqnarray}&&H_{SB}^{(i)}\;:=\;{{H}_{S}}\left( \frac{i-1}{n-1} \right)+V+{{H}_{B}}\end{eqnarray} \tag{ A.2 }$$

with $i=1,...,n$. This sequence of Hamiltonians will be used as a sequence of quenches on the equilibrated sub-system, as discussed in section appendix F. More precisely, consider a protocol in which, after the first quench from ${{H}^{(0)}}$ to ${{H}^{(1)}}$ described above, one applies a sequence of state thermalizations as (7) followed by quenches $H_{SB}^{(i)}\mapsto H_{SB}^{(i+1)}$ with $i=1,...,n-1$. One finds that

$$\begin{eqnarray}&&{{\langle W\rangle }^{1,n}}=\mathop{\sum }\limits_{i=1}^{n-1}{{\langle W\rangle }^{i,i+1}}=\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}\left( (H_{S}^{(i)}-H_{S}^{(i+1)})\otimes {{\mathbb{I}}_{B}}\omega (H_{SB}^{(i)}) \right).\end{eqnarray} \tag{ A.3 }$$

In the limit of n tending to infinity, the expected work cost of these sequence of quenches can be written as

$$\begin{eqnarray}&&\mathop{{\rm lim} }\limits_{n\to \infty }{{\langle W\rangle }^{1,n}}=-\int _{0}^{1}{\rm tr}\left( \frac{\partial {{H}_{S}}(\lambda )}{\partial \lambda }\otimes {{\mathbb{I}}_{B}}\omega ({{H}_{S}}(\lambda )+V+{{H}_{B}}) \right){\rm d}\lambda .\end{eqnarray} \tag{ A.4 }$$

Let us denote ${{H}_{SB}}(\lambda )={{H}_{S}}(\lambda )+V+{{H}_{B}}$, then

$$\begin{eqnarray}&&\begin{array}{l} \frac{\partial }{\partial \lambda }{\rm ln} ({\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})) \\ \quad =\frac{\frac{\partial }{\partial \lambda }{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})}=\frac{{\rm tr}(\frac{\partial }{\partial \lambda }{{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})} \\ \quad =\frac{{\rm tr}(\int _{0}^{1}{{{\rm e}}^{-\alpha \beta {{H}_{SB}}(\lambda )}}\frac{\partial }{\partial \lambda }(-\beta {{H}_{SB}}(\lambda )){{{\rm e}}^{-(1-\alpha )\beta {{H}_{SB}}(\lambda )}}){\rm d}\alpha }{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})} \\ \end{array}\end{eqnarray} \tag{ A.5 }$$

$$\begin{eqnarray}\begin{array}{clc} {} & = & \frac{\int _{0}^{1}{\rm tr}({{{\rm e}}^{-\alpha \beta {{H}_{SB}}(\lambda )}}\frac{\partial }{\partial \lambda }(-\beta {{H}_{SB}}(\lambda )){{{\rm e}}^{-(1-\alpha )\beta {{H}_{SB}}(\lambda )}}){\rm d}\alpha }{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})} \\ {} & = & \frac{{\rm tr}(\frac{\partial }{\partial \lambda }(-\beta ({{H}_{SB}}(\lambda )){{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}(\lambda )}})} \\ \end{array}\end{eqnarray} \tag{ A.6 }$$

$$\begin{eqnarray}&&=-\beta {\rm tr}\left( \frac{{\rm d}{{H}_{S}}(\lambda )}{{\rm d}\lambda }\otimes {{\mathbb{I}}_{B}}\omega ({{H}_{S}}(\lambda )+V+{{H}_{B}}) \right),\end{eqnarray} \tag{ A.7 }$$

where (A.5) follows from the Wilcox formula for matrix exponential derivatives [36]. By combining (A.7) with (A.4) and $F(\omega (H),H)=-{\rm ln} ({\rm tr}({{{\rm e}}^{-\beta H}}))/\beta$, one finds

$$\begin{eqnarray}&&\mathop{{\rm lim} }\limits_{n\to \infty }{{\langle W\rangle }^{1,n}}=F(\omega (H_{SB}^{(1)}),H_{SB}^{(1)}))-F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)})),\end{eqnarray} \tag{ A.8 }$$

hence, the total work extracted in the process is

$$\begin{eqnarray}\begin{array}{rcl} \mathop{{\rm lim} }\limits_{n\to \infty }{{\langle W\rangle }^{1,n}}({{H}^{(0)}},\rho _{S}^{(0)}) & = & {\rm t}{{{\rm r}}_{S}}((H_{S}^{(0)}-H_{S}^{(1)})\rho _{S}^{(0)})+F(\omega (H_{SB}^{(1)}),H_{SB}^{(1)})) \\ {} & {} & -F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)})) \\ \end{array}\end{eqnarray} \tag{ A.9 }$$

$$\begin{eqnarray}&&=\;F({{\tilde{\rho }}_{SB}},H_{SB}^{(0)})-F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)})-F({{\tilde{\rho }}_{SB}},H_{SB}^{(1)})-F(\omega (H_{SB}^{(1)}),H_{SB}^{(1)})\end{eqnarray} \tag{ A.10 }$$

$$\begin{eqnarray}&&=\;F({{\tilde{\rho }}_{SB}},H_{SB}^{(0)})-F(\omega (H_{SB}^{(0)}),H_{SB}^{(0)})-\mathop{{\rm min} }\limits_{{{\tilde{H}}_{S}}}\left[ F({{\tilde{\rho }}_{SB}},{{\tilde{H}}_{SB}})-F(\omega ({{\tilde{H}}_{SB}}),{{\tilde{H}}_{SB}}) \right],\end{eqnarray} \tag{ A.11 }$$

where (A.10) follows from calculations equivalent to the ones from the proof of theorem 1 and (A.11) is derived from the choice of $H_{SB}^{(1)}$.

In the weak-coupling limit the effect of the bath is to drive the system S to an equilibrium state $\omega ({{H}_{S}})$ that is a Gibbs state, since

$$\begin{eqnarray}&&\omega ({{H}_{SB}})\approx \omega ({{H}_{S}})\otimes \omega ({{H}_{B}}),\end{eqnarray} \tag{ B.1 }$$

then one can take ${{\tilde{\rho }}_{SB}}=\omega ({{\tilde{H}}_{S}})\otimes \omega ({{H}_{B}})$ and ${{\tilde{H}}_{S}}=-{\rm ln} (\rho _{S}^{(0)})/\beta$ and a simple calculation shows that in that case $\Delta {{F}_{{\rm irrev}}}=0$ and $\Delta {{F}_{{\rm rev}}}$ does not depend on H_B, so that

$$\begin{eqnarray}&&\mathop{{\rm max} }\limits_{\mathcal{P}}\langle W\rangle (\mathcal{P},{{H}^{(0)}},\rho _{S}^{(0)})\leqslant F(\rho _{S}^{(0)},H_{S}^{(0)})-F(\omega (H_{S}^{(0)}),H_{S}^{(0)}).\end{eqnarray} \tag{ B.2 }$$

Let us now comment on the precise role of the two terms $\Delta {{F}_{{\rm rev}}}$ and $\Delta {{F}_{{\rm irrev}}}$ appearing in the bound (9), as defined in (16) and (17). Consider first a modified scenario in which SB is treated as a larger working medium that we denote by $S^{\prime}$. In such scenario, one has full control over the Hamiltonian of $S^{\prime}$, that is H_SB, and furthermore, that $S^{\prime}$ can be driven to the Gibbs equilibrium state $\omega ({{H}_{SB}})$—this may be achieved by weak coupling with a bath $B^{\prime}$ that interacts with $SB$. In this case, similar analysis to the one leading to (B.2) shows that the maximum work extracted is precisely $\Delta {{F}_{{\rm rev}}}$. Hence, $\Delta {{F}_{{\rm irrev}}}$ should be understood as a work penalty due to our lack of control over ${{H}_{B}}$ and, therefore, through expression (7), over the equilibrium state of $S$.

From the proof of theorem 2 it is clear that one needs a specific initial state of the battery ${{\rho }_{W}}(0)=|{{\Psi }^{(0)}}\rangle \langle \Psi (0)|$ in order to guarantee that the state of R is not altered by the change of Hamiltonian. This is encapsulated in the following condition

$$\begin{eqnarray}&&K_{j,j^{\prime} }^{i,i^{\prime} }\;:=\;\langle {{\Psi }^{(0)}}|\Gamma ({{\Delta }_{i,j}}-{{\Delta }_{i^{\prime} ,j^{\prime} }})|{{\Psi }^{(0)}}\rangle =1\forall i,i^{\prime} ,j,j^{\prime} .\end{eqnarray} \tag{ C.1 }$$

This condition can be achieved by employing an initial state vector of the battery $|{{\Psi }^{(0)}}\rangle =|{{\Psi }_{\sqcap }}\rangle$ with

$$\begin{eqnarray}&&|{{\Psi }_{\sqcap }}\rangle =\frac{1}{N(E_{W}^{(0)},\Delta )}\sum W(0)E_{W}^{(0)}+\Delta w=E|w\rangle ,\end{eqnarray} \tag{ C.2 }$$

where $N(E_{W}^{(0)},\Delta )$ is the number of states with energy between $E_{W}^{(0)}$ and $E_{W}^{(0)}+\Delta$, according to the discretization chosen. Then,

$$\begin{eqnarray}&&K_{j,j^{\prime} }^{i,i^{\prime} }=\frac{\Delta -|E_{j^{\prime} }^{(1)}-E_{j}^{(1)}+E_{i}^{(0)}-E_{i^{\prime} }^{(0)}|}{\Delta }.\end{eqnarray} \tag{ C.3 }$$

Therefore, by assuming

$$\begin{eqnarray}&&\mathop{{\rm max} }\limits_{i,i^{\prime} ,j^{\prime} ,j}\frac{|E_{j^{\prime} }^{(1)}-E_{j}^{(1)}+E_{i}^{(0)}-E_{i^{\prime} }^{(0)}|}{\Delta }\leqslant \epsilon \end{eqnarray} \tag{ C.4 }$$

by taking Δ so that $\epsilon \gt 0$ is arbitrarily small we obtain ${{K}_{j,j^{\prime} }}$ arbitrarily close to one.

Let us now analyse how the state of the battery is changed after the quench from $H_{R}^{(0)}$ to $H_{R}^{(1)}$. Starting from an initial state

$$\begin{eqnarray}&&\rho _{RWQ}^{(0)}=\rho _{R}^{(0)}\otimes |{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}{{|}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}},\end{eqnarray} \tag{ C.5 }$$

applying (22) one finds that

$$\begin{eqnarray}\begin{array}{rcl} \rho _{RW}^{(1)} & = & \mathop{\sum }\limits_{i,i^{\prime} ,j,j^{\prime} }\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|j{{^{\prime} }^{(1)}}\rangle |{{j}^{(1)}}\rangle \\ {} & {} & \times \langle j{{^{\prime} }^{(1)}}{{|}_{R}}\Gamma ({{\Delta }_{i,j}})|{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}|{{\Gamma }^{\dagger }}({{\Delta }_{i^{\prime} ,j^{\prime} }}). \\ \end{array}\end{eqnarray} \tag{ C.6 }$$

Let us define the expected work extracted in the process as the mean energy difference between the initial and the final state of the battery. Then,

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle & := & {\rm t}{{{\rm r}}_{W}}({{H}_{W}}(\rho _{W}^{(1)}-\rho _{W}^{(0)})) \\ {} & = & {\rm t}{{{\rm r}}_{W}}({{H}_{W}}\mathop{\sum }\limits_{i,i^{\prime} ,j}\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|{{j}^{(1)}}\rangle \Gamma ({{\Delta }_{i,j}})|{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}|{{\Gamma }^{\dagger }}({{\Delta }_{i^{\prime} ,j}})) \\ {} & {} & -{\rm t}{{{\rm r}}_{W}}({{H}_{W}}|{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}|) \\ {} & = & \mathop{\sum }\limits_{i,i^{\prime} ,j}\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|{{j}^{(1)}}\rangle \frac{1}{N(E_{W}^{(0)},\Delta )} \\ {} & {} & \times \mathop{\sum }\limits_{E_{W}^{(0)}+\Delta }^{w,w^{\prime} =E_{W}^{(0)}}\mathop{\sum }\limits_{\infty }^{e=-\infty }e\langle e|{{\Delta }_{i,j}}+w\rangle \langle {{\Delta }_{i^{\prime} ,j}}+w^{\prime} |e\rangle -{\rm t}{{{\rm r}}_{W}}({{H}_{W}}|{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}|). \\ \end{array}\end{eqnarray} \tag{ C.7 }$$

From condition (C.4) in the limit $\epsilon \to 0$, we get

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle & = & \mathop{\sum }\limits_{i,i^{\prime} ,j}\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|{{j}^{(1)}}\rangle \frac{1}{N(E_{W}^{(0)},\Delta )}\mathop{\sum }\limits_{E_{W}^{(0)}+\Delta }^{w=E_{W}^{(0)}}({{\Delta }_{i,j}}+w) \\ {} & {} & -{\rm t}{{{\rm r}}_{W}}({{H}_{W}}|{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}|) \\ \end{array}\end{eqnarray} \tag{ C.8 }$$

$$\begin{eqnarray}&&=\mathop{\sum }\limits_{i}\langle {{i}^{(0)}}|\rho _{R}^{(0)}|{{i}^{(0)}}\rangle E_{i}^{(0)}-\mathop{\sum }\limits_{i,i^{\prime} ,j}\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|{{j}^{(1)}}\rangle E_{j}^{(1)}\end{eqnarray} \tag{ C.9 }$$

$$\begin{eqnarray}&&=\;{\rm t}{{{\rm r}}_{R}}((H_{R}^{(0)}-H_{R}^{(1)})\rho _{R}^{(0)}),\end{eqnarray} \tag{ C.10 }$$

where (C.9) follows from the fact that

$$\begin{eqnarray}&&{\rm t}{{{\rm r}}_{W}}({{H}_{W}}|{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}|)=\frac{1}{N(E_{W}^{(0)},\Delta )}\sum W(0)E_{W}^{(0)}w=Ew.\end{eqnarray} \tag{ C.11 }$$

In short, (C.10) formalizes the intuition that the expected energy provided (stored) by the battery is just the expected energy gained (lost) by the system R when the quench is applied. Indeed, (C.10) can be derived straightforwardly from the conservation of expected energy of RWQ and the fact that the state of R does not change. However, we derive it explicitly for a consistency check, and also as an illustrative example of how to deal with similar calculations that appear in later sections.

The unitary (22) is the transformation that changes the effective Hamiltonian acting on R, while leaving the state invariant. As shown in previous sections, a sufficiently coherent initial state of the battery is necessary to perform such a transformation. Here, we study what the effect of the unitary (22) is if the initial state of the battery is a classical state. We will show how the state of R is indeed disturbed when one implements that change of the Hamiltonian and how it relates to the work extracted by the battery in such a process. Let us consider an initial state

$$\begin{eqnarray}&&\rho _{RWQ}^{(0)}=\rho _{R}^{(0)}\otimes |0\rangle \langle 0{{|}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}}.\end{eqnarray} \tag{ D.1 }$$

We choose the battery to be initialized in the state $|0\rangle \langle 0{{|}_{W}}$ for ease of notation, but the extension to other pure initial states or convex mixtures of eigenstates of H_W is straightforward. The final state of RW after the quench is

$$\begin{eqnarray}&&\rho _{RW}^{(1)}=\mathop{\sum }\limits_{i,i^{\prime} ,j,j^{\prime} }\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|j{{^{\prime} }^{(1)}}\rangle |{{j}^{(1)}}\rangle \langle j{{^{\prime} }^{(1)}}{{|}_{R}}\otimes |{{\Delta }_{i,j}}\rangle \langle {{\Delta }_{i^{\prime} ,j^{\prime} }}{{|}_{W}}.\end{eqnarray} \tag{ D.2 }$$

The final state of the system R will depend heavily on the degeneracies of both $H_{R}^{(0)}$ and $H_{R}^{(1)}$, and also on the degeneracies of the energy differences ${{\Delta }_{i,j}}$. Let us assume that the initial state is diagonal in the eigenbasis of $H_{R}^{(0)}$. That is

$$\begin{eqnarray}&&\rho _{R}^{(0)}=\mathop{\sum }\limits_{i}\langle {{i}^{(0)}}|\rho _{R}^{(0)}|{{i}^{(0)}}\rangle |{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|.\end{eqnarray} \tag{ D.3 }$$

In this case

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle & = & {\rm t}{{{\rm r}}_{W}}({{H}_{W}}\rho _{W}^{(1)}) \\ {} & = & {\rm t}{{{\rm r}}_{W}}({{H}_{W}}\mathop{\sum }\limits_{i,i^{\prime} ,j}\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|{{j}^{(1)}}\rangle |{{\Delta }_{i,j}}\rangle \langle {{\Delta }_{i,j}}{{|}_{W}}) \\ {} & = & \mathop{\sum }\limits_{i,j}\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|{{j}^{(1)}}\rangle (E_{i}^{(0)}-E_{j}^{(0)}) \\ {} & = & \mathop{\sum }\limits_{i}\langle {{i}^{(0)}}|\rho _{R}^{(0)}|{{i}^{(0)}}\rangle E_{i}^{(0)}-\mathop{\sum }\limits_{j}\langle {{j}^{(1)}}|(\mathop{\sum }\limits_{i}\langle {{i}^{(0)}}|\rho _{R}^{(0)}|{{i}^{(0)}}\rangle |{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|)|{{j}^{(1)}}\rangle E_{j}^{(0)} \\ {} & = & {\rm t}{{{\rm r}}_{R}}((H_{R}^{(0)}-H_{R}^{(1)})\rho _{R}^{(0)}). \\ \end{array}\end{eqnarray} \tag{ D.4 }$$

Note that condition (D.4) is a necessary condition for the set of operations of the work-extracting protocol. Therefore, for classical states of the battery, the quench formalism can only be applied to extract work if the initial state $\rho _{R}^{(0)}$ is diagonal.

We now turn to the discussion of the physical mechanism that renders the thermalization map plausible. Indeed, it captures what one naturally would expect when bringing a small body into contact with a heat bath. In the above axiomatic approach we again leave the mechanism unspecified; here, we will explain why the above framework is indeed very meaningful and physically plausible. In in one way or another, the evolution to an equilibrium Gibbs state is essential in the functioning of any thermal machine. The precise setting considered, however, varies within recent approaches to the study of thermal machines. Within the formalism presented in [23, 25] a classical system is put in contact with a thermal bath. The system is classical in the sense that it is described a state ${{\sigma }_{S}}={{\sum }_{j}}{{\sigma }_{j}}|j\rangle \langle j|$ that is diagonal at all times, where $\{|j\rangle \}$ denotes the eigenvectors of a Hamiltonian H_S in a given state of the process. The evolution towards the Gibbs state in this formalism states that the probability distribution is modified and eventually reaches an equilibrium state given by

$$\begin{eqnarray}&&\omega ({{H}_{S}})=\frac{{{{\rm e}}^{-\beta {{H}_{S}}}}}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{S}}}})}.\end{eqnarray} \tag{ E.1 }$$

An alternative approach that has been employed successfully in the study of thermal machines is rooted in the framework of quantum mechanical resource theories. Within such resource theories, the allowed operations have to be specified, as well as the 'free resources'. Here, the role of the 'free resources' is assumed by Gibbs states with respect to some Hamiltonians and inverse temperature [28, 29, 31]. The work extraction process is described by a global unitary transformation on the sub-systems prepared in Gibbs states, a system S, as well as a battery. Within such an approach, actual evolution generated by Hamiltonians is not made explicit, and neither is the dynamics leading to equilibration and thermalization. Nevertheless, the allowed resource states are Gibbs states, which are, even if this is not made explicit, of course the result of some equilibration process, possibly involving a larger system. Again, the Gibbs states considered a resource are of the form as in (E.1), with the role of H_S taken over by the Hamiltonians of the sub-systems constituting the resources. In this sense, both approaches are similar in that they crucially rely on Gibbs states of Hamiltonians that are entirely non-interacting with any other part of the system.

However, such an assumption can be a rather implausible one in a number of situations. In fact, this assumption is often excessively restrictive, whenever the thermalizing sub-systems are not entirely decoupled from their environment. Gibbs states have been shown to emerge in small systems interacting very weakly with a large physical body under a number of standard assumptions on the density of states [11]. Such an approach is meaningful in a regime in which

$$\begin{eqnarray}&&\parallel V\parallel \ll {{\beta }^{-1}}.\end{eqnarray} \tag{ E.2 }$$

As $\parallel V\parallel$ is in general extensive, however, and $1/\beta$ is an intensive quantity, such a regime is only meaningful in spin chains or restricted forms of interactions [12]. One can certainly hope for better bounds that also extend to a wider range of physical situations. However, in systems with non-negligible interactions, one would not even expect the above to be a good approximation: one would not expect sub-systems to be well described by Gibbs states with respect to the Hamiltonians of the respective sub-systems. Thermalization then naturally means rather that the reduced states becomes locally indistinguishable from the reduced state of a global Gibbs state (see, e.g., [37–39]). Specifically, if one thinks of a local Hamiltonian H_SB that can for any region of the lattice S and its complement B be decomposed into

$$\begin{eqnarray}&&{{H}_{SB}}={{H}_{S}}+{{H}_{B}}+V,\end{eqnarray} \tag{ E.3 }$$

one would not expect ${{\mathbb{E}}_{t}}({{\rho }_{S}}(t))$ to be close to $\omega ({{H}_{S}})$: certainly the interaction captured by V will alter ${{\mathbb{E}}_{t}}({{\rho }_{S}}(t))$ significantly. In light of these considerations, it seems inadequate to ground the analysis of thermal machines on the existence of resource systems prepared in equilibrium Gibbs states in situations in which interactions can not be considered negligible.

Still, Gibbs states of course play an important role in the description of typical equilibrium reduced states of many-body systems, only it is the Gibbs states of larger systems that have to be taken into account. Consider again a system S and a system B that embodies a large number of degrees of freedom, evolving under the Hamiltonian ${{H}_{SB}}={{H}_{S}}+{{H}_{B}}+V$, where no assumption is made about the strength of the interaction term $V$. For typical local interactions and initial states, and in the absence of local conserved quantities, one expects that

$$\begin{eqnarray}&&{{\mathbb{E}}_{t}}\left\|{{\rho }_{S}}(t)-{\rm t}{{{\rm r}}_{B}}\left( \frac{{{{\rm e}}^{-\beta {{H}_{SB}}}}}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}}})} \right)\right \|{{}_{1}}\ll 1,\end{eqnarray} \tag{ E.4 }$$

where ${{\mathbb{E}}_{t}}$ denotes the expectation in time. This is a consequence of the sub-system being close in trace-norm for most times if the so-called effective dimension is large [12, 14], and the expectation that the time averaged state reduced to S is indistinguishable from ${\rm t}{{{\rm r}}_{B}}(\omega (H))$. That is, again, sub-systems are for most times expected to be operationally indistinguishable from the reduced state of the Gibbs state on a the entire system SB. This is precisely the kind of evolution towards equilibrium on which we base our description of thermal machines.

Assumption 1 Consider a system composed of a sub-system S, a bath B and a battery W. This assumption states that one can place an interaction $V$ between the sub-system and the bath such that the evolution under the Hamiltonian

$$\begin{eqnarray}&&{{H}_{SBW}}={{H}_{S}}+{{H}_{B}}+V+{{H}_{W}}={{H}_{SB}}+{{H}_{W}}\end{eqnarray} \tag{ E.5 }$$

for any initial state ${{\rho }_{SBW}}(t=0)$ and after an appropriately chosen relaxation time τ fulfils

$$\begin{eqnarray}&&{{\rho }_{S}}(t=\tau )={\rm t}{{{\rm r}}_{B}}\left( \frac{{{{\rm e}}^{-\beta {{H}_{SB}}}}}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}}})} \right).\end{eqnarray} \tag{ E.6 }$$

The time $\tau \gt 0$ may be well chosen probabilistically based on a suitable measure, and the statement can be weakened to be true with overwhelming probability. Certainly, one would expect ${{\rho }_{S}}$ to be locally close to the reduction of the time average for the overwhelming proportion of, but not all, times [12, 14]. However, precise error bounds for the equilibration time beyond free models [13] are still an area of active research. For the purposes of this work, therefore, we will take the pragmatic attitude that appropriate times τ can be taken such that the natural condition (E.6) holds true. In the framework of our formalism, this assumption will be taken as a physically plausible assumption and no attempts will be made to derive bounds to equilibration times.

Treating the thermalization map (7) as the result of an actual time evolution compels one to also apply a time evolution to the battery. As we discuss in section appendix H this will result in a loss of the coherence of the battery, which renders it in general impossible to perform further Hamiltonian quenches on SB. However, in realistic situations, the thermal machine SBW can be assumed to be weakly interacting with the surrounding environment. This will effectively produce decoherence—that is, damping the off diagonal terms in the Hamiltonian eigenbasis [14]. As there is no interaction between SB and W, however, both are weakly interacting with a local environment, decoherence is is expected to be most effective on the product eigenbasis of ${{H}_{SB}}+{{H}_{W}}$. This effect, as we show in (D.4) allows one to perform further quenches without the need for coherence.

Assumption 2 Consider a system composed of a system SB and a battery W, equipped with a non-interacting Hamiltonian

$$\begin{eqnarray}&&{{H}_{SBW}}={{H}_{SB}}+{{H}_{W}}.\end{eqnarray} \tag{ E.7 }$$

This assumption states that the evolution induced by the interaction of SBW with a suitable natural environment E is equivalent with the application of a decoherence map $\mathcal{E}$ described by

$$\begin{eqnarray}&&\mathcal{E}({{\rho }_{SBW}})=\mathop{\sum }\limits_{i,w}({{\rho }_{SBW}})_{i,i}^{w,w}|i\rangle \langle i|\otimes |w\rangle \langle w|,\end{eqnarray} \tag{ E.8 }$$

where ${{\rho }_{SBW}}={{\sum }_{i,i^{\prime} ,w,w^{\prime} }}({{\rho }_{SBW}})_{i,i^{\prime} }^{w,w^{\prime} }|i\rangle \langle i^{\prime} |\otimes |w\rangle \langle w^{\prime} |$, ${{H}_{SB}}={{\sum }_{i}}{{E}_{i}}|i\rangle \langle i|$ and ${{H}_{W}}={{\sum }_{w}}w|w\rangle \langle w|$.

We will now turn to analysing the formalism of quenches described in section ?? when the change of the Hamiltonian is implemented on a sub-system S in contact with a thermal bath B. Consider an initial global state $\rho _{SBW}^{(0)}(t=0)$ and an initial Hamiltonian for the thermal machine $H_{S}^{(0)}+{{H}_{B}}+V+{{H}_{W}}$. We then allow this system to equilibrate according to this Hamiltonian, so that the evolution fulfils assumptions 1 and 2. Hence, at a large enough time τ the state can be written as,

$$\begin{eqnarray}&&\rho _{SBW}^{(0)}(\tau )=\mathop{\sum }\limits_{i,w}(\rho _{SBW}^{(0)}(0))_{i}^{w}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\otimes |w\rangle \langle w|\end{eqnarray} \tag{ F.1 }$$

where $(\rho _{SBW}^{(0)}(0))_{i}^{w}=\langle {{i}^{(0)}}|\otimes \langle w|\rho _{SBW}^{(0)}(0)|{{i}^{(0)}}\rangle \otimes |w\rangle$ and

$$\begin{eqnarray}&&H_{SB}^{(0)}\;:=\;H_{S}^{(0)}+{{H}_{B}}+V=\mathop{\sum }\limits_{i}E_{i}^{(0)}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|.\end{eqnarray} \tag{ F.2 }$$

Also, the equilibrated state fulfils

$$\begin{eqnarray}&&\rho _{S}^{(0)}(\tau )={\rm t}{{{\rm r}}_{B}}\left( \frac{{{{\rm e}}^{-\beta {{H}_{SB}}}}}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}}})} \right).\end{eqnarray} \tag{ F.3 }$$

At time τ we perform a quench $H_{S}^{(0)}+{{H}_{B}}+V\mapsto H_{S}^{(1)}+{{H}_{B}}+V={{\sum }_{i}}E_{i}^{(1)}|{{i}^{(1)}}\rangle \langle {{i}^{(1)}}|$. The state after the quench $\rho _{SBW}^{(1)}(t=T)$ satisfies

$$\begin{eqnarray}&&\rho _{SBW}^{(1)}(\tau )=\mathop{\sum }\limits_{i,w}(\rho _{SBW}^{(0)}(0))_{i}^{w}U_{RW}^{{\rm on}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\otimes |w\rangle \langle w|U_{RW}^{{\rm on}\dagger },\end{eqnarray} \tag{ F.4 }$$

where ${{U}^{{\rm on}}}$ is the quench unitary as defined in (22). Hence, the work extracted at the battery is

$$\begin{eqnarray*}\begin{array}{rcl} \langle W\rangle & = & {\rm t}{{{\rm r}}_{W}}({{H}_{W}}\rho _{W}^{(1)}(\tau ))-{\rm t}{{{\rm r}}_{W}}({{H}_{W}}\rho _{W}^{(0)}(\tau )) \\ {} & = & \mathop{\sum }\limits_{i,w}(\rho _{SBW}^{(0)}(0))_{i}^{w}{\rm t}{{{\rm r}}_{W}}({{H}_{W}}{\rm t}{{{\rm r}}_{SB}}(U_{RW}^{{\rm on}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\otimes |w\rangle \langle w|U_{RW}^{{\rm on}\dagger })) \\ {} & {} & -\mathop{\sum }\limits_{i,w}(\rho _{SBW}^{(0)}(0))_{i}^{w}{\rm t}{{{\rm r}}_{W}}({{H}_{W}}{\rm t}{{{\rm r}}_{SB}}(|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\otimes |w\rangle \langle w|) \\ \end{array}\end{eqnarray*}$$

$$\begin{eqnarray}\begin{array}{cll} {} & = & {{\sum }_{i,w}}(\rho _{SBW}^{(0)}(0))_{i}^{w}({\rm t}{{{\rm r}}_{SB}}((H_{SB}^{(0)}-H_{SB}^{(1)})|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|)+w) \\ {} & {} & -{{\sum }_{i,w}}(\rho _{SBW}^{(0)}(0))_{i}^{w}w \\ \end{array}\end{eqnarray} \tag{ F.5 }$$

$$\begin{eqnarray}\begin{array}{clc} {} & = & {\rm t}{{{\rm r}}_{SB}}\left( (H_{SB}^{(0)}-H_{SB}^{(1)})\rho _{SB}^{(0)}(\tau ) \right) \\ {} & = & {\rm t}{{{\rm r}}_{SB}}\left( \left( (H_{S}^{(0)}-H_{S}^{(1)})\otimes {{\mathbb{I}}_{B}} \right)\frac{{{{\rm e}}^{-\beta {{H}_{SB}}}}}{{\rm tr}({{{\rm e}}^{-\beta {{H}_{SB}}}})} \right), \\ \end{array}\end{eqnarray} \tag{ F.6 }$$

where (F.5) follows from (D.4) and (F.6) from (F.3).

We now combine the statements of (C.10) and (F.6) in order to show that the work extraction protocol as defined in section 2 can be implemented.

(Physical implementation in a unitary framework).

Corollary 1 Given an initial state of the form ${{\rho }^{(0)}}=\rho _{SB}^{(0)}\otimes |{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}{{|}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}}$, with $|{{\Psi }_{\sqcap }}{{\rangle }_{W}}$ as defined in (C.2), and an arbitrary initial Hamiltonian ${{H}^{(0)}}$, assuming the validity of assumptions 1 and 2, any protocol $\mathcal{P}$ can be implemented with a unitary transformation acting on a system composed of the thermal machine SBW, the control qubit Q and an environment E.

Proof. This statement follows straightforwardly from assumptions 1 and 2, and (C.10) and (F.6).Given the initial state ${{\rho }^{(0)}}=\rho _{SB}^{(0)}\otimes |{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}{{|}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}}$, (C.10) shows that the quench unitary (22) performs the first Hamiltonian transformation of an arbitrary protocol $\mathcal{P}$—before the first state thermalization—so that it fulfils condition (4). Then, the unitary evolution under the composed system SBWQE satisfying assumptions 1 and 2 results in further quenches being fulfilled (F.6), which in turn implies that it fulfils (4) when applied on thermalized states as in (7).□

As we have discussed in appendix C, a coherent state of the battery allows one to perform a Hamiltonian quench. This can be easily seen from (C.6), if one applies a quench to an initial state of the form

$$\begin{eqnarray}&&\rho _{RWQ}^{(0)}=\rho _{R}^{(0)}\otimes |{{\Psi }_{\sqcap }}\rangle \langle {{\Psi }_{\sqcap }}{{|}_{W}}\otimes |0\rangle \langle 0{{|}_{Q}}\end{eqnarray} \tag{ H.1 }$$

—R plays the role of system plus bath—the reduced final state on R does not change, that is

$$\begin{eqnarray}&&\rho _{R}^{(1)}=\rho _{R}^{(0)}.\end{eqnarray} \tag{ H.2 }$$

Let us suppose that now we let the system RW undergo a time evolution under the Hamiltonian $H_{R}^{(1)}+{{H}_{W}}$—this is precisely what one does if R embodies both a system S and a bath, and the time evolution is intended to drive $\rho _{R}^{(1)}$ towards a thermalized state of the form (7). How does this time evolution affect the coherence in the state of the battery? Is the battery still coherent so that it can perform further quenches? Here we show that this is not the case. Coherence is a resource that is lost under such a time evolution. To see this, let us compute the time-evolved state after time t of $\rho _{R}^{(1)}$ which is given by

$$\begin{eqnarray*}\begin{array}{rcl} \rho _{RW}^{(1)}(t) & = & \mathop{\sum }\limits_{i,i^{\prime} ,j,j^{\prime} }\langle {{j}^{(1)}}|{{i}^{(0)}}\rangle \langle {{i}^{(0)}}|\rho _{R}^{(0)}|i{{^{\prime} }^{(0)}}\rangle \langle i{{^{\prime} }^{(0)}}|j{{^{\prime} }^{(1)}}\rangle {{{\rm e}}^{-{\rm i}(E_{j}^{(1)}-E_{j^{\prime} }^{(1)})t}}|{{j}^{(1)}}\rangle \langle j{{^{\prime} }^{(1)}}{{|}_{R}} \\ {} & {} & \otimes \frac{1}{N(E_{W}^{(0)},\Delta )}\mathop{\sum }\limits_{E_{W}^{(0)}+\Delta }^{w=E_{W}^{(0)}}{{{\rm e}}^{-{\rm i}({{\Delta }_{i,j}}+w)t}}|{{\Delta }_{i,j}}+w\rangle \mathop{\sum }\limits_{E_{W}^{(0)}+\Delta }^{w^{\prime} =E_{W}^{(0)}}{{{\rm e}}^{{\rm i}({{\Delta }_{i^{\prime} ,j^{\prime} }}+w)t}}\langle {{\Delta }_{i^{\prime} ,j^{\prime} }}+w^{\prime} |. \\ \end{array}\end{eqnarray*}$$

From this equation one can straightforwardly, but tediously, conclude that

$$\begin{eqnarray}&&\rho _{R}^{(1)}(t)={{{\rm e}}^{-{\rm i}H_{R}^{(1)}t}}\rho _{R}^{(0)}{{{\rm e}}^{{\rm i}H_{R}^{(1)}t}},\end{eqnarray} \tag{ H.3 }$$

that is, as one would expect, the initial state evolved under $H_{R}^{(1)}$ at time t. Now, if one intends to perform further quenches on this state—that is, a unitary of the form (22) changing $H_{R}^{(1)}$ to $H_{R}^{(2)}$ without altering the state on R—one finds that this is not possible, because the state of the battery has been changed by the evolution under H_W and it no longer serves as a coherent resource fulfilling (C.1). This can be shown by a tedious calculation applying the unitary (22) on (H.3). To avoid such a calculation and to merely grasp the intuition behind the mechanism, note that the state vector

$$\begin{eqnarray}&&|{{\Psi }_{\sqcap }}(t)\rangle =\frac{1}{N{{(E_{W}^{(0)},\Delta )}^{1/2}}}\sum W(0)E_{W}^{(0)}+\Delta w=E{{{\rm e}}^{-{\rm i}({{\Delta }_{i,j}}+w)t}}|{{\Delta }_{i,j}}+w\rangle \end{eqnarray} \tag{ H.4 }$$

no longer fulfils (C.1) when a new quench from $H_{R}^{(1)}$ to $H_{R}^{(2)}$—with energy gaps $\Delta _{i,j}^{(2)}$—is applied. Indeed, it is easy to see that for most times t

$$\begin{eqnarray}&&\langle {{\Psi }_{\sqcap }}(t)|\Gamma (\Delta _{i,j}^{(2)}-\Delta _{i^{\prime} ,j^{\prime} }^{(2)})|{{\Psi }_{\sqcap }}(t)\rangle \approx 0.\end{eqnarray} \tag{ H.5 }$$

In other words, the coherence of the battery is lost due to time evolution under its own Hamiltonian H_W, and this is an obstacle against performing further quenches in general. In the specific protocol leading to corollary 1, further quenches can be applied because coherence is no longer needed after the decoherence map specified in assumption 2 has been applied. We expect this decoherence map to reasonably represent plausible and realistic physical situations. However, it should be clear that alternative protocols in which, for instance, coherence is re-established in the battery by a certain operation—possibly employing a device playing the role of a source of coherence—are also of great interest. As a matter of fact, the role of coherence and how it should be accounted for as a resource in thermodynamics is an interesting open question which we leave open for future work. Note that the role played by coherence in this work is quite different from that in [31]. There, coherence is a catalytic resource, in the sense that it is not consumed in the protocol and can be re-used an arbitrary number of times. Our analysis points out that such a catalytic role of the coherence may be only an artefact of the specific framework of operations considered there, where time evolutions are not taken into account.

As far as work extraction is concerned, in our work, we follow the approach of, e.g., [30] and consider average work extraction. Our results hence apply to the expected work for individual systems. Note that we do not have to assume at any point—similarly to [29]—that we process N copies collectively in order to obtain (9). Due to the linearity of the work extraction process, it is implied by a basic argument of typicality that when processing N copies, the total work extracted per copy will be essentially deterministic in the limit of large N—the variance increases with $\sqrt{N}$ and the total work with N. However, it is still of interest to analyse the spread of the probability distribution of the energy in the battery for a single copy. This is relevant with generalizations to single-shot work extraction in the spirit of [23, 25, 28]. Note that such analysis is out of place within the abstract formalism defined in section 2: the operations just preserve the expected energy, thus transformations reducing arbitrarily the spread of the energy of the battery are allowed, similarly as in the formalism defined in [30]. Nonetheless, note that the unitary implementation of the protocol of corollary 2 does preserve the probability distribution of the entire machine SBW. This is the case because (i) the unitary defined in theorem 2 does not only preserve the mean total energy, but it also commutes with the total Hamiltonian and (ii) the dephasing map employed when the system relaxes to an equilibrium state, as defined in assumption 2, by definition preserves the probability distribution of energies of SBW. Therefore, one could restrict the set of operations defined in section 2 by substituting the assumption of mean energy conservation for a conservation of the probability distribution of total energy, and a protocol saturating (9) would still be attainable. In conclusion, the formalism itself, in contrast to the one in [30], can be easily modified to account for a possible generalization in therms of single-shot work extraction.

Nevertheless, there is another issue that prevents one from applying the findings of [23, 25, 28] in a straightforward manner. This is the impossibility of performing quenches with deterministic classical states of the battery. As detailed in section appendix C, one needs to employ initial state vectors of the battery $|{{\psi }_{\sqcap }}{{\rangle }_{W}}$. Therefore, the initial probability distribution of energies of the battery is already 'infinitely spread out'. As discussed in [31], a distinction between ordered work—as in the single-shot work extraction—and disordered work would need to take into account the energy carrier—in this case the battery—and how the initially spread distribution of the battery is affected by the protocol. We leave this as an interesting open question that lies outside the scope of this work.

We now turn to the discussion of the typicality of irreversibility and the relationship to an instance of the second law. The equivalence between optimality and reversibility in work extraction protocols is widely known in the context of phenomenological thermodynamics, the analysis of the Carnot engine being the most seminal example. More generally, Clausius' theorem states that overall heat flow vanishes over all reversible cyclic processes. That is,

$$\begin{eqnarray}&&{{\oint }_{{\rm rev}}}\frac{\delta Q}{T}=0\end{eqnarray} \tag{ J.1 }$$

where $\delta Q$ is the inexact differential of the heat Q and T is the temperature. This motivates the definition of the entropy state function as ${\rm d}S\;:=\;\delta Q/T$, T taking the role of the integrating factor. Furthermore, Clausius' inequality establishes that for general processes—not necessarily reversible or cyclic—it is true that

$$\begin{eqnarray}&&\frac{\Delta Q}{T}=\int _{i}^{f}\frac{\delta Q}{T}\leqslant \Delta S={{S}_{f}}-{{S}_{i}},\end{eqnarray} \tag{ J.2 }$$

where equality holds in the reversible case.

This theorem is formulated within the framework of phenomenological thermodynamics. However, similar expressions can be shown to hold within a statistical mindset with the von Neumann entropy taking over the role of thermodynamic entropy [26]. Indeed, in the weak-coupling setting, it is not difficult to show that (J.2) is indeed equivalent to the bounds on expected work extraction in terms of the free-energy difference, and also that optimal work extraction processes are reversible.

To see this, consider a protocol of work extraction by Hamiltonian quenches as defined in section 2. In the weak-coupling regime, the state thermalization map (7) is replaced by $\rho _{S}^{(i)}\mapsto \omega (H_{S}^{(i)})$. Equivalent to (12), the expected work extracted in a general protocol in the weak-coupling limit is given by

$$\begin{eqnarray}&&\langle {{W}^{{\rm wc}}}\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)})={\rm tr}(\rho _{S}^{(0)}(H_{S}^{(0)}-H_{S}^{(1)}))+\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}(\omega (H_{S}^{(i)})(H_{S}^{(i)}-H_{S}^{(i+1)})),\end{eqnarray} \tag{ J.3 }$$

which, recalling (B.2), fulfils

$$\begin{eqnarray}&&\mathop{{\rm max} }\limits_{\mathcal{P}}\langle {{W}^{{\rm wc}}}\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)})\leqslant F(\rho _{S}^{(0)},H_{S}^{(0)})-F(\omega (H_{S}^{(0)}),H_{S}^{(0)}).\end{eqnarray} \tag{ J.4 }$$

Equality is achieved here by a reversible process. Now let us see that a similar conclusion can be reached from (J.2). If we define the heat flow $\Delta Q$ as the energy lost by the bath—or equivalently, the energy gained by the system in the state thermalization process—one can see that

$$\begin{eqnarray}\begin{array}{rcl} \langle \Delta {{Q}^{{\rm wc}}}\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)}) & = & {\rm tr}((\omega (H_{S}^{(1)})-\rho _{S}^{(0)})H_{S}^{(1)}) \\ {} & {} & +\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}((\omega (H_{S}^{(i+1)})-\omega (H_{S}^{(i)}))H_{S}^{(i+1)})) \\ \end{array}\end{eqnarray} \tag{ J.5 }$$

$$\begin{eqnarray}&&=-{\rm tr}(\rho _{S}^{(0)}H_{S}^{(1)})+\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}(\omega (H_{S}^{(i)})(H_{S}^{(i)}-H_{S}^{(i+1)}))+{\rm tr}(\omega (H_{S}^{(n)})H_{S}^{(n)})\end{eqnarray} \tag{ J.6 }$$

$$\begin{eqnarray}&&=\;\langle {{W}^{{\rm wc}}}\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)})+{{\langle \Delta E\rangle }_{S}}\end{eqnarray} \tag{ J.7 }$$

where

$$\begin{eqnarray}&&{{\langle \Delta E\rangle }_{S}}\;:=\;{\rm tr}(\omega (H_{S}^{(n)})H_{S}^{(n)})-{\rm tr}(\rho _{S}^{(0)}H_{S}^{(0)})\end{eqnarray} \tag{ J.8 }$$

is the expected energy difference between the initial and final states. Therefore, identifying $\langle \Delta {{Q}^{{\rm wc}}}\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)})$ with the heat flow, and the von Neumann entropy with the thermodynamic entropy in (J.2), one obtains

$$\begin{eqnarray}\begin{array}{rcl} \langle {{W}^{{\rm wc}}}\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)}) & \leqslant & T\Delta {{S}_{S}}-{{\langle \Delta E\rangle }_{S}} \\ {} & = & F(\rho _{S}^{(0)},H_{S}^{(0)})-F(\omega (H_{S}^{(0)}),H_{S}^{(0)}), \\ \end{array}\end{eqnarray} \tag{ J.9 }$$

where, according to Clausius' theorem, equality again holds when the process is reversible. This equivalence between Clausius' theorem and the work extraction bounds means that indeed (J.4) may be understood as an alternative formulation of the second law of thermodynamics applied to expectation values. Also, the fact that there exists an optimal reversible protocol saturating (J.4) is to be understood as saturation of the second law.

Let us now investigate the situation where the interaction between the bath and system is not necessarily weak and the thermalization map is of the form (7). As anticipated in section 2, in general the coupling between bath and systems prevents one from saturating the second law in the form stated above and to perform reversible processes.

The first difference when analysing the strong-coupling case is that the very definition of heat is problematic. In a system evolving from $\rho _{SB}^{(i)}$ to $\rho _{SB}^{(i)}$ equipped with the Hamiltonian ${{H}_{S}}+{{H}_{B}}+{{V}_{SB}}$, it is not quite clear how much energy is lost by the bath—this is the canonical definition of heat—because the energy contribution of the interaction is not negligible, and it is not obvious which part corresponds to the bath and to the system. To find a way of circumventing this problem, let us consider an specific example. Let us partition the bath B into two regions ${{B}_{{\rm b}}}$ (the buffer) and ${{B}_{{\rm r}}}$ (the reservoir). The buffer represents the region of the bath that is surrounding the system S and the reservoir is the region that is not directly in contact with S. Let us suppose that ${{B}_{{\rm b}}}$ and ${{B}_{{\rm r}}}$ are weakly coupled, so that the operator norm ${{V}_{{{B}_{{\rm b}}},{{B}_{{\rm r}}}}}$ is much smaller than the energy gaps of their respective Hamiltonians. This would be the case if, for instance, S and ${{B}_{{\rm b}}}$ are parts of a conducting material, and ${{B}_{{\rm r}}}$ is just a surrounding gas that interacts weakly with ${{B}_{{\rm b}}}$. For such a set-up the equilibration towards the equilibrium of S will fulfil,

$$\begin{eqnarray}&&{{\rho }_{S,{{B}_{{\rm b}}}}}=\omega ({{H}_{S{{B}_{{\rm b}}}}}),\end{eqnarray} \tag{ J.10 }$$

where ${{H}_{S{{B}_{{\rm b}}}}}={{H}_{S}}+{{H}_{{{B}_{{\rm b}}}}}+{{V}_{S{{B}_{{\rm b}}}}}$, and ${{V}_{S{{B}_{{\rm b}}}}}$ is an arbitrarily strong interaction that only has support in ${{B}_{{\rm b}}}$ (but not in ${{B}_{{\rm r}}}$). In this case, the weak interaction between ${{B}_{{\rm b}}}$ and ${{B}_{{\rm r}}}$ establishes a clear cut that allows one to unambiguously define the energy that was lost by the the reservoir ${{B}_{{\rm r}}}$—in contrast to the energy that has flowed from ${{B}_{{\rm b}}}$ to S that is ambiguous due to the strong coupling in ${{V}_{{{B}_{{\rm b}}},S}}$. Hence, the definition of heat can be made unambiguous as the energy lost by the reservoir ${{B}_{{\rm r}}}$, or equivalently, the energy gained by $S{{B}_{{\rm b}}}$. Taking this as the definition of heat, one obtains

$$\begin{eqnarray}\begin{array}{rcl} \langle \Delta Q\rangle (\mathcal{P},H_{S{{B}_{{\rm b}}}}^{(0)},\rho _{S{{B}_{{\rm b}}}}^{(0)}) & = & {\rm tr}\left( \left( \omega (H_{S{{B}_{{\rm b}}}}^{(1)})-\rho _{S{{B}_{{\rm b}}}}^{(0)} \right)H_{S{{B}_{{\rm b}}}}^{(1)} \right) \\ {} & {} & +\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}\left( \left( \omega (H_{S{{B}_{{\rm b}}}}^{(i+1)})-\omega (H_{S{{B}_{{\rm b}}}}^{(i)}) \right)H_{S{{B}_{{\rm b}}}}^{(i+1)} \right) \\ \end{array}\end{eqnarray} \tag{ J.11 }$$

$$\begin{eqnarray}\begin{array}{clc} {} & = & -{\rm tr}\left( \rho _{S{{B}_{{\rm b}}}}^{(0)}H_{S{{B}_{{\rm b}}}}^{(1)} \right)+\mathop{\sum }\limits_{i=1}^{n-1}{\rm tr}\left( \omega (H_{S{{B}_{{\rm b}}}}^{(i)})(H_{S{{B}_{{\rm b}}}}^{(i)}-H_{S{{B}_{{\rm b}}}}^{(i+1)}) \right)+{\rm tr}\left( \omega (H_{S{{B}_{{\rm b}}}}^{(n)})H_{S{{B}_{{\rm b}}}}^{(n)} \right) \\ {} & = & \langle W\rangle (\mathcal{P},H_{S{{B}_{{\rm b}}}}^{(0)},\rho _{S{{B}_{{\rm b}}}}^{(0)})+{{\langle \Delta E\rangle }_{S{{B}_{{\rm b}}}}}, \\ \end{array}\end{eqnarray} \tag{ J.12 }$$

where

$$\begin{eqnarray}&&{{\langle \Delta E\rangle }_{S{{B}_{{\rm b}}}}}\;:=\;{\rm tr}\left( \omega (H_{S{{B}_{{\rm b}}}}^{(n)})H_{S{{B}_{{\rm b}}}}^{(n)} \right)-{\rm tr}\left( \rho _{S{{B}_{{\rm b}}}}^{(0)}H_{S{{B}_{{\rm b}}}}^{(0)} \right).\end{eqnarray} \tag{ J.13 }$$

Using (J.2) and identifying $\Delta S=S(\omega (H_{S{{B}_{{\rm b}}}}^{(n)}))-S(\rho _{S{{B}_{{\rm b}}}}^{(0)})$, one obtains

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle (\mathcal{P},H_{S{{B}_{{\rm b}}}}^{(0)},\rho _{S{{B}_{{\rm b}}}}^{(0)}) & \leqslant & T\Delta {{S}_{S{{B}_{{\rm b}}}}}-{{\langle \Delta E\rangle }_{S{{B}_{{\rm b}}}}} \\ {} & = & F\left( \rho _{S{{B}_{{\rm b}}}}^{(0)},H_{S{{B}_{{\rm b}}}}^{(0)} \right)-F\left( \omega (H_{S{{B}_{{\rm b}}}}^{(0)}),H_{S{{B}_{{\rm b}}}}^{(0)} \right) \\ {} & = & -\Delta {{F}_{{\rm rev}}}, \\ \end{array}\end{eqnarray} \tag{ J.14 }$$

where (J.14) is a consequence of theorem 1 and taking ${{\tilde{\rho }}_{S{{B}_{{\rm b}}}}}=\rho _{S{{B}_{{\rm b}}}}^{(0)}$.

Lastly, in the case of the thermalization map of the form (7), where no assumption is made about a cut between the buffer and the reservoir, the entire bath has to be considered as being the buffer ${{B}_{{\rm b}}}$ and the reservoir is not present. Then, in analogy to (J.10), if we strengthen condition (7) by assuming that the equilibrium state fulfils

$$\begin{eqnarray}&&{{\rho }_{SB}}=\omega ({{H}_{SB}}),\end{eqnarray} \tag{ J.15 }$$

one can define heat unambiguously as the energy gained by the whole machine—which vanishes by an argument based on the conservation of energy. Indeed, we find

$$\begin{eqnarray}&&\langle \Delta Q\rangle (\mathcal{P},H_{SB}^{(0)},\rho _{SB}^{(0)})=\langle W\rangle (\mathcal{P},H_{S}^{(0)},\rho _{S}^{(0)})+{{\langle \Delta E\rangle }_{SB}}=0,\end{eqnarray} \tag{ J.16 }$$

where the final equality follows simply from the expected energy conservation. Therefore, in a scenario based on a a thermalization map of the kind considered in(7), the second law can be written simply as

$$\begin{eqnarray}&&0\leqslant \Delta {{S}_{SB}}\end{eqnarray} \tag{ J.17 }$$

where equality is fulfilled by a reversible process. This together with (J.16) gives again

$$\begin{eqnarray}\begin{array}{rcl} \langle W\rangle (\mathcal{P},H_{SB}^{(0)},\rho _{S{{B}_{{\rm b}}}}^{(0)}) & \leqslant & T\Delta {{S}_{SB}}-{{\langle \Delta E\rangle }_{SB}} \\ {} & = & F\left( \rho _{SB}^{(0)},H_{SB}^{(0)} \right)-F\left( \omega (H_{SB}^{(0)}),H_{SB}^{(0)} \right) \\ {} & = & -\Delta {{F}_{{\rm rev}}}, \\ \end{array}\end{eqnarray} \tag{ J.18 }$$

where the equality is satisfied by reversible protocols of work extraction. The bound of theorem 1 thus usually imposes a limitation, quantifiable by $\Delta {{F}_{{\rm irrev}}}$, against saturating the second law of thermodynamics (J.17). The reason, as the very formulation of the second law by Clausius' theorem already takes into account, is that the process is not reversible. This can be easily seen from (A.10, A.11). The optimal protocol specifies a Hamiltonian $H_{SB}^{(1)}$, and a parametrized curve of Hamiltonians describing a trajectory from $H_{SB}^{(1)}$ to $H_{SB}^{(0)}$. Now one can reverse the protocol, that is, given $\mathcal{P}$ by $\{H_{SB}^{(i)}\}_{i=1}^{n-1}$ and ${\bf k}$, we define the inverse protocol ${{\mathcal{P}}^{-1}}$ by $\{H_{SB}^{(n-i)}\}_{i=1}^{n-1}$ and ${{{\bf k}}^{-1}}(i)\;:=\;{\bf k}(n-i)$, and a simple calculation shows

$$\begin{eqnarray}&&\begin{array}{l} \langle W\rangle ({{\mathcal{P}}^{(-1)}},H_{SB}^{0},{{\omega }_{S}}(H_{SB}^{(0)})) \\ \quad ={\rm tr}\left( \omega (H_{SB}^{(1)})(H_{S}^{(1)}-H_{S}^{(0)}) \right)+F\left( \omega (H_{SB}^{(0)}),H_{SB}^{(0)} \right)-F\left( \omega (H_{SB}^{(1)}),H_{SB}^{(1)} \right) \\ \quad =F\left( \omega (H_{SB}^{(0)}),H_{SB}^{(0)} \right)-F\left( \omega (H_{SB}^{(1)}),H_{SB}^{(0)} \right) \\ \quad =\langle W\rangle (\mathcal{P},H_{SB}^{0},\rho _{SB}^{(0)}))+\Delta {{F}_{{\rm irrev}}}. \\ \end{array}\end{eqnarray} \tag{ J.19 }$$

That is, the work difference between the optimal protocol and its reversed protocol is precisely $\Delta {{F}_{{\rm irrev}}}$. This quantity is exactly the amount by which the work extraction bounds differ from the maximum ones allowed by the second law stated in the form $0\leqslant \Delta {{S}_{SB}}$. Altogether, this suggests that theorem 1 may be viewed as a generalization of the second law of thermodynamics which accounts for strong couplings and the unavoidable irreversibility that it induces.

The irreversibility of the optimal process may result in a tension with theorem 1, where it is shown that a global unitary evolution performs the optimal protocol, and therefore it must be reversible. This apparent paradox is resolved by noting that being reversible at the level of abstract protocols—that is, as we define ${{\mathcal{P}}^{-1}}$—is not equivalent to being reversible in the sense of time-reversed implementation. Note that the time-reversed evolution can take equilibrium states to states out of equilibrium, however, a reversed protocol in the sense of ${{\mathcal{P}}^{-1}}$ does not allow for such passages from equilibrium to non-equilibrium states. This is precisely the case, for example, in the first step of the optimal protocol detailed in the proof of theorem 1. There, the initial Hamiltonian $H_{SB}^{(0)}$ is quenched to $H_{SB}^{(1)}$, and then the state of SB is driven to equilibrium, so that $\rho _{S}^{(0)}\mapsto {{\omega }_{S}}(H_{SB}^{(1)})$. Clearly, this equilibration is eventually due to some unitary evolution of the composed system SB, and indeed could in principle be reversed if one had control over the exact time that we waited until

$$\begin{eqnarray}&&\rho _{S}^{(0)}\to {{\omega }_{S}}(H_{SB}^{(1)})\end{eqnarray} \tag{ J.20 }$$

has converged. However, at the abstract level mainly considered here, where work extraction protocols $\mathcal{P}$ are being defined, the protocols neither explicitly take time into account nor any other dynamical analysis of the state thermalization. Therefore, a reversed protocol of the previous example would just amount to a quench from $H_{SB}^{(1)}$ to $H_{SB}^{(0)}$ on the state ${{\omega }_{S}}(H_{SB}^{(1)})$. The use of the abstract map (7) is grounded precisely in typicality arguments, as explained in appendix E. In other words, the irreversibility exhibited by the optimal protocols should also be understood as a feature of typicality: given the precise times that one has waited in each equilibration process, $t=({{\tau }_{1}},{{\tau }_{2}},...,{{\tau }_{l}})$, for most times, with overwhelmingly high probability, the optimal protocol extracts $-(\Delta {{F}_{{\rm rev}}}-\Delta {{F}_{{\rm irrev}}})$. If one applies the reversed protocol, with suitable times for equilibration, to most times and all initial states, with overwhelming probability, the work extracted in the inverse protocol would be $-\Delta {{F}_{{\rm rev}}}$. Therefore, the optimal protocol is typically irreversible.