Quantum fluctuation theorems for arbitrary environments: adiabatic and non-adiabatic entropy production

We analyze the production of entropy along non-equilibrium processes in quantum systems coupled to generic environments. First, we show that the entropy production due to final measurements and the loss of correlations obeys a fluctuation theorem in detailed and integral forms. Second, we discuss the decomposition of the entropy production into two positive contributions, adiabatic and non-adiabatic, based on the existence of invariant states of the local dynamics. Fluctuation theorems for both contributions hold only for evolutions verifying a specific condition of quantum origin. We illustrate our results with three relevant examples of quantum thermodynamic processes far from equilibrium.


I. INTRODUCTION
Classical thermodynamics and statistical mechanics provide a systematic approach to the phenomenology of a system immersed in a large environment. Within these frameworks, two complementary strategies are employed. The first is to explicitly model the environment -often an equilibrium thermal reservoir-to obtain an effective reduced dynamics for the system alone, which then can be analyzed. The second is to derive fundamental constraints in the form of inequalities using the second law of thermodynamics and magnitudes like entropy, entropy production and free energy. The recent introduction of an entropy for stochastic trajectories [1] allows one to refine these inequalities with exact equalities for arbitrary nonequilibrium processes, results generically known as fluctuation theorems (FT's) [2,3].
These two strategies have also been successfully applied to quantum systems. Open quantum system dynamics -the determination and analysis of the system's reduced dynamics-is a well-developed and active field [4,5]. Complementing this approach, a variety of quantum FT's have been derived [6][7][8][9][10][11] to asses the statistics of the relevant quantities. Different proposals to obtain these statistics in the laboratory have been reported, using techniques related to quantum tomography [12][13][14][15][16][17], and some of them have already been used to carry out experimental verifications of FT's [18,19]. However, most of the research on quantum FT's is only valid for equilibrium reservoirs with a focus on the energy exchange between the system and the environment in the form of heat and work. By contrast, classical FT's have * gmanzano@ucm.es been formulated more generally for generic Markov systems [20][21][22][23][24] using the entropy production instead of heat and work, which are only meaningful in physical situations where a system exchanges energy with equilibrium reservoirs.
The task of deriving FT's for generic quantum dynamics also implies a more detailed characterization of entropy production in nonequilibrium quantum contexts, a problem that has attracted a growing interest in recent years [8,30,31,33,[50][51][52][53][54][55][56][57]. In Ref. [33] we derived a FT for a class of completely-positive trace-preserving (CPTP) quantum maps, which model a variety of quantum processes. Through this analysis, we identified a quantity that coincides with the entropy production for thermalization processes and resembles the nonadiabatic entropy production introduced in the classical context [21][22][23]. The purpose of this paper is to clarify and extend those results considering together the system and its surroundings. By tracing over the environment, we can then recover the quantum map for the reduced system dynamics. This setup allows us to unveil the origin of entropy production in arbitrary processes from coarse-graining, and derive corresponding FT's. We also split the entropy production into an adiabatic and a nonadiabatic arXiv:1710.00054v2 [quant-ph] 1 Jul 2018 contribution, exactly as in classical stochastic thermodynamics. However, contrary to what happens in classical systems, the split is not always possible. A condition, derived in Ref. [33], is necessary. We will explore the similarities and differences between classical and quantum FT's in concrete examples.
The paper is organized as follows. In Section II, we introduce a thermodynamic process for a generic bipartite system that models a system and its environment. We will define in this section the entropy production along the process and the concomitant reduced system dynamics. We develop a FT for this entropy production in Section III using a time-reversed or backward thermodynamic process. In Section IV, FT's for the adiabatic and nonadiabatic entropy production are derived. Our results are also extended both to the case of concatenations of CPTP maps and multipartite environments. This is applied to the specific case of quantum trajectories unraveled from Lindblad Master Equations in V. Finally, relevant examples to illustrate our results are given in Section VI while concluding in Section VII with some final remarks.

II. QUANTUM OPERATIONS AND ENTROPY PRODUCTION
Along the paper we consider an isolated quantum system composed of two parts, system and environment (or ancilla), with Hilbert space H ≡ H S ⊗ H E , where H S and H E are the local Hilbert spaces of the system and the environment respectively. We focus our attention on the entropy production along the generic process depicted in Fig. II, consisting of initial and final local projective measurements that bracket a unitary evolution. The outcomes of the measurements constitute a quantum trajectory, which plays a crucial role in the formulation of FT's.

A. The (forward) process
The process begins with the global system in an uncorrelated product state ρ SE = ρ S ⊗ ρ E . The spectral decomposition of the local states reads where p n and q ν are the eigenvalues, and {P n } and {Q ν } are orthogonal projectors onto their respective eigensubpaces (for simplicity we assume they are rank-1 projectors). At t = 0 an initial projective measurement on the system and environment is performed using the eigenprojectors in Eq. (1) and outcomes n and ν are obtained. This measurement projects the system and environment onto pure states |ψ n ψ n | S ≡ P n and |φ ν φ ν | E ≡ Q ν , without modifying the average or unconditional state of the global system ([P n Q ν , ρ SE ] = 0). Subsequently, we evolve the compound system during the time interval [0, τ ]. The corresponding unitary operator U Λ is generated by the Hamiltonian H(t) = H(λ t ), which depends on time through an external parameter λ t that we vary according to a prescribed protocol Λ = {λ t : 0 ≤ t ≤ τ }: where T + denotes the time-ordering operator. As a result, the compound system at time t = τ is described by the new density matrix which in general contains classical and quantum correlations. The reduced (or local) states of the system and the environment can be obtained by partial tracing: To complete the process, a second local projective measurement is performed at time t = τ on both the system and environment. The measurement operators are arbitrary (rank-1) orthogonal projectors, denoted as {P * m } and {Q * µ }. Unlike in the first measurement, in this case the average global state is disturbed, transforming into Notice also that this is not a product state: the final local measurement does not eliminate the classical correlations contained in ρ * SE [58]. However, the measurement collapses the local states of the system and environment into pure states |ψ * m ψ * m | S ≡ P * m and |φ * µ φ * µ | E ≡ Q * µ . Thus, the spectral decomposition of the reduced states after the final measurement is where p * m = µ FIG. 1. Schematic picture of the forward process presented in the main text. System and environment start from an uncorrelated state ρS ⊗ ρE. A local measurement of observables with projectors {Pn, Qν } is carried out, which does not alter the density matrix in the average evolution but selects a pure state |ψn ⊗|φν at the trajectory level. System and environment then interact with each other according to the unitary evolution UΛ, ending in an entangled state ρ SE . Finally, we measure again now using projectors {P * m , Q * µ .}. In the last measurement quantum correlations in state ρ SE are erased, while the final state ρ * SE may still have in general non-zero classical correlations. The reduced evolution of the system conditioned to the measurement in the environment is described through the quantum operation Eµν (shaded area).
action of a quantum CPTP map E that admits a Kraus representation [59] with a set of Kraus operators M µν satisfying µ,ν M † µν M µν = I. There exist many Kraus representations {M µν } that reproduce the reduced dynamics on the system. We choose This specific representation retains all the details of the evolution of the environment, relating unequivocally each Kraus operator M µν with a transition |φ ν E → |φ * µ E in the environment. This is a key point in order to characterize the thermodynamics of the process at the trajectory level, as we will see shortly. Let us finally define the quantum operation: which describes the conditioned evolution of the system when the environment starts in the pure state |φ ν E and ends in the state |φ * µ E after measurement [60].

C. Average entropy production
We now discuss the entropy change along the process. We analyze here the von Neumann entropy, S(ρ) = −Tr[ρ ln ρ] of the global system. Recall that the von Neumann entropy coincides with the thermodynamic entropy for equilibrium states (setting the Boltzmann constant k = 1). For non-equilibrium states there are some situations where the von Neumann entropy still can be interpreted as a thermodynamic entropy [61]. However, in this paper we will refrain from identifying S(ρ) with a thermodynamic entropy and refer to it simply as the entropy or the quantum entropy of state ρ.
Along the process described above, the quantum entropy of the global system changes as This quantity is the quantum entropy production along the process. We will refer to ∆ i S inc as the inclusive entropy production to distinguish it from the entropy production when the system and the environment are separated at the end of the process and the final classical correlations are lost (see below). The inclusive entropy production is always non-negative, since von Neumann entropy cannot decrease in a projective measurement and stays constant along any unitary evolution, i.e., S(ρ SE ) = S(ρ SE ) ≤ S(ρ * SE ). Notice also that S(ρ) equals the classical Shannon entropy of the probability distribution of pure states in the eigenbasis of ρ. In particular we have To express the entropy of the global state in terms of local entropies and correlations, one can use the mutual information. For an arbitrary state σ SE with reduced states σ S and σ E , the mutual information is defined as (13) Here we have introduced the quantum relative entropy, S(ρ||σ) ≡ Tr[ρ(ln ρ − ln σ)], a non-symmetric and nonnegative measure of the distinguishability between states ρ and σ, which vanishes if and only if ρ = σ [62]. This property implies that mutual information becomes zero only for product (uncorrelated) states σ SE = σ S ⊗ σ E . Using mutual information, the inclusive entropy production can be rewritten as: where we have taken into account that the initial state is uncorrelated and, therefore, I(ρ SE ) = 0. The second equality shows that there are two sources of entropy production. The first one is the measurement disturbance of the final local states ρ S → ρ * S and ρ E → ρ * E , which can be avoided only by measuring in the eigenbasis of the reduced states ρ S and ρ E . The second source, captured by the term I(ρ SE ) − I(ρ * SE ) ≥ 0, is the erasure of quantum correlations in the state ρ SE . This is due to the local character of the measurements, being zero only if the global interaction U Λ does not generate quantum correlations [63,64].
In most situations the classical correlations remaining after the final measurement are irreversibly lost, with an entropic cost equal to the mutual information I(ρ * SE ). This is the case if we separate system and environment after the process and all subsequent manipulations are local and do not incorporate any feedback based on the outcomes of the final measurements. The entropy production in those situations is We will refer to ∆ i S as the non-inclusive entropy production or simply entropy production. The positivity of the non-inclusive entropy production in Eq. (15) has been already identified with the second law [48] and the existence of a thermodynamic arrow of time [65,66]. Notice that ∆ i S ≥ ∆ i S inc ≥ 0, since the mutual information I(ρ * SE ) is always non-negative. The differences between inclusive and non-inclusive entropy production will be illustrated in a specific example in Sec. VI A.

A. Forward and backward trajectories
We now extend the previous analysis to stochastic entropy changes at the level of individual quantum trajectories. A trajectory γ of the process introduced in the previous section (hereafter, we will call it the forward process) is simply given by the outcome of the four measurements, i.e., γ = {n, ν, µ, m}. This trajectory corresponds to the following transition between pure states Notice that, in virtue of our choice of the Kraus representation for the reduced dynamics [Eq. (8)], a trajectory γ is also a trajectory of the reduced dynamics, where the pair (ν, µ) now indicates the Kraus operation affecting the system instead of the initial and final states of the environment (which is otherwise hidden in the reduced dynamics). The probability to observe that trajectory γ is given by To introduce the backward process, we make use of the anti-unitary time-reversal operator in quantum mechanics, Θ, satisfying ΘΘ † = Θ † Θ = I and Θi = −iΘ. This operator changes the sign of odd variables under time reversal, like linear and angular momenta or magnetic fields [6,67]. We will consider the separate time reversal operators for system, Θ S , and environment, Θ E , as well as the one for the total bipartite system Θ = Θ S ⊗ Θ E .
The backward process is defined as follows. We start with a generic initial state of the form As in the forward process, the first step at time t = 0 is a local measurement of the family of projectors According to Eq. (18), the outcomes m and µ are obtained with probability˜ mµ . We then let the global system evolve under the Hamiltonian ΘH(λ t )Θ † inverting the time-dependent protocol as Λ ≡ {λ t | 0 ≤ t ≤ τ } withλ t = λ τ −t . This evolution is given by the unitary transformation Finally, at time t = τ we perform new local measurements on the system and environment using projectors and the corresponding backward trajectoryγ = {m, µ, ν, n} occurs with probabilitỹ

B. Fluctuation theorem
The unitary transformations corresponding to the forward and the backward process satisfy the so-called micro-freversibility principle for non-autonomous systems [6,68]: This is the key property that relates the probabilities of trajectories γ andγ in a quantum FT. By comparing the probabilities (17) and (21), using micro-reversibility (22) and the cyclic property of the trace, we immediately get where we have defined the quantities The terms in Eq. (24) are related to entropy changes per trajectory in the system and the environment, whereas I m,µ in Eq. (25) corresponds to the stochastic version of the mutual information [25,69] in the initial state of the backward process (18). From the detailed FT in Eq. (23), we immediately have the integral version and, using Jensen's inequality e x ≥ e x , one obtains a second-law-like expression The interpretation of ∆ i s γ depends on the choice of ρ SE , the initial global state of the backward process. If we setρ SE = Θρ * SE Θ † , then˜ mµ = * mµ and ∆ i s γ becomes the inclusive entropy production per trajectory. Its average equals the inclusive entropy production defined in Eq. (10). If the initial condition for the backward process is the uncorrelated stateρ SE = Θ(ρ * S ⊗ ρ * E )Θ † , theñ mµ = p * m q * µ and ∆ i s γ is the non-inclusive entropy production per trajectory, whose average yields the entropy production defined in Eq. (15) A third choice sets the environment in the (inverted) initial state of the forward process,ρ SE = Θ(ρ * S ⊗ ρ E )Θ † , which yields˜ mµ = p * m q µ . In this case both initial and final local measurements in the environment are performed in the same basis Q * µ = Q µ , and we obtain which includes an extra contribution measuring the disturbance on the environment during the process. The term S(ρ * E ||ρ E ), unlike S(ρ * E ) − S(ρ E ), is negligible when the environmental state is modified only infinitesimally (see appendix A), as is the case e.g. of a large reservoir. Moreover, when ρ E is a Gibbs state, Eq. (29) is the entropy production proposed in Ref. [50], and S(ρ * E ||ρ E ) corresponds to the thermodynamic entropy production due to irreversibly reseting the ancilla back to ρ E in contact with an equilibrium reservoir at the same temperature. Finally, we stress that for equilibrium canonical initial conditions both in the forward and in the backward processes, the trajectory entropy production (23) equals the stochastic dissipative work and one recovers the celebrated Crooks work theorem and the original Jarzynski equality [6,7].

IV. DUAL PROCESSES: ADIABATIC AND NON-ADIABATIC ENTROPY PRODUCTION
We now focus on the reduced dynamics. Our aim is to obtain FT's involving only the quantum trajectory defined in Sec. III and the initial and final states of the system. To do that, we follow our previous work [33], where we derived a FT for CPTP maps based on the dual dynamics introduced by Crooks in Ref. [70]. Remarkably, the resulting FT goes beyond the one that we have obtained considering the global dynamics, Eq. (23), as it will reveal an interesting split of the total entropy production into two terms: the adiabatic entropy production, which accounts for the irreversibility of the stationary regime, and the non-adiabatic entropy production, which measures how far the system is from that stationary state.
We apply the formalism in Ref. [33] to E, the map governing the reduced dynamics of the process, as well as to the map corresponding to the backward dynamics. Therefore, we first need to introduce the reduced dynamics in the backward process, which will be described by a new CPTP mapẼ. To do that, it is necessary that the system and the environment start the backward process in an uncorrelated stateρ SE =ρ S ⊗ρ E , i.e., we have to imposeĨ mµ = 0 [see Eq. (25)]. Otherwise the CPTP map of the backward reduced dynamics would depend on the initial state of the system. In that case, similarly to our choice (8) for the forward process, a useful representation ofẼ isẼ where the backward Kraus operators are given bỹ Notice that here we have swapped subscripts with respect to the definition of the forward operators given by Eq. (8). This can be done since the pair (µ, ν) is just a label of the Kraus operator. The choice in Eq. (31) means that the operationẼ νµ is equivalent to obtaining µ in the initial measurement of the backward process and ν in the final one. Now, micro-reversibility (22) implies an intimate relationship between the forward and backward Kraus operators: It is important to notice that the FT for the total entropy production (23) can be derived directly from the above equation. In other words, Eq. (32) expresses the fundamental symmetry under time reversal yielding the FT.
A. The dual-reverse process and non-adiabatic entropy production FT In order to go beyond the FT for the total entropy production, we proceed as in Refs. [33,70]. These works, inspired by classical stochastic thermodynamics, introduce a quantum dual dynamics that reveals the irreversibility associated to a map at the steady state. In the following we denote π an invariant state of the forward map, E(π) = π, which we indeed assume to be positive definite. The dual dynamics, which we will call here dual-reverse following the criterion used for classical systems [21][22][23], is defined as a mapD(ρ) such thatπ ≡ Θ S πΘ † S is an invariant state, i.e.,D(π) =π. Furthermore, when the map is applied several times starting from the stationary stateπ, it generates trajectoriesγ distributed asP D (γ|π) = P (γ|π). Here the trajectories are γ = {n, (ν 1 , µ 1 ), . . . , (ν N , µ N ), m} and γ = {m, (µ N , ν N ), . . . , (µ 1 , ν 1 ), n}, corresponding to N applications of the map.
Summarizing, in the stationary regime the dual-reverse generates the same ensemble of trajectories as the forward process, but reversed in time. For instance, if the map describes the dynamics of a system in contact with a single thermal bath (thermalization), then the forward process generates reversible trajectories (indistinguishable from their reversal) and the dual-reverse coincides with the forward map. In nonequilibrium situations, the dual generically inverts flows. For instance, for a system in contact with two thermal baths at different temperatures, the dual-reverse is usually obtained by swapping the temperatures of the baths, hence inverting the flow of heat.
In any case, one can prove that a Kraus representation of the dual-reverse map is given by the operators [33,70]: Finally, the dual-reverse process is the dual-reverse map complemented by a specific choice of the initial condition for the system (the environment does not appear explicitly in the dual map, which acts only on the system). The appropriate initial condition for the dual-reverse process isρ S , i.e., the same as the backward process. We now have three processes: the forward E, the backwardẼ, and the dual-reverseD, each one inducing an evolution in the system characterized by trajectories γ = {n, ν, µ, m}. We can compute the probability of observing a trajectory γ or its reverseγ = {m, µ, ν, n} in each of those evolutions. With a self-explanatory nota-tion, these probabilities read To obtain FT's from these expressions we need a condition of proportionality between operators M † µν , andD νµ , similar to the relationship (32) between M † µν andM νµ . In [33] inspired by [71], we found that a necessary and sufficient condition for that proportionality is the following. We first define the nonequilibrium potential Φ = − ln π, from the invariant state π. Its spectral decomposition reads: where φ i = − ln π i , and π i and {|π i } are, respectively, the eigenvalues and eigenstates of the invariant density matrix π. Now we require that each Kraus operator M µν is unambiguously related to a nonequilibrium potential change ∆φ µν (note however that the converse statement is not necessarily true, i.e. we may have for different values of µ and ν the same value of ∆φ µν ). In the invariant state eigenbasis this condition reads: with m µν ij = 0 whenever φ j − φ i = ∆φ µν . As pointed in [33], this condition does not imply single jumps between pairs of π eigenstates, but it could account for any set of correlated transitions between different pairs with same associated ∆φ µν . An extreme example are unital maps, where π is proportional to the identity matrix. In that case, ∆φ µν = 0 and any complex coefficients m µν ij satisfy Eq. (38). It is not hard to show that condition (38) is equivalent to [33] [Φ, This alternative formulation of (38) indicates that, when ∆φ µν = 0, M µν can be interpreted as ladder operators in the eigenbasis of the invariant state π. For thermalization or Gibbs preserving maps, with π = e −β(H−F ) , β = (kT ) −1 being the inverse temperature and F the equilibrium free energy, the potential is Φ = β(H − F ) and kT ∆Φ is the energy transfer between the system and the environment, i.e., the heat. Condition (38) in this case implies that the Kraus operators produce jumps between levels with the same energy spacing or, equivalently, jumps with a well defined value of the heat.
Introducing condition (38) in Eq. (33), one easily derives the following relationship between the forward and the dual-reverse Kraus operators [33]: and, using Eq. (32), one gets: Finally, inserting (40) and (41) into the expressions for the probability of trajectories (34-36) we obtain the following FT's: We will call ∆ i s a µν the adiabatic entropy production and ∆ i s na µν the non-adiabatic entropy production, following the terminology used in classical stochastic thermodynamics [21][22][23]. They contribute to the total entropy production per trajectory, ∆ i s γ = ∆ i s a µν + ∆ i s na γ , as defined in Eq. (23). Below we discuss the averages of the adiabatic and non-adiabatic entropy production in some cases, clarifying the origin of the terms.

B. The dual process and adiabatic entropy production FT
Notice that (43) is not a proper FT for the forward process. In particular, we cannot derive a Jarzynski-like equality for exp(∆ i s a µν ) averaged over forward trajectories, P (γ). To achieve this goal we need a further assumption that will allow us to apply the results of Ref. [33] to the backward process. In this way, we will obtain the dual-reverse of the backward process, which we simply call the dual map D. If condition (38) is satisfied, then, by virtue of (32), the backward Kraus operators can be written as: withm νµ ij ≡ e −σ E µν /2 (m µν ji ) * . We observe that, setting ∆φ νµ = −∆φ µν , condition (38) is recovered for the backward process. However, a requirement to apply the theoretical framework developed in Ref. [33] is that |π ≡ Θ S |π is an invariant state of the backward map E. This is not warranted by the definition ofẼ, not even when the Kraus operators are of the form (38). Therefore, we have to add this extra assumption. In particular, it is satisfied when the driving protocol associated to the map is time-symmetric, the Hamiltonian of the environment is invariant under time reversal, and we perform the same measurements at the beginning and the end of the process on the environment. This is the case of the infinitesimal maps that govern the dynamics of a quantum Markov process since, even in the case of arbitrary driving, each map is generated by a constant Hamiltonian.
We now obtain the dual operators D µν , applying transformation (33) to the backward Kraus operatorsM νµ (with the role of Θ S and Θ † S swapped [33]). Similarly to (40), condition (38) on the backward operators imply and, using Eq. (32), The dual process is given by the dual map with initial condition ρ S . The trajectories generated by this process are distributed as Combining Eqs. (34) and (47) and using condition (32), we get a new FT for the adiabatic entropy production:

C. Integral fluctuation theorems
We can now derive integral FT's for the adiabatic and non-adiabatic entropy productions: which follow from the detailed versions by averaging over trajectories γ. Finally, convexity of the exponential function provides the following two second-law-like inequalities as a corollary ∆ i s na γ ≥ 0 and ∆ i s a γ ≥ 0. As for the FT for the total entropy production (23), the meaning of these average entropies becomes clearer if the initial condition of the backward process is specified. Setting ρ = Θ(ρ * S ⊗ ρ * E )Θ † , the average of the adiabatic and nonadiabatic entropy production defined by (40) and (41) reads (51) and the sum equals the total non-inclusive average entropy production ∆ i S [see Eq. (15)]. It is interesting to notice that the average change of the potential, can be alternatively written in terms of averages over the states of the system, ρ S and ρ S if condition (38) is fulfilled. That condition implies [Φ, M µν ] = M µν ∆φ µν (see also Ref. [33]). Introducing the commutator in (52), where we have used the cyclic property of the trace and Eqs. (7) and (8). Therefore, the average potential change ∆φ can be expressed as the change in the expected value of the operator Φ due to the map. Recall that the operator Φ acts on the Hilbert space of the system H S , i.e., is a local observable on the system. If the final measurement does not alter the state of the system, i.e., if ρ * S = ρ S , or if the final measurement is skipped, as it is the case when we concatenate maps and the system is measured only after the whole concatenation (see Sec. IV E below), we can write the average non-adiabatic entropy production in an appealing form: where we have used the definition Φ = − ln π of the potential operator in terms of the invariant state π. Here we see that the non-adiabatic entropy production is related to the distance between the state of the system and the invariant state π. During the evolution, the state of the system can only approximate the invariant state and the non-adiabatic entropy production is a measure of the irreversibility associated to such convergence. In fact, inequality in Eq. (54) follows from direct application of Ulhman's inequality (monotonicity of quantum relative entropy) holding for general CPTP evolutions [51,62].

D. Multipartite environments
The results obtained so far can be also applied to multipartite environment. The corresponding Hilbert space is decomposed as H E = R r=1 H r , corresponding to R independent ancillas or reservoirs interacting with the open system. We assume an uncorrelated initial state of the environment, ρ E = ρ 1 ⊗ ... ⊗ ρ R , and that the measurements are performed locally in each environmental ancilla.
In this case the adiabatic entropy production per trajectory and its average read (see details in appendix B):

E. Concatenation of CPTP maps
Up to now, we have considered a single interaction of duration τ between the system and the environment [see Eq.
(2)]. The CPTP map E describes the evolution of the system when the environment is measured before and after interaction. This framework can be extended to concatenations of CPTP maps, where the system interacts sequentially with the environment. Each single interaction in an time interval [t, t + τ ] is described by a single CPTP map like E. The map describing the reduced dynamical evolution for N interactions, from t = 0 to t = N τ , is: where, in particular, each map E (l) may have a different (positive-definite) invariant state π (l) . We assume that the system interacts from time t l−1 ≡ (l − 1)τ to time t l ≡ lτ with a "fresh" (uncorrelated) environment in a generic state ρ α , and, as in the single map case, the environment is measured before and after interaction with the system by projective measurements. On the other hand, the system is only measured at the beginning and end of the the whole concatenation (57), as depicted in Fig. 2.
In this case, trajectories are specified by the set of outcomes γ = {n, (ν 1 , µ 1 ), ..., (ν N , µ N ), m}, which can be compared to the backward trajectoriesγ = {m, (ν N , µ N ), ..., (ν 1 , µ 1 ), n} generated by the reverse se- . We find that all our above results apply as well to the concatenations setup (see appendix C) yielding the following three detailed fluctuation theorems: where σ S nm is given by Eq. (24), σ E µ l ν l = ln q (l) is the entropy change in the environment due to the lth map, and ∆φ ν is the change in non-equilibrium potential for the l-th map.

V. LINDBLAD MASTER EQUATIONS
The results of the last section can be applied to Lindblad master equations [31,52]. Consider the following master equation in Lindblad form [4, 60,72], depending on an external parameter λ t : where H(λ t ) is the system Hamiltonian in the selected picture and L k (λ t ) are positive Lindblad operators, which generally depend on the control parameter λ t and describe jumps in some (possibly time-dependent) basis. We assume that there exists an instantaneous invariant state π λ , which is the steady state of Eq. (61) when the external control parameter is frozen: L λ π λ = 0 [5]. Λ depending on the protocol Λ l . The ancilla is measured before and after interaction generating outcomes ν l and µ l respectively.
The Lindblad equation (61) can be written as a concatenation of CPTP maps with the Kraus representation Recall that this Kraus representation is not unique [60]. As before, the representation (63-64) is related to a specific detection scheme for the jumps, that is, it implies a specific choice of the initial state and the local observables to be monitored in the environment (the set of projectors {Q ν } and {Q * µ }). The Kraus representation (63-64) is based on a family of operations M k with k = 1, . . . , K that induce jumps in the state of the system and occur with probabilities of order dt, and a single operation M 0 that induces a smooth nonunitary evolution and occurs with probability of order 1. This implies that a trajectory γ consists of a large number of 0's punctuated by a few jumps M k with k = 1, . . . , K. An alternative way of describing the trajectory is to specify the jumps k j and the times t j where they occur, i.e., γ = {n, (k 1 , t 1 ), ..., (k j , t j ), ..., (k N , t N ), m}, where, as before, n and m denote the outcomes of the initial and final measurements in the system at times t = 0 and t = t f . Jump k is given by the operation E k (ρ) ≡ M k ρM † k , whereas between two consecutive jumps at t j and t j+1 the evolution is given by the repeated application of the operation corresponding to the Kraus operator M 0 (λ t ) in (63). This results in a smooth evolution given by the operator: with an effective non-hermitian Hamiltonian that reads In this representation, the probability of a trajectory γ = {n, (k 1 , t 1 ), ..., (k j , t j ), ..., (k N , t N ), m} is with U tj+1,tj (ρ) = U eff (t j+1 , t j )ρ U † eff (t j+1 , t j ).

A. Backward, dual and dual-reverse dynamics
Consider now the backward dynamics. The timeinversion of the evolution of the global system corresponds to a time-reversed version of the Lindblad master equation (61). As in the previous section, the backward process is generated by inverting the sequence of operations together with time-inversion of each operation in the sequence. The map corresponding to an infinitesimal time-step in the time-reversed dynamics, ρ t+dt =Ẽ(ρ t ), admits a Kraus representation with Kraus operatorsM k (λ t ). To obtain the backward map, we would need to know details about the environment that induces the Markovian dynamics given by the Lindblad equation (61). However, in the previous sections we have derived a relationship between the forward and backward CPTP maps, namely Eq. (32): Imposing the backward maps to be trace-preserving, that is,M † 0M 0 + kM † kM k = I, we obtain σ E 0 = 0, and the consistency condition Any set of numbers {σ E k } satisfying Eq. (70) defines, through (68), an admissible backward process. The existence of such set is warranted, since any Lindblad equation can be derived from the interaction between the system and an ancilla.
For any trajectory γ = {n, (k 1 , t 1 ), ..., (k N , t N ), m} generated in the forward process with probability P (γ), there exist a backward trajectoryγ = {m, (k N , t N ), ..., (k 1 , t 1 ), n} occurring in the backward process with probabilityP (γ). The backward trajectory can also be identified by the times of successive jumps. In this representation, the probability of trajectoryγ can be written as: whereẼ k (ρ) =M kρM † k . The smooth evolution between jumpsŨ t ,t (ρ t ) =Ũ eff (t , t)ρ tŨ † eff (t , t) is given by the operator where {λ t } again corresponds to the inverse sequence of values for the control parameter. It can be shown that this smooth evolution obeys the micro-reversibility relationship Θ † SŨ eff (t , t)Θ S = U eff (t , t) † . Let us discuss now the dual and dual-reverse dynamics. The condition (39), necessary to define the dual-reverse process, reads [31,52]: These commutation relationships indicate that the Lindblad operators L k (λ t ) promote jumps between the eigenstates of π λt . Furthermore, as the condition must be fulfilled for the operator M 0 in Eq. (63) as well, we need [H, k L † k L k ] = [H, Φ] = 0, which in turn implies ∆φ 0 = 0. This means that the instantaneous steady state of the dynamics must be diagonal in the basis of the Hamiltonian term appearing in Eq. (61). This condition is fulfilled by equilibrium Lindblad equations and in situations in which the operator H becomes the identity operator in an appropriate interaction picture (see e.g. Refs. [45,73]). However, the condition can be broken in nonequilibrium situations, a genuine quantum effect.
In section VI C we present an example of a periodically driven cavity mode where the adiabatic entropy production can be negative. Finally, as discussed in section IV B, we recall that the fluctuation theorem for the adiabatic entropy production can be stated when the backward mapsẼ admitsπ λ ≡ Θ S π λ Θ † S as an invariant state. If these conditions are fulfilled, the dual process is defined by the dual operations D k (·) = D k (·)D † k with Kraus operators {D k } as defined in Eq. (46), whereas the dualreverse process is given by operationsD k =D k (·)D † k with Kraus operators {D k } defined in Eq. (40) (see also Eqs. (C7-C8) in App. C). The probability of a trajectory γ in the dual process, P D (γ), can be calculated from Eq. (67) by using the same map U t ,t for the no-jump time evolution intervals, and replacing the operations E k by the dual operations D k . Analogously, for the dual-reverse process the probability of trajectoryγ,P D (γ), can be constructed from Eq. (71) withŨ t ,t for the no-jump evolution, and dual-reverse operationsD k . We further notice that in general D k = M k , andD k =M k , that is, σ E k = −∆φ k . In many applications, the Lindblad operators come in pairs and the corresponding pair of terms in the sum (70) cancel. That occurs if, for a specific pair of operators and Γ j (λ t ) being positive transition rates, and L(λ t ) some arbitrary (possibly time-dependent) system operator. Then, condition (70) implies (cf. [74]) and the (inverted) Kraus operators of the backward map are also operators of the forward map: where we have used the detailed-balance relation (32). Moreover,π λ ≡ Θπ λ Θ † is invariant under the backward map:Ẽ B. Entropy production rates The above considerations lead us to reproduce the three different detailed FT's in Eqs. (60) for quantum trajectories generated by Lindblad master equations. From the integral fluctuation theorems we can derive secondlaw-like inequalities analogous to Eqs. (50) and (51) for the entropy production rates [31]: whereṠ = −Tr[ρ t ln ρ t ] is the derivative of the von Neumann entropy of the system,φ ≡ Tr[ρ t Φ(λ t )] = −Tr[ρ t ln π λt ] is the nonequilibrium potential change rate, and σ E (λ t ) dt = k l Tr[E k l (ρ t )]σ E k l (λ t ) the entropy change in the monitored environment during dt [52]. The three above equations guarantee the monotonicity of the average entropy production, ∆ i S, and the adiabatic and non-adiabatic contributions, ∆ i S na and ∆ i S a , during the whole evolution.
The physical interpretation of the adiabatic and nonadiabatic entropy production now becomes clear. The non-adiabatic part can be written aṡ which is the continuous time version of Eq. (54). If the control parameter changes quasi-statically, we have ρ t π λt and, therefore, the non-adiabatic entropy production vanishes. This is analogous to the classical non-adiabatic entropy production introduced in Refs. [20][21][22][23]. On the other hand, the adiabatic contributionṠ a is in general different from zero even if the driving is extremely slow. In a physical system, this term accounts for the entropy production required to keep the system out of equilibrium when λ is fixed, and the associated dissipated energy is usually referred to as housekeeping heat [20]. At this point, it is worth it to remark an important difference between classical and quantum systems. In classical systems, the split of the entropy production in two terms, adiabatic and non-adiabatic, can always be done at the level of trajectories, and both terms obey fluctuations theorem that ensure the positivity of their respective averages. This is possible for quantum systems only if (39) [or (73) for Lindblad operators] is met. One can still use (79) as a definition for the average non-adiabatic entropy productionṠ na andṠ a =Ṡ i −Ṡ na for the average adiabatic entropy production rate. However, these definitions cannot be extended to single trajectories and, furthermore, they do not obey a fluctuation theorem. In the next section, we discuss a specific example where the condition is not fulfilled and, as a consequence, the average adiabatic entropy production rate can be negative.
Finally, it is also important to notice thatṠ and σ E in Eqs. (77) are exact differentials, i.e., can be written as the time derivative of the system and the environment entropy, respectively. On the other hand, the termφ ≡ Tr[ρ t Φ(λ t )], as well as the adiabatic and nonadiabatic entropy production rates in Eq. (78), cannot be expressed, in general, as a time derivative. One important exception is the case of a constant invariant state π λt = π like, for instance, in a relaxation in the absence of driving. In that case, all the quantities in Eqs. (77)(78) are exact differentials. In particular, the non-adiabatic entropy production when the system relaxes from ρ 0 to ρ t is given by ∆S na = −S(ρ τ ||π) + S(ρ 0 ||π) ≥ 0 (80) which equals ∆S na = S(ρ 0 ||π) for a full relaxation to ρ τ = π. The later coincides with the entropy production introduced by Spohn [75].

VI. EXAMPLES
We illustrate our findings with three paradigmatic examples. In the first one, we consider a two-qubit CNOT gate as a simple process with a finite size environment to illustrate the differences between the inclusive and noninclusive entropy production introduced in section II C. The second and third examples correspond to two representative examples of non-equilibrium quantum Markov systems. The second example is an autonomous system coupled to several thermal baths. In this case, the nonadiabatic entropy production is zero except during the transient relaxation to the steady state. However, it provides an intuitive picture of how entropy is produced in non-equilibrium setups. The third example is a driven system that does not fulfill condition (38) and, consequently, does not admit the splitting of the entropy into adiabatic and non-adiabatic contributions with positive averages.

A. A two-qubit CNOT gate
The difference between the inclusive and non-inclusive entropy production introduced in section II C becomes especially relevant for processes where the system of interest repeatedly interacts with a finite-size reservoir. As an extreme case we consider both system and environment to be qubits with the same energy spacing . Their Hamiltonians are given by H S = |1 1| S and H E = |1 1| E . We assume the initial state of the system to be partially coherent, ρ S = (I + α σ x )/2 with 0 ≤ α ≤ 1, and the environmental qubit starting in a thermal state ρ E = e −βH E /Z E at inverse temperature β = 1/k B T ≥ 0, Z E = 1 + e −β being the partition function. The initial state can be written as where κ ≡ tanh(β /2), σ j with j = x, y, z are the Pauli matrices, and we take the standard qubit basis {|0 , |1 } for both the system and the environment. The eigenbasis of ρ SE determines the projectors of the initial measurements {P n , Q ν }, which are in this case P ± = |ψ ± ψ ± | with |ψ ± = (|0 ±|1 )/ √ 2 and Q ν = |ν ν| with ν = 0, 1. System and environment interact through a CNOT gate, U CNOT , where the system acts as the control qubit [62]. The interaction leads to the following global systemenvironment state Notice that ρ SE has maximally mixed reduced states both in system and environment. As a consequence, for any choice of the final projectors {P * m , Q * µ }, we have ρ S = ρ * S = ρ E = ρ * E = I/2. In contrast, the global state ρ * SE depends on the final projectors. The average work done during the interaction is W = Tr[(H S + H E )(ρ SE − ρ SE )] = (1/2 − e −β /Z E ) > 0, while there is no further energy contributions from local measurements.
The inclusive entropy production in Eq. (14) is just given by the erasure of quantum correlations in the final measurements, ∆ i S inc = I(ρ SE ) − I(ρ * SE ). This is the so-called mutual induced disturbance introduced by Luo [63]. Moreover if, following Refs. [76,77], we maximize I(ρ * SE ) over {P * m , Q * µ }, then the inclusive entropy production is equal to the (symmetric) quantum discord [64,78] of the state ρ SE . On the other hand, the non-inclusive entropy production in Eq. (15) is given by the total correlations in state ρ SE , that is, ∆ i S = ∆ i S inc + I(ρ * SE ) = I(ρ SE ), and it is independent of the choice of the local projectors of the final measurements {P * m , Q * µ }. The entropy production per trajectory ∆ i s γ can be calculated as explained in Sec. III. Recall that we may obtain both the inclusive or non-inclusive entropy production depending on our choice for the initial state of the backward process, and that the two quantities verify the integral fluctuation theorem (26).
In Fig. 3 we show the probability distribution of the entropy production P (∆ i s γ ) for β = 2.5 and α = 0.8. Blue solid bars correspond to the non-inclusive version and purple dashed bars correspond to the inclusive one. The latter does depend on the final measurements. Here we have taken as final projectors the local energy eigenbasis, {P * m = |m m| S , Q * µ = |µ µ| E } for m, µ = 0, 1. The different types of average entropy production are plotted in the inset figure as functions of β for the same value of α = 0.8. There, black and blue solid lines correspond to the average non-inclusive entropy production with and without including the term S(ρ * E ||ρ E ) due to local disturbance of the environment [see Eq. (29)], respectively. Dashed and dotted lines show to the average inclusive entropy production for different choices of the local projectors in the final measurement {P * m , Q * µ }. The purple dashed line is obtained when the final projectors are given by the local energy eigenbasis {P * m = |m m| S , Q * µ = |µ µ| E } for m, µ = 0, 1. The orange dotted line is the symmetric quantum discord, obtained when maximizing I(ρ * SE ). As mentioned in section II C, inclusive and noninclusive entropy production apply to different physical situations depending on how the system and the environment are manipulated after the process. If system and environment are separated and every further manipulation is local, then we do not make use of the classical correlations given by the mutual information I(ρ * SE ); in this case the non-inclusive entropy production is the magnitude that adequately describes the increase of entropy. On the other hand, global operations on the whole system+environment can make use of those correlations and, for instance, extract more energy from a thermal bath. We illustrate this possibility in our simple example by considering a second CNOT interaction after the final local measurements. For simplicity we perform the final measurements in the local energy basis, {P * m = |m m| S , Q * µ = |µ µ| E } for m, µ = 0, 1. Applying these projectors to state ρ SE in Eq. (82), one obtains the final global state Applying the second CNOT to this state, one gets where ρ E is the initial thermal state of the environment. As we can see, in this second process system and environment become completely decorrelated after interaction, while a work W ext = Tr[(H S + H E )(ρ SE − ρ * SE )] = (1/2 − e −β /Z E ) is extracted when performing the second gate. Notice that the extracted work equals the work performed in the first gate. This work extraction is impossible if we only have local operations at our disposal, for which the final state ρ * SE is completely equivalent to the uncorrelated state ρ * S ⊗ ρ * E . This simple example highlights the importance of distinguishing between inclusive and non-inclusive entropy production in a small finite-size environment. Similar conclusions can be applied for the term S(ρ * E ||ρ E ).
FIG. 4. Schematic diagram of a three-level thermal machine acting as a refrigerator. The three transitions of the machine are weakly coupled to thermal reservoirs at temperatures β1, β2 and β3, inducing jumps between the machine energy levels (double arrows). In a refrigeration cycle the machine performs a sequence of three jumps |g → |eA → |eB → |g , where it absorbs a quantum of energy ω1 from the cold reservoir, together with a quantum ω2 from the hot one, while emitting a quantum ω3 into the reservoir at intermediate temperature.

B. Autonomous thermal machine
Consider an autonomous three-level thermal machine powered by three thermal reservoirs at different temperatures, as depicted in Fig. 4 [44,[79][80][81][82]. Each bath mediates a different transition between the energy levels, {|g , |e A , |e B }. The Hamiltonian of the system is that is, the three possible transitions g ↔ e A , e A ↔ e B and g ↔ e B have frequency gaps ω 1 , ω 2 , and ω 3 ≡ ω 1 + ω 2 , respectively. Each transition is weakly coupled to a bosonic thermal reservoir in equilibrium at inverse temperature β r = 1/kT r with r = 1, 2, 3, where we assume β 1 ≥ β 3 ≥ β 2 for concreteness. The dynamics of the three-level thermal machine can be described by a Lindblad master equation obtained in the weak coupling limit by standard techniques from open quantum systems theory [4, 5,60]. It readṡ where ρ t is the density operator of the three level system and Lamb-Stark shifts have been neglected. The three dissipative terms in the above equation describe the irreversible dynamical contributions induced by each of the three thermal reservoirs: r = 1, 2, 3 where a 1 = |g e A |, a 2 = |e A e B | and a 3 = |g e B | are the ladder operators of the three-level system. Equation (87) describes the emission and absorption of excitations of energy ω r to or from the reservoir r, at rates Γ (r) ↓ = γ r (n th r + 1) and Γ (r) ↑ = γ r n th r , fulfilling detailed balance Γ (r) ↓ = e βr ωr Γ (r) ↑ . Here n th r = (e βr ωr − 1) −1 is the mean number of excitations of energy ω r in reservoir r, and γ r ω r ∀r, r = 1, 2, 3 are the spontaneous emission decay rates associated to each transition. The heat fluxes entering from the reservoirs associated to the imbalance in emission and absorption can be defined asQ r = Tr[H S L r (ρ t )] for r = 1, 2, 3, and the first law of thermodynamics readsU = Tr[H SρS ] = Q 1 +Q 2 +Q 3 . Therefore, in our example we have six Lindblad operators (r = 1, 2, 3): that define the infinitesimal CPTP map (62) with the Kraus representation given by Eqs. (63)(64). In particular: Here the stochastic jumps during the evolution correspond to simple transitions between the energy levels {|g , |e A , |e B }. Therefore, the stochastic dynamics is completely equivalent to a classical Markov process if the initial state of the machine is diagonal in the Hamiltonian eigenbasis. In particular, the stationary state reads: In appendix D, we explicitly calculate the occupation probabilities π g , π A , and π B . Nevertheless, the transient dynamics could exhibit some quantum effects when the initial state exhibit coherences in the Hamiltonian eigenbasis. For instance it has been recently pointed that initial coherence can be used to reach lower temperatures during the transient dynamics [83,84]. t N ), . . . , (k 1 , t 1 ), n} is defined by the inverse sequence of events with respect to γ, occurring in the backward process. We consider the initial state of the backward process the inverted final state of the forward process, Θ S ρ t f Θ † S , while the thermal reservoirs have the same state as in the forward process. We further assume the simplest form for the time inversion operator Θ S , namely, the complex conjugation, i.e., Θ S ψ = ψ * , which commutes with any matrix with real entries, as the Hamiltonian and the jump operators a r , a † r . The Lindblad operators in this case come in pairs L (r) ↓ = e βr ωr/2 L (r) † ↑ . Hence the stochastic entropy change in the environment σ E k for each operator L k is given by Eq. (74), where the label k takes on the six possible values k = (↑, r) and k = (↓, r) with r = 1, 2, 3: This is as expected, since the upward jump r induced by the operator L (r) ↓ in the forward trajectory γ dissipates a heat ω r to the reservoir at inverse temperature β r . Equivalently, in the downward jump r a heat ω r is extracted from the thermal bath reducing its entropy by an amount β r ω r . The Kraus operators of the backward map are given by Eq. (75) for the effective evolution operator describing the dynamics between jumps in the backward process. From the above equations we see that the backward mapẼ is obtained from the forward map E inverting the jumps. We also notice that, consequently, the backward mapẼ admits the time-reversed steady statẽ π = Θ S πΘ † S = π as an invariant state. We next construct the dual and dual-reverse processes for the model. The condition for the Lindblad operators to be of the form in Eq. (38) is fulfilled here. Indeed, the non-equilibrium potential, Φ = − ln π, obeys [Φ, H S ] = 0 and where the nonequilibrium potential changes associated to each jump in the trajectory read ∆φ 0 = 0, ∆φ r↓ = −β r ω r , ∆φ r↑ = β r ω r . (94) Here we have introduced the quantities β 1 = ln( π0 π1 )/ ω 1 , β 2 = ln( π1 π2 )/ ω 2 and β 3 = ln( π1 π2 )/ ω 3 , which can be seen as the local inverse temperature (or virtual temperature [85][86][87]) of each transition in the steady state π. As shown in appendix D, they determine the direction of the heat flow in the stationary regime, i.e., if β r > β r , then the temperature of reservoir r is higher than the local temperature of the machine and the heatQ r is positive (energy flows from the reservoir to the machine), and viceversa. Moreover, for β r β r the heat flow is proportional to β r − β r ; therefore, the difference β r − β r can be considered as a thermodynamic force for the heat flow between the reservoir r and the system. In Fig. 5 we plot the local inverse temperatures β r compared to the reservoir temperatures β r for a specific choice of β 2 = 0.5 and β 3 = 4 and as a function of β 1 , the inverse temperature of the coldest bath (we use units of ( ω 1 ) −1 ). There is a point, around β 1 = 9.3, where β r = β r and all the heat flows in the stationary regime vanish. Below that point, the steady heat flow from the coldest reservoir at inverse temperature β 1 is positive, i.e., the machine acts as a refrigerator that extract energy from the coldest bath 1. On the other hand, for β 1 > 9.3, heat flows from the machine to the hottest bah at inverse temperature β 2 , so In the refrigerator regime, the transition g ↔ eA is at an effective temperature colder than the coldest reservoir, β 1 ≥ β1, inducing heat extraction from it, while the other transitions induce dissipation of heat to the reservoir at intermediate temperature, β2 ≥ β 2 , and absorption of heat in the hotter one β 2 ≥ β2. In the heat pump regime the three heat flows change its directions as the previous inequalities become inverted.
the machine acts as a heat pump capable to heat up the hottest reservoir 2.
The Kraus operators for dual and dual-reverse maps, D andD, can be obtained from Eqs. (46) and (41) respectively, by using Eqs. (92) and (94). They read: We see that the dual and dual-reversed dynamics induce similar jumps in the three-level system, but with modified rates depending on the differences β r − β r . Only when β r = β r for each r, the dual process becomes equal to the forward process, and hence the dual-reverse process equals the backward process (see Fig. 5).
Notice that Eq. (93), together with the backward map havingπ = π as an invariant state, are sufficient conditions to ensure the existence of the three fluctuation theorems for the adiabatic, non-adiabatic and total entropy productions during trajectory γ. They explicitly read: where σ S nm is the stochastic entropy increase in the system, and q ↓ ) is the stochastic heat entering the system from reservoir r, n (r) ↑↓ being the total number of upward or downward jumps in transition r. The expression for the adiabatic entropy production is particularly interesting, since it is equal to the entropy generated by the heat transfer between reservoirs at inverse temperatures β r and β r . In particular, the adiabatic entropy production is identically zero when β r = β r even though it is possible to have transient flows of heat.
We can now calculate the average rates of nonequilibrium potential and reservoirs entropy changes: where we split in three parts the nonequilibrium potential flowΦ =Φ 1 +Φ 2 +Φ 3 = −Tr[ρ S ln π]. The entropy production rates hence read: showing the same structure as the trajectory entropies in Eqs. (99-101). Since there is no driving in this example, the non-adiabatic entropy production reads as in Eq. (80), and equals ∆S na = S(ρ 0 ||π) for a full relaxation to the steady state π.
In the steady state we haveṠ na = 0, and the first law becomes rQ ss r = 0. This implies that the only contribution to the entropy production rate is the adiabatic one, which can be written as: This equation can be used to bound the efficiency of the machine in the different regimes of operation. For instance, the efficiency of the machine acting as a refrigerator is given by ... where C is the so-called Carnot efficiency of a refrigerator [85].

C. Periodically driven cavity mode
Our third example consists of a single electromagnetic field mode with frequency ω in a microwave cavity with slight losses in one of the two mirrors. The losses of the cavity are produced by the weak coupling of the cavity mode to a bosonic thermal reservoir in equilibrium at inverse temperature β = 1/kT . In addition, an external laser of the same frequency ω and weak intensity drives the cavity mode producing excitations. The Hamiltonian of the system can be expressed as H S (t) = H 0 + V S (t) consisting of two terms: the first one is the Hamiltonian of the undriven mode H 0 = ωa † a and describes the effect of the classical resonant laser field with complex amplitude = | |e iϕ . Here the subscript S stands for the Schrödinger picture, whereas operators and density matrices without any subscript will correspond to the interaction picture with respect to H 0 (H 0 is of course the same in the two pictures). Figure 6 shows a schematic picture of the setup. The reduced evolution of the cavity mode can be described by the following Lindblad master equation [60] with the dissipative part This term accounts for emission and absorption of photons by the cavity mode from the equilibrium reservoir at respective rates Γ ↓ = γ 0 (n th + 1) and Γ ↑ = γ 0 n th . Here again n th = (e −β ω − 1) −1 , and γ 0 is the spontaneous emission decay rate in absence of driving. The dissipative term L(ρ) does not depend on the driving: it induces jumps in the eigenbasis of H 0 and is also invariant under the change of picture. Notice that this is an approximation. For slow driving, for instance, the bath induces jumps between the instantaneous eigenstates of the Hamiltonian H S (t). The dissipator (111) is valid for weak driving and weak coupling with the thermal bath, that is γ 0 ∼ | | ω [88]. In the interaction picture with respect to H 0 , the Lindblad equation (110) reads [60] where V = i ( a † − * a) is the driving Hamiltonian in the interaction picture, which turns to be constant. Before discussing the FT applied to this example, let us calculate the energetics of the system from the Lindblad equation. For this purpose, it is more convenient to express the internal energy in the Schrödinger picture: U (t) = Tr[H S (t)ρ S (t)]. The first law readṡ U (t) =Ẇ (t) +Q(t), where the average work is given byẆ and the average energy change not accounted for by work we denotė although it is not necessarily equal to the heat, i.e., the energy reversibly exchanged with a thermal reservoir that accounts for the reservoir's entropy change [89]. Below we discuss in detail the physical nature of this energy transfer.
In contrast to the undriven case, here the cavity does not reach equilibrium with the reservoir: coherences in the energy basis do not decay to zero due to the work performed by the external laser. Notice also that the state π defines a limit cycle (unitary orbit) in the Schrödinger picture. In the stationary regime π S (t) = e −iH0t/ π e iH0t/ , i.e., the mode rotates in optical phase space, according to the free evolutionπ S = The energetics in this stationary regime is rather simple. The internal energy is constant, even though the state π S (t) depends on time: U ss = Tr[H S π S ] = Tr[(H 0 + V )π] = Tr[H 0 π] = ω(n th + |α| 2 ), bigger than the thermal average energy ωn th . The laser introduces energy at a rate: which is dissipated as heat to the thermal bathQ ss = −Ẇ ss . The FT can be applied both to the Srchrödinger or the interaction picture. Here it is is more convenient to determine the forward and backward processes in the interaction picture, where there is no driving. The Kraus operators for the map E in Eq. (62) read in this case: for the no-jump evolution, and for the downward and upward jumps corresponding to emission and absorption of photons. The trajectory γ = {n, (k 1 , t 1 ), ..., (k N , t N ), m} is then constructed as in the previous example by counting the jumps induced by the reservoir and registering the times at which they occur.
Since the forward dynamics is governed by a single pair of Lindblad operators {L ↓ = Γ ↓ a, L ↑ = Γ ↑ a † }, condition (70) allows us to obtain the stochastic entropy changes in the reservoir: That is, in a downward (upward) jump, the entropy in the environment increases (decreases) by β ω, corresponding to a transfer of energy ω. In average, this transfer of energy equals Tr[H 0 L(ρ(t))], whereas the energy not accounted for by work is given by Eq. (114), i.e., bẏ ]. The origin of this discrepancy is the choice of a dissipator (111) independent of the driving. As already mentioned, the dissipator is valid for weak driving [88], when the term Tr[V L(ρ(t))] ∼ O(γ| 0 |) is negligible. However, it is worth it to notice that our approach does not depend on the physical nature of the environment and its interaction with the system. As shown in section V, once a Lindblad equation like (110) with a given set of Lindblad operators for its unraveling has been specified, no matter how it has been derived, induces an entropy change in the environment given by Eq. (117). This is a direct consequence of micro-reversibility that yields condition (70) on the Lindblad operators. When these operators come in pairs, as it is the case in our example, the condition determines the entropy change in the environment [see Eq. (74)].
Therefore, if one could conceive physical situations where the Lindblad equation (110) is exact, then the entropy production would be given by Eq. (117) and the energy transfer Tr[V L(ρ(t))] would not contribute to the entropy of the environment. A clue on the nature of this energy transfer is provided by Ref. [55]. In that paper, Elouard et. al. introduce a driven two-level system in an engineered thermal bath where excitations occur through a third level with a very short lifetime and transitions are monitored by measuring emitted photons. The resulting Lindblad equation is the analog of Eq. (110) for a twolevel system and the entropy change in the environment is precisely (117). They show that the energy transfer Tr[V L(ρ(t))] is due to the collapse of a coherent state induced by the photon detection. This energy transfer does not change the entropy of the universe and has been categorized either as "measurement work" [31,74] or as "quantum heat" [54,55] due to measurement.
The Kraus operators of the backward map are given by Eq. (75) implying again that the forward and the backward map are equivalent, i.e., the jumps up (down) in the forward process are transformed in jumps down (up) in the backward process.
The main feature of this example is that the key condition (38) is not fulfilled. Recall that this condition is needed to define the dual and dual-reverse dynamics as well as the stochastic adiabatic and non-adiabatic entropy production at the trajectory level. Using the expression for the stationary state Eq. (115), we can calculate the nonequilibrium potential in the interaction picture where we have introduced the field quadrature The nonequilibrium potential Φ in Eq. (121) does not obey Eq. (73), because the Lindblad operators appearing in the dynamics (112) promote jumps among the eigenstates of the unperturbed Hamiltonian H 0 , instead of the eigenstates of the steady density matrix π. This implies that we cannot associate a single change in the nonequilibrium potential to each of the Lindblad jump operators, nor to M 0 . As a consequence, the entropy production per trajectory cannot be decomposed in an adiabatic and a non-adiabatic contributions and the corresponding fluctuation theorems do not apply. However, the conditions in Eq. (73) can be recovered in some cases by properly including the effect of the driving on the Lindblad operators [74].
Even though the adiabatic and non-adiabatic entropy production cannot be defined at the trajectory level, we can calculate their averages [75] using, for instance, Eqs. (78): Here we have usedṠ i =Ṡ − βQ andU (t) = Ḣ 0 = Q +Ẇ , and introducedẊ ϕ ≡ Tr[x ϕρ (t)], the rate at which the cavity field mode is displaced by the laser (with phase ϕ) until the steady state is reached. Since there are no fluctuation theorems for these quantities, in principle they could be negative. The non-adiabatic entropy production, however, still obeys (79) and, since the steady state π is constant in the interaction picture, it can be written as the change of the relative entropy between the state ρ(t) and π, which is always positive: S na = −Ṡ(ρ(t)||π) ≥ 0. This is not the case of the adiabatic entropy productionṠ a , which indeed can take on negative values. The expression forṠ a in Eqs. (123) does equal the entropy production in the steady state, S a → βW ss = −βQ ss ≥ 0. However, it can be negative in the transient regime, as shown in Fig. 7 (see also appendix E). In Fig.7(a) we depict the evolution of the three entropy production rates when the cavity mode starts in a Gibbs thermal state in equilibrium with the reservoir temperature, ρ 0 = exp(−βĤ 0 )/Z 0 . We find that the entropy of the mode is kept constant during the evolution,Ṡ = 0 ∀t, which impliesṠ i = −βQ ≥ 0, anḋ S na = β( ω|α|Ẋ ϕ −U ) ≥ 0. On the other hand, the adiabatic entropy production rateṠ a = β(Ẇ − µẊ ϕ ) is negative for times t < t n ≡ 2 ln 2/γ 0 . It is worth mentioning that for this initial condition the term Tr[V L(ρ(t))] in the energetics vanishes at any time t.
To explore the origin of this purely quantum effect, we plot the energetics of the relaxation process Fig. 7(b). The cavity mode absorbs energy at constant entropy from the external laser until the periodic steady state is reached,U =Ẇ e −γ0t/2 , whereẆ =Ẇ ss (1−e −γ0t/2 ) ≥ 0 is the input power and, consequently, heat is dissipated at a rate −Q =Ẇ (1 − e −γ0t/2 ) ≥ 0. When the relaxation is completed, the input laser power is fully dissipated into the reservoir, i.e.Q ss = −Ẇ ss . The energy absorbed by the cavity mode during the evolution is fully employed to generate the unitary displacement α, that is, ∆U = ω|α|∆X ϕ = ω|α| 2 . However, the transient dynamics ruling this process is far from trivial. The cavity mode is always displaced, i.e. gaining coherence, at a higher rate than energy,U = ω|α|Ẋ S (1 − e −γ0t/2 ), in accordance with the positive non-adiabatic entropy production rate. In addition, by comparing Figs. 7(a) and 7(b) the energetic meaning of the adiabatic entropy production rate can be clarified. In the initial transient whereṠ a < 0 the coherence gain surpass the input power, , and total (Ṡi) entropy production rates represented by color solid lines, and (b) input power (Ẇ ), rate at which the cavity mode absorbs energy (U ), and displacement rate (Ẋϕ). The cavity mode starts in equilibrium with the thermal reservoir, ρ0 = e −βĤ 0 Z, and the laser driving is switched on at t = 0. The dashed line in (a) corresponds to βẆss, and we used vertical dotted lines to highlight the instant of time at which the adiabatic entropy production rate changes its sign (tn). We used parameters = 0.02ω, γ0 = 0.01ω, and temperature kT = 10 ω, for ω = 1 units.
i.e. ω|α|Ẋ ϕ >Ẇ , which in turn implies that the rate at which the cavity mode gains energy speeds up in this periodÜ > 0. At time t n , whenṠ a = 0, we havė W = ω|α|Ẋ ϕ =Ẇ ss /2, andU peaks at its maximum. After this time, the adiabatic entropy production rate is positiveṠ a > 0, implying ω|α|Ẋ ϕ <Ẇ , andU decreases until it becomes zero in the long time run, when the periodic steady state is reached. In conclusion, we obtained that the sign of the adiabatic entropy production rate spotlights the acceleration in the internal energy changes of the cavity mode.

VII. CONCLUSIONS
In this paper we have analyzed the production of entropy in general processes embedded in a two measurement protocol, with local measurements performed in both system and environment. Our first main result is the fluctuation theorem (23), which compares the probability of forward and backward trajectories. Particularizing this expression to certain initial conditions of the backward process, one can obtain FT's for the change of inclusive (14) and exclusive (15) entropy production, i.e., for the entropy production that results from keeping or neglecting the classical correlations generated between the system and the environment during the process. These differences have been illustrated for the case of two qubits interacting through a CNOT quantum gate.
We have also explored whether it is possible to split the total entropy production into adiabatic and nonadiabatic contributions, as it is customary in classical systems far from equilibrium [22,23]. For that purpose, we have introduced a dual dynamics for the reduced evolution of the system, which in turn allowed us to clarify the interpretation of previous FT's derived for quantum CPTP maps [33]. We have shown that the aforementioned decomposition is possible only if the reduced dynamics satisfies certain condition, namely Eq. (38). In fact, we give an explicit example where that condition is not fulfilled and the adiabatic entropy production takes on negative values. This is a pure quantum feature whose consequences, we believe, are worth to be further explored.
Our results can be applied to a broad range of quantum processes including multipartite environments and concatenations of CPTP maps. In particular, we developed in detail the application to processes described by Lindblad master equations. We have introduced a general method to identify the environmental entropy change during the trajectories induced by quantum jumps [see Eq. (70) and below], which allows us to recover the FT's. The meaning of the terms adiabatic and non-adiabatic become clear in this situation, since the non-adiabatic contribution tends to zero for quasi-static driving.
We have finally studied the decomposition of the total entropy production in two specific situations of interest: an autonomous three-level thermal machine and a dissipative cavity mode resonantly driven by a classical field.
Summarizing, our results provide an exhaustive characterization of the entropy production in open quantum systems undergoing arbitrary processes. This includes: systems in contact with non-thermal or finite-size reservoirs, configurations with several equilibrium baths with different temperatures or chemical potentials, driven systems, etc. In all those cases, one should be able to assess the entropy production and characterize its fluctuations within the theoretical framework presented in this paper. Therefore, our results clarify the origin of entropy production from coarse-graining and its link to thermodynamical notions when particular choices for the environment are made.

ACKNOWLEDGMENTS
We thank Gili Bisker for comments. The authors acknowledge funding from MINECO (Proyecto TerMic, FIS2014-52486-R). G. M. acknowledges support from FPI grant no. BES-2012-054025. This work has been partially supported by COST Action MP1209 "Thermodynamics in the quantum regime". In this appendix we show that the term S(ρ * E ||ρ E ) appearing in Eq. (29) is negligible for infinitesimal changes in the environment density matrix. Let us assume the change in the environment density operator in the following general form: where Tr[∆ρ E ] = 0, and ≥ 0 is a small real number.
Using the definition of the quantum relative entropy, we can then write: In the following we show that if 1, then , and consequently S(ρ * E ||ρ E ) → 0. This can be done by applying perturbation theory. Let the eigenvalues and eigenstates of ρ * E , the set {q * µ , |φ * µ }, be expanded to second order in : where the zeroth order contributions obey ρ E |φ µ = q µ |φ µ , and we have for the first order terms: We now calculate the entropy change up to second order in : (1) µ ln q µ − 2 q The above equations (A7) and (A8) are equal up to first order, differing in O( 2 ). Therefore, using Eq. (A2), we conclude that: and when 1 we can consider S(ρ * E ||ρ E ) → 0 up to first order, in contrast to the change in entropy [Eq. (A7)].
for r = 1, ..., R. (B3) That is, the Kraus operators of the forward process are given by and analogously for the Kraus operators of the backward process (31) we havẽ The key relation (32) necessary to obtain the fluctuation theorem for the total entropy production (23) hence follows as well in this case, with a decomposition of the environment boundary term being σ (r) µ (r) ν (r) ≡ − lnq (r) µ (r) + ln q (r) ν (r) . The application of the above formalism introducing the dual and dual-reverse processes follows immediately in the same manner, leading to the fluctuation theorems for the adiabatic and non-adiabatic entropy production in detailed and integral versions, Eqs. (42), (48) and (49). The adiabatic entropy production per trajectory and its average then read in this case: ∆ i S a = R r=1 S(ρ * r ) − S(ρ r ) + ∆φ ≥ 0, (B8) where in the averaged version we set again (uncorrelated) reversible boundaries,ρ SE = Θ(ρ * S ⊗ ρ * 1 ⊗ ... ⊗ ρ * R )Θ † .

Appendix C: Concatenations of CPTP maps
In the following we focus in the derivation of FT's for concatenations of CPTP maps reported in Sec. IV E.
Here we assume the environment is a single reservoir or ancilla. However, the extension to multiple reservoirs follows in the same manner than in the one map case (see Appendix B).
Consider the maps concatenation Ω in Eq. (57). For any map E (l) in the sequence the environmental ancilla starts in a generic state and it is measured at the beginning and at the end of the interaction with the system, generating outcomes labeled as ν l and µ l , respectively. The measurements are specified by the rank-one projective operators {Q (2). Here we consider always the same total time-dependent Hamiltonian H(t), following an arbitrary driving protocol Λ = {λ t | 0 ≤ t ≤ N τ }. For convenience the latter can be also split into N intervals; hence the partial protocol Λ l = {λ t | t l−1 ≤ t ≤ t l } generates the unitary operator U (l) Λ . A quantum trajectory in this context is defined as follows. At time t = 0 we start with our system in ρ S , which is measured with eigenprojectors {P n }, obtaining outcome n. Then the sequence of maps Ω defined in Eq. (57) is applied, obtaining outcomes {µ l , ν l } from each of the l = 1, ..., N pairs of measurements in the environment. Finally at time t = N τ the system is measured again with arbitrary (rank-one) projectors {P * m } giving outcome m.
Dual and dual-reverse maps and operations also follow from its definitions in Sec. IV when conditions (38) and E (l) (π (l) ) =π (l) are met for each map in the sequence. The corresponding probabilities for trajectory γ in the dual process, and trajectoryγ in the dual-reverse are: where in the dual-reverse trajectories we took again the sequence of maps in inverted order, that is, we applied D where the nonequilibrium potential changes are defined with respect to the invariant state π (l) of each map E (l) as in the single map case: ∆φ ν . The set of equations (C7-C9) immediately implies the detailed FT's for concatenations in Eqs. (58)(59)(60). Its corresponding integral versions and second-law-like inequalities follow immediately as a corollary.
Finally, it is interesting to consider the expression of the average nonequilibrium potential change during the whole sequence. By denoting ρ S (t l ) the reduced state of the system at time t l , we have: Tr[E (l) µ l ν l (ρ S (t l−1 ))]∆φ (l) Tr[(ρ S (t l ) − ρ S (t l−1 ))Φ l ], where Φ l = − ln π (l) . The above expression can be decomposed into the following boundary and path contributions: When all the maps in the concatenation have the same invariant state, Φ l+1 = Φ l ≡ Φ ∀l, we obtain ∆Φ p = 0, while ∆Φ b = Tr[(ρ S − ρ S )Φ] and we recover the expression for the single map case, c.f. Eq. (53). In the other hand the boundary term only vanishes for cyclic processes, such that ρ S = ρ S , implemented by cyclic concatenations with Φ N = Φ 1 . In this case ∆Φ b = 0 while ∆Φ p gives in general a non-zero contribution. The dynamical versions of these boundary and path terms read: which are also analogous to their classical counterparts [21][22][23].