Joint Fluctuation Theorems for Sequential Heat Exchange

We study the statistics of heat exchange of a quantum system that collides sequentially with an arbitrary number of ancillas. This can describe, for instance, an accelerated particle going through a bubble chamber. Unlike other approaches in the literature, our focus is on the joint probability distribution that heat Q1 is exchanged with ancilla 1, heat Q2 is exchanged with ancilla 2, and so on. This allows us to address questions concerning the correlations between the collisional events. For instance, if in a given realization a large amount of heat is exchanged with the first ancilla, then there is a natural tendency for the second exchange to be smaller. The joint distribution is found to satisfy a Fluctuation theorem of the Jarzynski–Wójcik type. Rather surprisingly, this fluctuation theorem links the statistics of multiple collisions with that of independent single collisions, even though the heat exchanges are statistically correlated.


Introduction
Fluctuations of thermodynamic quantities, which are usually negligible in macroscopic systems, are known to play a dominant role in the micro-and mesoscopic domain. These fluctuations are embodied in the so-called fluctuation theorems (FT) [1][2][3][4], a collection of predictions for systems evolving under nonequilibrium conditions valid beyond linear response. They can be summarized as [5,6] where P(Σ) denotes the probability that an amount of entropy Σ is produced in a certain process and P(Σ) denotes the corresponding probability for the time-reversed process.
Of the many scenarios which present FTs, one which is particularly interesting is that of heat exchange between a system S, prepared in equilibrium with a temperature T s , and an environment E, prepared in a different temperature T e . In this case, as first shown by Jarzynski and Wójcik in Ref. [7], the distribution P(Q) of the heat exchanged between them, satisfies P(+Q) P(−Q) = e ∆βQ , where ∆β = β e − β s (with β = 1/T and k B = 1). Here, and throughout the paper, Q denotes the net heat transfer from the system to the environment. Quite surprising, in this case it turns out that P(Q) = P(Q), meaning the statistics of the forward and backward processes are the same. Equation (2) was subsequently generalized to allow for the exchange of both energy and particles between several interacting systems initially at different temperatures and chemical potentials [6,8,9].
Here we consider a generalization of this scenario, where the system interacts sequentially with multiple parts of the environment, exchanging heat with each part. One can imagine, for instance, an accelerated particle crossing a bubble chamber. In this case, the system will leave a trail on E, represented by the heat exchanged in each point. In the microscopic domain this process will be stochastic, with a random amount of heat exchanged in each interaction.
The key idea that we will explore in this paper is to look at the joint probability distribution for the heat exchanged with each part, P(Q 1 , Q 2 , Q 3 , . . .). This allows us to understand the correlations between the different heat exchanges.
For instance, in a situation where all the ancillas have the same temperature, from a stochastic perspective a large exchange in the first collision increases the probability that the second collision exchanges less. This feature is fully captured by the joint distribution. This happens because thermal operations have the property of bringing the system closer to its thermal equilibrium state, σ eq , i.e., [10] where D(ρ ||ρ) = Tr(ρ ln ρ − ρ ln ρ) is the quantum relative entropy. If in the first interaction the system exchange a large quantity of heat, the system gets a lot closer to its steady state. So in the next interaction, the system should exchange less heat.
To formalize this idea, we split the environment into a set of ancillas A i , with which the system interacts sequentially, producing a collisional model [11][12][13][14]. The process is schematically illustrated in Figure 1 and the formal framework is developed in Section 2. In Section 3 we then show that P(Q 1 , Q 2 , Q 3 , . . .) satisfies a fluctuation theorem that generalizes (2). Moreover, we show how this fluctuation theorem relates the joint distribution to the statistics of a single collision, even though the events are statistically correlated.
Forward processes Figure 1. Schematic representation of a system S interacting sequentially with a series of ancillas. The system starts in the state σ 0 and the ancillas in an initial states ρ i , which are assumed to be thermal but at possibly different temperatures. Each SA i interaction is also governed by a possibly different unitary U i .

Formal Framework
We consider a quantum system S, with Hamiltonian H s , prepared in a thermal state σ 0 = e −β s H s /Z s , with temperature T s . The system is put to interact sequentially with a series of N ancillas A i , as depicted in Figure 1. The ancillas are not necessarily identical. Each has Hamiltonian H i and is prepared in a thermal state ρ i = e −β i H i /Z i , with possibly different temperatures T i . Each collision is described by a unitary operator U i acting only on SA i , which may also differ from one interaction to another.
In order to comply with the scenario of Ref. [7], we assume that the U i satisfy the strong energy-preservation condition Or, what is equivalent, that each collision is a thermal operation [10,15]. This implies that all energy that leaves S enters A i , so nothing is stuck in the interaction. As a consequence, there is no work involved and all the change in energy of the system can be unambiguously identified as heat flowing to the ancillas [13].
We label the eigenvalues and eigenvectors of the system as H s |α = E s α |α . For concreteness, we assume these levels are non-degenerate. Time is labeled discretely by i = 1, 2, 3, . . ., representing which collisions already took place. For instance, the initial state is decomposed as σ 0 = ∑ α 0 p 0 (α 0 )|α 0 α 0 |, with p 0 (α 0 ) = e −β s E s α 0 /Z s and we use α 0 to emphasize that this is before the first collision. Similarly, the eigenvalues and eigenvectors of the ancillas are labeled as H i |n i = E e n i |n i . The initial state of each A i is thus decomposed as ρ i = ∑ n i q i (n i )|n i n i | where q i (n i ) = e −β i E e n i /Z i . The dynamics depicted in Figure 1 generates a stroboscopic map for the system. The joint state of SA i after the interaction is given by Taking the partial trace over A i then leads to the updated state σ i . Conversely, tracing over the system leads to the reduced state ρ i of the ancilla after the interaction (Figure 1).
The fact that the unitary is energy preserving (Equation (4)), together with the assumption that the energy levels are non-degenerate, means that it is possible to construct quantum trajectories for the system in two equivalent ways. The first is to assume a two-point measurement scheme in S at each step [16,17]. Equation (4) implies that the system will remain diagonal in the energy basis, so that measurements in this basis are non-invasive (that is, have no additional entropy production associated to it). Measuring S in the energy basis after each collision then leads to the trajectory The heat associated with each collision is then readily defined as Alternatively, one can construct a quantum trajectory by measuring the ancillas, before and after each collision, plus a single measurement of the system before the process starts. That is, one can consider instead a quantum trajectory of the form γ e = {α 0 , n 1 , n 1 , n 2 , n 2 , . . . , n N , n N }. (8) This, in a sense, is much more natural since the ancillas are only used once and thus may be experimentally more easily accessible. Furthermore, as far as heat exchange is concerned, this turns out to be equivalent to the trajectory (6). The reason is that Equation (4) implies the restriction where δ(x) is the Kronecker delta. In addition, since the energy values are taken to be non-degenerate, energies uniquely label states. Thus, for instance, if we know α 0 , n 1 , n 1 we can uniquely determine α 1 , and so on. The converse, however, is not true: from α 0 and α 1 we cannot specify n 1 and n 1 (which is somewhat evident given that the number of points in Equation (6) is smaller than that in Equation (8)). This, however, is not a problem if one is interested only in the heat exchanged, which can also be defined from the trajectory (8) as Due to Equation (9), this must coincide with Equation (7); i.e., The assumption in Equation (4) may at first seem somewhat artificial. However, this is not the case. This assumption is a way to bypass the idea of weak coupling, which is one of the conditions used in [7]. Moreover, the interesting thing about the present analysis is that it establishes under which conditions Equations (6) and (8) are equivalent. Naively one would expect that this is often the case. However, as the above arguments show, several assumptions are necessary for this to be the case. This reflects some of the challenges that appear in describing thermodynamics in the quantum regime.

Path Probabilities from Measurements in S
Thermal operations imply that the probability that, after the i-th collision, the system is in a given eigenstate |α i depends only on the probabilities in the previous time. That is, the dynamics of populations and coherences completely decouple [18]. Indeed, Equation (5) together with Equation (4) imply that where The populations therefore evolve as a classical Markov chain, with M i (α i |α i−1 ) representing the transition probability of going from α i−1 to α i . Moreover, Equation (9) together with the fact that the ancillas are initially thermal, imply that M i (α i |α i−1 ) satisfies detailed balance where, notice, what appears here is the temperature β i of ancilla A i . The path probability associated with γ s in Equation (6) will then be which is nothing but the joint distribution of a Markov chain. We call attention to the clear causal structure of this expression: marginalizing over future events has no influence on past ones. For instance, summing over α N leads to a distribution of the exact same form. Conversely, marginalizing over past variables completely changes the distribution. The joint distribution of heat can then be constructed from Equation (14) in the usual way: This is the basic object that we will explore in this paper.
Instead, from P [γ se ] one can now obtain P [γ e ] by marginalizing over α 1 , . . . , α N ; viz., The above analysis puts in evidence the Hidden Markov nature of the dynamics in Figure 1. When measurements are done in the ancilla, the system plays the role of the hidden layer, which is not directly accessible. Instead, predictions about the system must be made from the visible layer (i.e., the ancillas).
This Hidden Markov nature manifests itself on the fact that even though the system obeys a Markov chain [Equation (14)], the same is not true for the ancillas. In symbols, this is manifested by the fact that n i depends not only on n i and n i−1 , but on the entire past history (n 1 , n 1 , . . . , n i−1 , n i−1 , n i ). This is intuitive in a certain sense: the amount of heat exchanged at the i-th collision will depend on the heat exchanged in all past events.
With P [γ e ], the distribution of heat, Equation (15) can be equivalently defined using Equation (10). One then finds The reason why this is equivalent to Equation (15) becomes clear from the way we derived P [γ e ] above: we can expand the summation to γ se and then use the fact that Q i [γ s ] = Q i [γ e ].

Backward Process
To construct the fluctuation theorem, we must now establish the backward process. As shown in [19], however, there is an arbitrariness in the choice of the initial state of the backward process; different choices lead to different definitions of the entropy production. Here we are interested specifically in heat and the generalization of the Jarzynski-Wójcik fluctuation theorem [7]. Hence, we assume that in the backward process both system and ancillas are fully reset back to their thermal states. As usual, the time-reversed interaction between SA i now takes place by means of the unitary U † i . However, the order of the interactions must now be flipped around, as shown in Figure 2. More about the choice of backward process can be found in [20,21] and its relation to the notion of recovery maps is discussed in [22].  In the backward process, the system will therefore evolve according tõ where we index the states as α N−i instead of α i just so that the trajectory γ s can remain the same as in the forward process. The path probabilityP [γ s ] associated to this process will then bẽ which is similar to that used in the original Crooks fluctuation theorem [23]. The corresponding heat distribution isP where Q i continues to be the heat exchanged with A i (which is now different from the heat exchanged at collision i).

Joint Fluctuation Theorem for Heat Exchange
We are now ready to construct the fluctuation theorem. The detailed balance condition (13) immediately implies that Equations (15) and (20) will be related by This is a theorem for the joint distribution of the heat exchanged between multiple ancillas. It thus represents a generalization of Ref. [7] to the case where the system interacts sequentially with multiple reservoirs. This result has several features which are noteworthy. First, note that the temperature β i of the ancillas are not necessarily the same. Second, note how after the first collision the state of the system is no longer thermal. However, still, this does not affect the fluctuation theorem. All that matters is that before the first collision the system is in equilibrium.
It is also important to point out that any Markov chain satisfying the detailed balance relation also satisfies a fluctuation theorem [24]. This fact can be used to obtain Equation (21) when properly choosing the rates of the Markovian evolution. Beyond that, a generalization of the detailed FT to multiple reservoirs has also being obtained before, e.g., in Ref. [25].

Causal Order and Relation to Single Collisions
The causal order of the process plays a crucial role here. Marginalizing over future events has no effect on the fluctuation theorem. That is, from (21) one could very well construct a similar relation for P(Q 1 , . . . , Q N−1 ), by simply summing over Q N . This is not possible, however, for marginalization over past events. That is, P(Q 2 , . . . , Q N ), for instance, does not satisfy a fluctuation theorem.
The right-hand side of Equation (21) is very similar to what appears in the original FT (2). We can make this more rigorous as follows. Let us consider a different process, consisting of a single collision between the system thermalized in β s and an ancilla thermalized in β i (Figure 3). The associated heat distribution P sc (Q i ) will then satisfy Equation (2); viz., where, recall, in this case of a single collision the backward process coincides with the forward one, so that the distributionP sc in the denominator is simply P sc . It is very important to emphasize, however, that P sc (Q i ) is not the marginal of P(Q 1 , . . . , Q N ) (with the exception of Q 1 ). Notwithstanding, comparing with Equation (21), we see that the full process in Figure 1 is related to the single-collision processes according to This result is noteworthy, for the right-hand side is a product whereas the left-hand side is not. The full distribution P(Q 1 , . . . , Q N ) cannot be expressed as a product because the heat exchanges are, in general, not statistically independent. Notwithstanding, the ratio on the left-hand side of (23) does factor into a product. The point, though, is that this is not the product of the marginals, but of another distribution P sc .
One can also write a formula of the form (23), but for only some of the heat exchanges. For instance, it is true that This kind of decomposition, however, depends crucially on the causal structure since it can only be done for future exchanges. For instance, we cannot write something involving P(Q 2 , . . . , Q N ). The reason is that P(Q 1 , . . . , Q N−1 ) satisfies the fluctuation theorem (21), but P(Q 2 , . . . , Q N ) does not (since, after the first collision the system is no longer in a thermal state).

Information-Theoretic Formulation of the Entropy Production
We define the entropy production associated with Equation (21) as The second equality is obtained using the detailed balance relation (13). We emphasize that this is the entropy production associated with the choice of backward protocol used in Section 2.3, which may differ from other definitions in the literature [18,26]. As discussed in [19], the interpretation of the entropy production depends on the choice of the initial state of the backwards process. For instance, if we have chosen the initial state as the final state of the forward process, i.e., the state N (see Equation (5)), we would have a a contribution related to the correlations between the system and the ancillas. This type of entropy production was called the inclusive entropy production in Ref. [19]. This happens because this state carries the information about the correlations. Here we have choose a initial state for the backward process that does not have this contributions.
In [7], Jarzynski and Wójcik calculated an upper bound on the probability of observing a violation of the second law, i.e., the passage of heat from a colder to a hotter body. We can apply the same reasoning to Equation (25). Let us assume that all ancillas start in the same thermal state with temperature T a and β a − β s > 0. The probability that the heat transfer from the system to i-th ancilla will fall below a specified value q i in each interaction through the whole process, obeys the inequality which is the multiple-exchange extension of the result obtained in [7]. Equation (26) shows that observing a positive total transference of heat from the hot system to the cold ancillas dies exponentially with q 1 + · · · + q N . Alternatively, we can consider the entropy production from the perspective of the global trajectory γ se in Equation (16). Using also that Q i [γ s ] = Q i [γ e ], we can then write Σ[γ se ] as The average entropy production may then be written as where S(ρ) = −Tr(ρ ln ρ) is the von Neumann entropy. Here σ N is the final state of the system after the N collisions. In the Equation (28), we can identify where ∆I se is the change in the mutual information between the system and the ancillas. This way we can have a more clear meaning of the expression (28). One term is proportional to the total correlations built between system and ancillas and the other two relative entropy terms measure the disturbance on the environment and the system during the process. The important aspect of this result is that it depends only on local changes in the ancillas. That is, all quantities refer to the local states ρ i of each ancilla after the interaction. In reality, because the ancillas all interact with the system, they actually become indirectly correlated. These correlations are still represented indirectly in Σ[γ se ], but they do not appear explicitly. This, ultimately, is a consequence of the choice of backward process that is used in the Jarzynski-Wójcik scenario [7].

Initially Correlated Ancillas
One possible extension of our formalism is to consider the case of initially correlated system-ancillas. In this case, we could explore how the correlation between the system and the ancillas affect the XFT. This problem was studied for a single heat exchange in [27] and in our case, the same approach yields where the ∆I (γ se ) = I * − I with I * = ln p(α n , n 1 , . . . , n N ) p 0 (α N )q 1 (n 1 ) . . . q N (n N ) (31) where we define p(α 0 , n 1 , . . . , n N ) = α 0 , n 1 , . . . , n N |ρ SE |α 0 , n 1 , . . . , n N . Here ρ SE is the initial state for the system-ancillas. This result is similar to the one found in [27]. Because in our case, we are working with thermal operations, we can write Equation (30) as By taking the above equation and sum over all trajectories, to obtain the nonequilibrium equality for an initially correlated state and then using Jensen's inequality we have that So it is possible to obtain a type of Clausius relation where now the entropy production has a new lower bound.

Conclusions
To summarize, we have considered here the sequential heat exchange between a system and a series of ancillas. We assume all entities start in thermal state, but at possibly different temperatures. Moreover, all interactions are assumed to be described by thermal operations, which makes the identification of heat unambiguous. The main object of our study was the joint probability of heat exchange P(Q 1 , . . . , Q N ) for a set of N collisions. This object contemplates the correlations between heat exchange, a concept which to the best of our knowledge, has not been explored in the quantum thermodynamics community. We showed that P(Q 1 , . . . , Q N ) satisfies a fluctuation theorem, which relates this joint distribution with single collision events. This result, we believe, could serve to highlight the interesting prospect of analyzing thermodynamic quantities in time-series and other sequential models.
Author Contributions: J.S., A.T. and G.L. equally contributed to conceptualization, investigation and writing of the paper. All authors have read and agreed to the published version of the manuscript.