# Heat fluctuations and fluctuation theorems in the case of multiple reservoirs

###### Abstract

We consider heat fluctuations and fluctuation theorems for systems driven by multiple reservoirs. We establish a fundamental symmetry obeyed by the joint probability distribution for the heat transfers and system coordinates. The symmetry leads to a generalisation of the asymptotic fluctuation theorem for large deviations at large times. As a result the presence of multiple reservoirs influence the tails in the heat distribution. The symmetry, moreover, allows for a simple derivation of a recent exact fluctuation theorem valid at all times. Including a time dependent work protocol we also present a derivation of the integral fluctuation theorem.

###### pacs:

05.40.-a,05.70.Ln## I Introduction

There is a current interest in the thermodynamics and statistical mechanics of fluctuating systems in contact with heat reservoirs and driven by external forces. The current focus stems from the recent possibility of direct manipulation of nano systems and biomolecules. These techniques permit direct experimental access to the probability distribution functions for the work or for the heat exchanged with the environment Trepagnier et al. (2004); Collin et al. (2005); Tietz et al. (2006); Blickle et al. (2006); Wang et al. (2002); Imparato et al. (2007); Imparato et al. (2008a); Douarche et al. (2006); Garnier and Ciliberto (2007); Imparato et al. (2008b); Ciliberto et al. (2013, 2014); Koski et al. (2013). These methods have also yielded access to the experimental verification of the recent fluctuation theorems which relate the probability of observing entropy-generated trajectories with that of observing entropy-consuming trajectories Jarzynski (1997); Kurchan (1998); Gallavotti (1996); Crooks (1999, 2000); Seifert (2005a, b); Evans et al. (1993); Evans and Searles (1994); Gallavotti and Cohen (1995); Lebowitz and Spohn (1999); Gaspard (2004); Imparato and Peliti (2006); van Zon and Cohen (2003a); van Zon et al. (2004); van Zon and Cohen (2003b, 2004); Speck and Seifert (2005); Imparato et al. (2007); Esposito and den Broeck (2010).

In the present paper we address the issue of heat fluctuations and fluctuation theorems in the case of a system coupled to multiple reservoirs. The case of a system driven by two heat reservoirs, both in the case of one degree of freedom and in the case of many degrees of freedom, has been discussed extensively, see e.g. Fogedby and Imparato (2011, 2012), whereas the case of multiple reservoirs seems to have received less attention. In the case where a system is coupled to many reservoirs the heat flows exhibit a more complicated pattern which will influence the heat distribution.

We characterise the heat transfers to the system by the vector quantity with components , referring to the n-th heat reservoir maintained at temperature ; we set Boltzmann’s constant . We, moreover, assume that the reservoirs are independent. Noting that the transfer of heat also influences the internal state of the system characterised by the coordinate , the central quantity is the joint distribution . For simplicity we consider a single degree of freedom ; the generalisation to many degrees of freedom is straightforward and discussed at the end of the paper.

Assuming that the system is initially in the state at time and that the sampled heat vanishes at , the joint distribution describes the transition of the system from the initial state at time to the final state at time in contact with multiple heat reservoirs transferring the heat to the system during the time span . We, moreover, assume that the system moves in the time independent potential . Hence, there are no external forces acting on the system and the non equilibrium state is entirely driven by the heat transfers from the reservoirs.

It follows from the structure of Fokker-Planck equation that the joint distribution obeys the fundamental symmetry

(1) |

as discussed in Sec. II.2. Here we have singled out the k-th reservoir at temperature and, moreover, introduced the notation

(2) |

The symmetry (1) which is valid at all times is quite general for systems driven stochastically by heat reservoirs. The symmetry basically establishes a connection between the heat transfers during a transition from to and minus the heat transfers during the reverse transition from to . The symmetry incorporates the general features of a reservoir-driven non equilibrium transition. As will become clear, the symmetry implies both the asymptotic fluctuation theorem for the cumulant generating function valid at long times, see e.g. Gallavotti and Cohen (1995); Lebowitz and Spohn (1999), the integral fluctuation theorem Seifert (2005a), and the exact fluctuation theorem proposed in Cuetara et al. (2014) valid at all times.

In the present paper we discuss in detail the implications of the symmetry (1). For the purpose of a simple discussion we introduce in Sec. II the model of a single over damped particle moving in a time independent potential and at the same time driven by multiple heat reservoirs. This model allows for a simple derivation of the symmetry in (1). In Sec. III we consider the asymptotic long time regime and discuss the application of the symmetry to the cumulant generating function. Limiting the discussion to the harmonic potential we consider the branch cut structure of the cumulant generating function and the implications for the tails in the associated heat distribution. In Sec. IV we discuss implications of the symmetry valid at all times. We make contact with the trajectory formulation and associated entropy considerations. We also present a derivations of the integral fluctuation theorem and the exact fluctuation theorem. In Appendix E we summarise the extension of the formalism to the case of a time dependent potential modeling an external work protocol. This extension allows us to demonstrate the integral fluctuation theorem in the general case. In Sec. V we consider the case of a deterministic Hamiltonian system driven by heat reservoirs and present the corresponding symmetry for the joint distribution. In Sec. VI we present a summery and a conclusion. Technical matters are deferred to a series of appendices.

## Ii System coupled to multiple heat reservoirs - one degree of freedom

### ii.1 General

Here we consider for simplicity a single degree of freedom moving in the potential and at the same time driven by multiple heat reservoirs maintained at temperatures . We restrict our discussion to the autonomous case where the potential is time independent. Associating the damping with the n-th reservoir, the corresponding Langevin equation has the form, see e.g. Imparato et al. (2007); Seifert (2005a),

(3) |

where we have introduced the total damping, the total noise, and the noise correlations, i.e.,

(4) | |||

(5) | |||

(6) |

here a prime denotes the spatial derivative . Correspondingly, the heat flux from the n-th reservoir is given by , where is the force on the system originating from the n-th reservoir, i.e.,

(7) |

The configuration is depicted in Fig. 1.

Since only a single degree of freedom is coupled to multiple reservoirs, the system is maintained in a state described by a Boltzmann-like probability distribution with an effective temperature. This is in contrast to for example a harmonic chain coupled to two reservoirs, where a genuine non equilibrium situation is established, see e.g. Fogedby and Imparato (2012). By inspection of the Langevin equation (3) for we infer the effective equilibrium temperature

(8) |

and the stationary distribution . However, the individual heat fluxes between the reservoirs via the particle constitute a non equilibrium problem. From the equations of motion (3) and (7) we also infer

(9) |

expressing global energy conservation; here denotes the initial configuration. We assume in the following that , i.e., we start sampling the heat at . In other words, whereas the individual heat component is fluctuating, the sum locks onto the nonfluctuating potential difference .

The stochastic equations (3) and (7) define the problem we wish to study. However, for the present purposes it is more convenient to adhere to a Fokker-Planck description Risken (1989); van Kampen (1992). Since the transfer of heat changes the state of the system and thus couples to the coordinate , the distribution for the coordinate and heat is characterised by the joint distribution satisfying a Fokker-Planck equation with Liouville operator van Kampen (1992); Imparato et al. (2007), i.e.,

(10) | |||

(11) |

Since is linear in it is convenient to introduce the characteristic function defined according to

(12) |

where . then satisfy the Fokker-Planck equation with Liouville operator below

(13) | |||

(14) |

We note that setting corresponds to integrating over and we have , where is the distribution function for the coordinate in global equilibrium with the reservoirs at temperature given by (8). thus satisfies the Fokker-Planck equation , where , with stationary distribution .

### ii.2 Symmetry

Incorporating the initial condition at time in the notation by setting , imposing the boundary condition , and introducing the matrix notation

(15) |

the Fokker-Planck equation (13) takes the form

(16) |

In analogy with corresponding manipulations for the evolution operator in quantum mechanics Zinn-Justin (1989), we obtain the formal solution

(17) |

here the the n-th term in the expansion of the exponential is defined according to .

Choosing a particular reservoir at temperature inspection of the Liouville operator (14) in the form (15) yields the symmetry

(18) |

for details see Appendix B. From the solution (17) we then readily infer the corresponding symmetry for , i.e.,

(19) |

which holds for any choice of . From the inverse Laplace transform (12) Lebedev (1972),

(20) |

we finally obtain (1), i.e.,

(21) |

that again holds for all . In the special case we obtain in particular the simpler symmetries

(22) | |||

(23) |

This concludes the demonstration of the symmetries for the characteristic function and the joint heat distribution. In the next section we consider the implications for the cumulant generating function.

## Iii Cumulant generating function

### iii.1 General

At long time the characteristic function is dominated by the largest eigenvalue of the Liouville operator . This follows heuristically from the solution (17) and can be made more precise by considering the spectral representation of ; this discussion is deferred to Appendix A. We have at long times compared to

(24) |

where is the cumulant generating function. The prefactor is expressed in terms of the eigenfunctions associated with the eigenvalue , for details see Appendix A. From the definition (12) we also have, denoting the constrained average by the subscript ,

(25) |

and it follows that yields the mean heat and higher cumulants, i.e. , , etc.

Applying the symmetry (22) to (24) we readily obtain the usual asymptotic fluctuation theorem Gallavotti and Cohen (1995); Lebowitz and Spohn (1999) generalised to multiple reservoirs

(26) |

Likewise, we infer, applying the full symmetry (19) to (24), the generalised fluctuation theorem , for fixed . Here we note that since the exponential in (19) is subdominant in time compared with at long times. Moreover, since this symmetry holds for any choice of we conclude that only depends on the difference variable

(27) |

and we have

(28) |

This relation generalises the usual asymptotic fluctuation theorem to the case of multiple reservoirs.

From (20) the heat distribution at long times is given by

(29) |

Since the integral at long times is controlled by the exponential and depends on the difference , it follows that the integral is invariant under the shift , where is arbitrary. Consequently, and we obtain the constraint

(30) |

note that this constraint on the heat variables in the heat distribution at long times is not in conflict with (9), i.e. which expresses energy conservation for the fluctuating heat. At long times the dependence of on is confined to the sub manifold determined by (30).

The condition (30) also follows from the symmetry (1) noting that is subdominant in time and that can be chosen arbitrarily. We have in abbreviated notation

(31) |

In other words, for reservoirs there are only independent heat transfers as a result of the constraint in (30). For two reservoirs we have for example , where is the heat distribution for reservoir 1 and we infer that , i.e., the distribution for the outgoing heat from reservoir 1 is identical to the distribution for incoming heat to reservoir 2, expressing conservation of energy.

By the usual steepest descent argument Touchette (2009); Lebowitz and Spohn (1999); Fogedby and Imparato (2012) the prefactor locks onto , where is a solution of and we have the scaling form

(32) | |||

(33) |

where , , is the large deviation function in the case of multiple reservoirs. The symmetry (1) applied to then yields the fluctuation theorem for the large deviation function

(34) |

### iii.2 Branch points and tails

It follows from general principles Touchette (2009); den Hollander (2000) that the cumulant generating function is bounded and possesses a branch cut structure. For large the Laplace transform in (29), when closing the contour in its evaluation, samples the edges of the branch cuts closest to the origin determining the tails of the heat distribution. We note that in principle the singular structure of the prefactor can also influence the tails. This issue has been discussed in detail for the case of a particle driven by two heat reservoirs Visco (2006); Farago (2002). In the present context we only analyse the contribution arising from the singular structure of . In general is a downward convex function passing through the origin due to normalzation. From the symmetry (28) it then also follows that . Moreover, if has a branch point at there will also be a branch point present at .

In order to be more specific we return to the model in Sec. II and note that the similarity transformation Imparato et al. (2007), where

(35) |

maps the Liouville operator in (14) to a Hermitian Schrödinger form, i.e., , where

(36) |

Here

(37) | |||

(38) |

Consulting Appendix C it follows that and have identical spectra, in particular the same largest eigenvalue . From the structure of it, moreover, follows that depends on the form of and parametrically on the temperature dependent parameter . We also note that the function is invariant under the transformation , in accordance with (28)

In order to obtain an explicit expression for we consider in the following a harmonic oscillator potential . Since and the operator (36) corresponds to a quantum oscillator with mass moving in the potential . Referring to Appendix D the discrete spectrum is given by , , and we obtain in particular the cumulant generating function . In more detail inserting (37)

(39) |

#### iii.2.1 One reservoir

The situation is simple in the case of a single reservoir. We have , i.e., and from (37) . locks onto and we have for all , i.e., a vanishing cumulant generating function. This is consistent with the observation that a single degree of freedom coupled to a single heat reservoir is maintained in equilibrium at temperature with a bounded heat distribution Imparato et al. (2007); Fogedby and Imparato (2011).

#### iii.2.2 Two reservoirs

The case of two reservoirs with dampings and and temperatures and was considered by Derrida et al. Derrida and Brunet (2005), see also Fogedby and Imparato (2011); Visco (2006). Here we have , , , and the cumulant generating function is given by the explicit expression

(40) |

The cumulant generating function is a downward convex function passing through and . It, moreover, obeys the symmetry . Closing off reservoir 2 by setting we have since a single particle driven by one reservoir is in equilibrium. The branch points are given by Fogedby and Imparato (2011)

(41) |

In the equal temperature case for , i.e., we have in particular the branch points . For the branch points are independent of the damping and located at

(42) |

#### iii.2.3 Multiple reservoirs

For multiple reservoirs the analysis is based on the general expression in (39) originating from the harmonic potential case. We, moreover, for simplicity consider the case of reservoirs with identical damping constants, i.e., and thus . We find

(43) |

and the branch point condition reads . This is still a complex expression to analyse; however, focussing on reservoir 1 by setting and , for we have for and we obtain the branch points

(44) |

For we recover (42); for we have . As the number of reservoirs increase the branch points recede to infinity as .

### iii.3 Heat distribution

For small and assuming small we can replace by a parabolic approximation consistent with the symmetry (28) and the boundary condition i.e.,

(45) |

Here is symmetric and related to the mean heat flux according to . From the steepest descent calculation we then obtain after some algebra the following expression for the large deviation function

(46) |

where

(47) | |||

(48) |

Note that is quadratic in yielding a Gaussian heat distribution for small ; for the corresponding analysis in the case of two reservoirs, see Fogedby and Imparato (2012).

## Iv Fluctuation theorems

### iv.1 Exact fluctuation theorem

Here we discuss further consequences of the symmetry (1) valid at all times. Multiplying both sides with , integrating over and , and exchanging and in the integral on the left hand side, we obtain the identity

(51) |

where we have defined the average

(52) |

Here is the joint heat distribution associated with a non equilibrium transition from at time to at time with the system in equilibrium with the -th reservoir at the initial time . The expression (51) is a fluctuation theorem valid at all times and is a particular case of a more general expression first derived by Cuetara et al. Cuetara et al. (2014) for a system in contact with several energy and particle reservoirs, on the basis of entropic considerations at the level of single trajectories.

In the case of two reservoirs setting we obtain in particular for the joint heat distribution or integrating over the fluctuation theorem

(53) |

The interpretation is straightforward, see Cuetara et al. (2014). Preparing the system in equilibrium at time with the reservoir at temperature and subsequently monitoring the heat transferred from the reservoir at temperature , the heat distribution obeys the fluctuation theorem in (53) valid at all times.

In the case of multiple reservoirs we prepare the system in equilibrium with the k-th reservoir at temperature at time and the fluctuation theorem (51) applies to the joint heat distribution for the other reservoirs. Note that whereas the asymptotic fluctuation theorem (28) for the cumulant generating function is generic in the sense that it is independent of the potential , the fluctuation theorem in the present context requires that the system is in equilibrium with distribution with one of the reservoirs at the initial time.

### iv.2 Trajectory interpretation

#### iv.2.1 Time independent potential - no work protocol

In order to establish contact with the trajectory point of view pursued by Seifert Seifert (2005a), based on the stochastic thermodynamics scheme by Sekimoto Sekimoto (1998), and the role of entropy changes in the course of the non equilibrium time evolution, we express the definition of the characteristic function in the form

(54) |

Here is a solution of the Langevin equation (3) for a specific noise realisation defining a forward trajectory in configuration space from the initial configuration at time to the final configuration at time . Likewise, solving for , i.e., , we identity a backward trajectory from at time to at time . The symmetry (22) then reads

(55) |

For a particular noise realisation defining a trajectory in configuration space we obtain

(56) |

Note that is a solution of the Langevin equation (7). Since defines a forward trajectory the heat transfer is associated with this trajectory; likewise, is associated with a backward trajectory. For the special choice we obtain in particular

(57) |

This relationship is completely equivalent to the analysis in Seifert (2005a) based on the path integral formulation of the solution to the Fokker-Planck for the distribution .

Setting in (55) and integrating over the initial position we obtain

(58) |

where we have used the normalisation condition or, equivalently, , i.e., the conservation of probability.

Introducing the total increase of the heat bath entropy

(59) |

we can also express (58) in the form of a fluctuation theorem for the bath entropy

(60) |

This fluctuation theorem states that integrated over the initial configuration with weight one, the constrained average of along forward trajectories equals unity. The fluctuation theorem relates to the total reservoir entropy production in the course of a non equilibrium transition. The fluctuation theorem is a consequence of a basic symmetry of the joint heat-position distribution combined with the normalisation of the distribution for the position. We also note that the joint distribution incorporates on the Fokker-Planck level the stochastic thermodynamics scheme developed within a Langevin formulation Sekimoto (1998); Seifert (2005a).

#### iv.2.2 Time dependent potential - with work protocol

In the case of a work protocol characterised by a time dependent potential modelling work performed on the system, the symmetries in (19) and (22) are not available since the Liouville operator acquires an explicit time dependence. The extension of the scheme to the time dependent case is summarised in Appendix E. Whereas the fundamental symmetry still applies to the Liouville operator, i.e.,

(61) |

the characteristic function is given by the time ordered product Zinn-Justin (1989),

(62) |

Considering for example the second order term in the expansion we obtain applying the symmetry

(63) |

and we note that the time dependent protocol as exemplified by the time dependence of becomes entirely entangled in the time integrations. On the contrary, in the time independent case the factors in the expansion of the exponential in (62) cancel and we recover the symmetry in (19).

We also note in passing that in the case of a moving harmonic potential of the form , as discussed in Imparato et al. (2007), the time ordered expression reduces to integration over Gaussians (yielding error functions) and can possible be reduced. However, we abstain from such an analysis in the present context; we refer to Imparato et al. (2007) for a steepest descent analysis.

With respect to the bath-fluctuation theorem (60) we note that setting we have in operator form, see Appendix E,

(64) |

and in particular setting

(65) |

It follows from the form of the Liouville operator in (115) that , expressing conservation of probability in the Fokker-Planck equation for , . Correspondingly, for the transposed Liouvillian we have and it follows from (65) that also in the time dependent case do we obtain (58), yielding the fluctuation theorem (60).

The above analysis might appear cumbersome and possibly superfluous. General physical arguments imply that the fluctuation theorem for the bath entropy only monitors the transfer of heat and thus does not depend on the applied work protocol.

### iv.3 Integral fluctuation theorem

Seifert has proposed an integral fluctuation theorem valid at all times based on entropy production and consumption on the trajectory level Seifert (2005a). Within the present context the integral fluctuation theorem follows readily from (60) and is therefore a consequence of the basic symmetry of the Liouville operator.

The bath fluctuation theorem was formulated as an integral condition for the constrained average of the fluctuating heat or entropy for an initial configuration at at time and a final configuration at time . In order to make contact with the integral fluctuation theorem one must consider fluctuating initial and final configuration characterised by the normalised distributions and .

Based on the path integral formulation of the distribution it was shown in Seifert (2005a) that one can define entropy production and consumption for an individual trajectory from at time to at time based on the Gibbs expression . By inspection of the entropy associated with the forward and backward trajectories one easily extracts the entropy associated with the contact with the heat baths, . The additional entropy associated with the non equilibrium transition itself on the trajectory level is then given in terms of the initial and final distribution as .

Averaging over the initial distribution we infer from (60)

(66) |

where the distribution is balanced by the entropy contribution . Finally, integrating over the final state with the normalised weight , i.e., , we obtain the integral fluctuation theorem proposed by Seifert Seifert (2005a)

(67) |

or

(68) |

where the total entropy is given by

(69) | |||

(70) |

## V Reservoir-driven Hamiltonian system

The last issue to be considered in this section is the case of a system with many degrees of freedom. We consider an over damped system coupled to multiple reservoirs described by the Hamiltonian . For the equations of motion for , and we obtain

(71) | |||

(72) | |||