Thermodynamics of Micro- and Nano-Systems Driven by Periodic Temperature Variations

We introduce a general framework for analyzing the thermodynamics of small systems that are driven by both a periodic temperature variation and some external parameter modulating their energy. This set-up covers, in particular, periodic micro and nano-heat engines. In a first step, we show how to express total entropy production by properly identified time-independent affinities and currents without making a linear response assumption. In linear response, kinetic coefficients akin to Onsager coefficients can be identified. Specializing to a Fokker-Planck type dynamics, we show that these coefficients can be expressed as a sum of an adiabatic contribution and one reminiscent of a Green-Kubo expression that contains deviations from adiabaticity. Furthermore, we show that the generalized kinetic coefficients fulfill an Onsager-Casimir type symmetry tracing back to microscopic reversibility. This symmetry allows for non-identical off-diagonal coefficients if the driving protocols are not symmetric under time-reversal. We then derive a novel constraint on the kinetic coefficients that is sharper than the second law and provides an efficiency-dependent bound on power. As one consequence, we can prove that the power vanishes at least linearly when approaching Carnot efficiency. We illustrate our general framework by explicitly working out the paradigmatic case of a Brownian heat engine realized by a colloidal particle in a time-dependent harmonic trap subject to a periodic temperature profile. This case study reveals inter alia that our new general bound on power is asymptotically tight.


I. INTRODUCTION
Thermodynamic processes on the micro-and nanoscale in systems driven out of equilibrium by periodically changing control parameters like an external force or the temperature of their environment, can be scrutinized under the microscope by virtue of precise measurements of characteristic quantities such as applied work or exchanged heat [1][2][3][4][5][6]. Despite such groundbreaking results, so far, no general theoretical framework for the thermodynamic description of periodically driven systems beyond the quasi-static limit is available.
In principle, the concepts of irreversible thermodynamics, a phenomenological but powerful theory, which, building on the principle of local equilibrium, furnishes non-equilibrium steady states with a universal thermodynamics structure [7], can be transferred to periodic systems [8][9][10]. These results are, however, crucially tied to specific models and require rather involved and nonintuitive definitions for currents and affinities.
In this paper, we overcome these limitations. Starting from the first law formulated for an arbitrary system driven out of equilibrium by both, periodic perturbations of its Hamiltonian and the temperature of a surrounding heat bath, we obtain natural and perfectly general identifications of fluxes and affinities. Since they are defined on the level of cycle averages, these quantities are timeindependent. Nevertheless, our new formalism captures essential finite-time properties of the driven system and permits a discussion of quantities like power, which are out of reach for the laws of classical thermodynamics due to its notorious lack of time scales.
Moreover, our approach provides a universal prescription for the construction of kinetic coefficients, which fully characterize the system in the linear response regime. By using a well established and rather general stochastic approach based on a Fokker-Planck equation to describe the underlying time-dependent dynamics [11], we prove that these coefficients are interrelated by a remarkable symmetry, which, just like Onsager's celebrated reciprocity relations [12,13], can be traced back to microscopic reversibility.
Ever since James Watt's steam engine, the urge to explore the fundamental principles governing the performance of cyclic heat engines, was one of the major quests in thermodynamics. Two key figures of merit are particularly important in this context. While efficiency, the first one, is universally bounded by the Carnot value as a direct consequence of the second law, similar constraints on power, the second one, could, so far, be obtained only within specific setups, see for example [14][15][16][17][18].
Thermoelectric heat engines, which, in contrast to cyclic ones, work in a steady state [19,20] are likewise subject of substantial research efforts concerning efficiency [21,22], power [23,24] and the relation between these quantities [25,26]. At least in the linear response regime, however, much more general results are currently available for this type of heat engines than for their periodic counterparts. Specifically, it can be shown that, in the linear regime, the power of such devices is bounded by a quadratic function of efficiency that vanishes at the Carnot value η C and attains its maximum at the Curzon-Ahlborn value η C 2 [27][28][29][30]. This result follows from a quite general analysis within the framework of linear irreversible thermodynamics. It must, however, be reconsidered in the presence of a magnetic field, which breaks the Onsager symmetry, an issue which is currently under active discussion [20,[31][32][33][34][35][36][37][38][39][40].
Here, by applying our new formalism, we establish the aforementioned quadratic bound on power for cyclic heat engines in linear response. In particular, by using a novel method, which does not require any additional assumptions, we prove that this bound still holds if the matrix of kinetic coefficients is not symmetric, which is the generic case for periodically driven systems. We emphasize that, by now, an analogous result for thermoelectric heat engines is not available on such a general level. To complete our analysis, we show that this bound on power is tight within a paradigmatic model of a Brownian heat engine, which was originally proposed in [14] and recently realized in a landmark experiment [2].
The rest of this paper is structured in three major parts. In section II, we develop our general formalism and prove a generalized reciprocity relation for the kinetic coefficients of periodic systems. Section III is devoted to the discussion of cyclic heat engines in the linear response regime and the derivation of a general bound on power. We illustrate our results by considering a simple model system in section IV. Finally, we conclude in section V.

II. FRAMEWORK
In this section, we will demonstrate that the notions of irreversible thermodynamics can be transferred to periodically driven systems.

A. Nonlinear Regime
We begin with a brief review of the energetics of driven systems in thermal contact with a heat bath [30]. Specifically, we consider a classical system with degrees of freedom x ≡ (x 1 , . . . , x n ) and time-dependent Hamiltonian which is immersed in a heat bath, whose temperature T (t) oscillates between the two values T c and T h > T c . Here, g w (x, t) denotes an externally controlled dimensionless function of order 1 and ∆H quantifies the strength of this time-dependent perturbation. The power extracted by the controller thus readṡ where p(x, t) denotes the normalized probability density to find the system in the state x at the time t and dots indicate time derivatives throughout the paper. To compensate for this loss in internal energy FIG. 1. Non-equilibrium periodic state of a system with one degree of freedom x in a sinusoidally shifted harmonic potential, represented by gray parabolas, which is embedded in an environment with periodically changing temperature as indicated by the periodic color gradient. The solid line in the x-t plane shows the motion of the center of the Gaussian phase space distribution, which lags behind position of the potential minimum shown as dashed line. At any time t, the width of the colored region equals the width of the phase space distribution, which varies according to the temperature.
the system picks up the heaṫ from the environment as stipulated by the first laẇ We will now pass from time-dependent to constant variables by exploiting the periodic boundary conditions where T is the length of one operation cycle. Quite naturally, we invoke the assumption that, given these conditions, the time evolution of the probability density p(x, t) eventually converges to a periodic limit as illustrated in Fig. 1 for a simple system. Once this periodic state is reached, the average entropy production per cycle arises solely due to heat exchange between the system and the environment, since, due to the periodicity of the the distribution p c (x, t), no net entropy is generated in the system in a full cycle, i.e., we havė We now focus on the linear regime with respect to the temporal gradients ∆H and ∆T ≡ T h − T c . By expanding the fluxes (12) and (13), we obtain with linearized affinities F w = ∆H T c and F q = ∆T T 2 c + O ∆ 2 (18) and kinetic coefficients The entropy production (8) thus reduces tȯ To guarantee that this expression is nonnegative for any F α as stipulated by the second law, the kinetic coefficients must obey the constraints which we prove explicitly in Sec. III B for a large class of systems. It is, however, not evident at this stage whether a reciprocity relation relating L αβ with L βα or any further constraints exist.

Adiabatic Limit
As a first step, we investigate the adiabatic regime, which is characterized by the Hamiltonian H(x, t) and the temperature T (t) changing slowly enough in time such that the system effectively passes through a sequence of equilibrium states, i.e., with and k B denoting Boltmann's constant. Expanding (22) to linear order in ∆H and ∆T and inserting the result into (12), (15) and (19) gives the universal expression for the adiabatic kinetic coefficients. Here, we introduced the notations and for any quantity A(x, t) and the equilibrium distribution of the unperturbed system with Z 0 denoting the canonical partition function. Notably, the coefficients (24) are fully antisymmetric, i.e., L ad αβ = −L ad βα .
As might be expected, this property, which can be proven by a simple integration by parts with respect to t, implies vanishing entropy production (20) in the adiabatic limit. This avoidance of dissipation can, however, be only achieved for an infinite cycle duration T and therefore inevitably comes with vanishing fluxes J α .

Stochastic Dynamics
For further investigations of the kinetic coefficients, we have to specify the dynamics, which governs the time evolution of the probability density p(x, t). Having in mind, in particular, mesoscopic systems surrounded by a fluctuating medium, a suitable choice is given by the Fokker-Planck equation [41] with where summation over identical indices is understood and natural boundary conditions are assumed. The Hamiltonian H(x, t) and the temperature T (t) enter via the drift and diffusion coefficients D i (x, H, T ) and D ij (x, H, T ), which thus become implicitly time-dependent [42].
We will now formulate a set of conditions on the general Fokker-Planck operator (30) to adapt it to the physical situation that we wish to discuss here. Since microreversibility plays a crucial role in linear irreversible thermodynamics, we have to ensure that our theory complies with this fundamental principle. To this end, first, at any time t any possible state x of the system must be associated with the same energy as the time-reversed state εx ≡ (ε 1 x 1 , . . . , ε n x n ) with ε i = 1 for even and ε i = −1 for odd variables, i.e., where, throughout the paper, the transformation x → εx is meant to include the reversal of external magnetic fields. Second, the unperturbed Fokker-Planck opera- for the canonical distribution (27), which uniquely satisfies The dagger showing up in (32), from here onwards, designates the adjoint of the respective operator. Note that, while in (33) the operator L 0 (x) acts on the function p eq (x), (32) is to be read as an operator identity becoming meaningful when applied to a specific function of x. Physically, the relation (32) means that, once the system has reached its equilibrium state, the rate of transitions from the state x to the state x ′ is balanced by the rate of transitions in the reverse direction. The equilibrium Fokker-Planck operator L 0 (x) can be naturally decomposed in a reversible and an irreversible contribution which are characterized by their respective behavior under time reversal. While the irreversible part accounts for dissipative effects induced by the presence of the heat bath, the reversible part describes the intrinsic coupling of the system's degrees of freedom, which is not directly affected by the fluctuating environment. Since this autonomous part of the dynamics should preserve the internal energy of the system, we have to impose the condition We note that this consideration does not play a role in the overdamped limit, within which the entire time evolution of the system is effectively irreversible due to strongly dominating friction forces. The notion of detailed balance can not be immediately generalized to situations with external driving and timedependent temperature. However, in analogy with (33), the full Fokker-Planck operator L(x, t) can still be characterized by the weaker property which is naturally obeyed in the absence of nonconservative forces and guarantees that the system follows the correct thermal equilibrium state if the Hamiltonian and the temperature are varied quasi-statically.

Finite-Time Kinetic Coefficients
In the linear response regime, the Fokker-Planck operator L(t) showing up in (29) can be replaced by the expansion where, for simplicity, from (38) onwards, we notationally suppress the dependence of any operator on x, whenever there is no need to indicate it explicitly. The Fokker-Planck equation (29) can then be solved with due consideration of the boundary condition (6) by treating L H (t) and L T (t) as first order perturbations. The result of this standard procedure [41,43] reads After some algebra involving condition (37), which we relegate to appendix A for convenience, this solution leads to the compact expression for the kinetic coefficients, where the generalized equilibrium correlation function is defined as for arbitrary quantities A(x, t) and B(x, t). We recall the definition (26) of the δ-notation.
The expression (40) admits an illuminating physical interpretation. It shows that the kinetic coefficients of periodically driven can be decomposed into an adiabatic contribution independently identified in (24) and a finitetime correction, which has the form of an equilibrium correlation function. This result might therefore be regarded as a generalization of the well established Green-Kubo relations, which relate linear transport coefficients like electric or thermal conductivity to equilibrium correlation functions of the corresponding currents [43]. In our case, the role of the currents is played by the fluctuation variables δġ α (x, t) and the ensemble average is augmented by a temporal average over one operation cycle. We note that a similar expression has been obtained for the special case of the effective diffusion constant of a periodically rocked Brownian motor in [44].

Reciprocity Relations
Time-reversal symmetry of microscopic dynamics appears as the detailed balance condition on the level of the Fokker-Planck equation (32). By using this relation and the Green-Kubo type formula (40), it is straightforward to derive the generalized reciprocity relation where, the Onsager coefficients are considered as functionals of the time-dependent Hamiltonian and temperature and an external magnetic field B. The technical details of the derivation leading to the relation (42), which constitutes our second main result, can be found in appendix B. Here, we emphasize that, although the setup of this paper differs significantly from the one Onsager dealt with in his pioneering work [12,13], the symmetry (42), which constitutes our second main result, and the original Onsager relations share microscopic reversibility as the common physical origin. Since, in the presence of time-dependent driving, full time reversal especially includes reversal of the driving protocols, naturally, these reversed protocols show up in (42).
The symmetry (42) holds individually for both, the adiabatic kinetic coefficients and the finite-time correction showing up in (40). Given the general relation (28), it follows that the L ad αβ must vanish if the driving protocols are symmetric under time-reversal.
The additional relation can be proven if the driving g w (x, t) introduced in (1) factorizes according to Hence, the off-diagonal kinetic coefficients change place if the magnetic field is reversed and the respective protocols determining the time dependence of the Hamiltonian and the temperature are interchanged. This symmetry does not involve the reversed protocols. It is, however, less universal than (42), since it requires the special structure (44), see appendix B for details.

III. CYCLIC STOCHASTIC HEAT ENGINES IN LINEAR RESPONSE
As a key application of our new approach we will discuss the performance of stochastic heat engines.

A. Power and Efficiency
The main two benchmark parameters here are, first, the average power output per operation cycle (45) and, second, the efficiency which is bounded by the Carnot value η C = 1 − T c T h as a direct consequence of the second lawṠ ≥ 0. The latter figure, which is naturally suggested by the representation (11) of the entropy production per cycle, should be regarded as a generalization of the conventional thermodynamic efficiency defined for a heat engine operating between two reservoirs of respectively constant temperature. Our formalism includes this scenario as the special case where γ q (t) is chosen as a step function with 0 < T 1 < T such that the system is in contact with the hot temperature T h in the first part of the cycle and the cold temperature T c in the second one. Under linear response conditions, which we will assume for the rest of this section, the fluxes J α can be eliminated using (17) such that the expressions (45) and (46) reduce to and respectively. Clearly, these figures are crucially determined by the kinetic coefficients L αβ . In contrast to the thermoelectric case, where the reciprocity relation L αβ = L βα holds without magnetic fields, the analysis of the preceding section has revealed that for cyclic heat engines this symmetry is typically broken if the driving protocols are not invariant under time reversal. As pointed out by Benenti et al. [31], a non symmetric matrix of kinetic coefficients leads to profound consequences for the performance of thermoelectric devices including the option of Carnot efficiency at finite power. These results apply similarly to the systems considered here, since our theoretical framework is structurally equivalent to the standard theory of linear irreversible thermodynamics used as a starting point in [31]. To demonstrate this correspondence explicitly, following [31], we define the dimensionless parameters which, due to the second law (21), are related by the inequalities By optimizing (48) and (49), respectively, as functions of F w the expressions Lqq as a function of the normalized efficiencyη ≡ η η C for selected values of the asymmetry parameter x ≥ 1 in the upper and x < −1 in the lower panel. For x → ±∞, convergence to P P0 =η is observed. The dotted black lines correspond to the maxima of (57) for 0 ≤ x < ∞ in the first and −∞ < x < −1 in the second plot thus indicating the relation between maximum power and efficiency at maximum power. The apparent divergence for negative x occurs in the limit x → −1 [47].
for maximum efficiency and for efficiency at maximum power [30] are obtained, where η C ≈ ∆T T c = T c F q denotes the Carnot efficiency in the linear regime. For y = h(x) and x ≥ 0, the maximum efficiency equals η C and the efficiency at maximum power can exceed the Curzon-Ahlborn value η C 2 [27-30, 45, 46], reaching even η C in the limit x → ±∞.
For a quantitative assessment of the relation between power and efficiency, following [39], we consider the max-imum power at given efficiency as a joint benchmark parameter, which is found by eliminating F w in (45) in favor of using (49). To ensure that the Carnot value η C is included in the range of accessible efficiencies, y must be replaced by its bound h(x) in (55). The resulting function is plotted in Fig. 2. While P (η, x) reduces to in the symmetric case x = 1 and thus vanishes linearly in the limit η → η C , strikingly, we observe that Carnot efficiency might be reached at finite power [47] P This result is a priori surprising, since this analysis has fully incorporated the constraints imposed by the second law.

B. A New Constraint
We will now prove the existence of an additional constraint on the kinetic coefficients (19), which, so far, has been missing in our considerations. To this end, we introduce the symmetric matrix where will play the role of a normalization constant and ⟪•⟫ was defined in (25). The matrix A has the nontrivial property of being positive semidefinite such that the determinant of any of its principal submatrices must be nonnegative, as we show in appendix C by using only the rather general assumptions of Sec. II B. Two important implications follow immediately from this insight. First, by taking the determinant of the lower right 2 × 2-submatrix, we recover the inequality (21), which we inferred from the second law on the phenomenological level and now have proven explicitly. Second, by evaluating the determinant of A, we get the new constraint which, in contrast to the bare second law, leads to a bound on power. Specifically, this bound is found by bounding L qq in (55) using (63) and then maximizing the resulting function with respect to y [39]. This procedure yields the simple result where we define the standard power and y becomes through the optimization. Remarkably, for any x ≥ 1 the bound (64), which is our third main result, restores the quadratic relation between power and efficiency (58), which we found in our first analysis not invoking the new bound (63) only for the symmetric case x = 1. Consequently, we have shown that, in the linear response regime, the power of any cyclic heat engine comprised by our theoretical framework must vanish at least linearly as its efficiency approaches the Carnot value. We emphasize that this quite natural result can neither be derived from the laws of thermodynamics nor from micro-reversibility, which appears in form of the reciprocity relation (42). Instead, it relies on the additional constraint (62), which is beyond both of these principles.

A. Model and Kinetic Coefficients
A particularly simple setup for a stochastic heat engine consists of a Brownian particle in one spatial dimension confined in a harmonic potential of variable strength κ(t) and immersed in a heat bath of time-dependent temperature T (t) as schematically shown in Fig. 3. Originally proposed in [14], this model has recently been realized in a remarkable experiment [2] and can be used to illustrate various aspects of stochastic thermodynamics like Here, by applying our general theory developed in the last sections we will calculate the kinetic coefficients for this stochastic heat engine and optimize the protocol κ(t) controlling the trap strength to obtain maximum power for given efficiency.
In the overdamped limit, due to the absence of kinetic energy, the Hamiltonian of the system depends only on the position x of the particle. Here, κ 0 is the equilibrium strength of the trap, κ quantifies the strength of the time-dependent driving, γ w (t) denotes the driving protocol and x 0 ≡ 2k B T c κ 0 is the characteristic length scale of the system. The time evolution of the probability density p(x, t) for finding the particle at the position x is generated by the Fokker-Planck operator Additionally, the optimal protocol γ * w (t) has been rescaled by the dimensionless factor κ0η C κ.
which, in equilibrium reduces to Here, µ denotes the mobility and the temperature Formula (40) can be easily evaluated by using the detailed balance relation (32) to transform L 0 into L 0 † , since it is readily seen that the function x 2 − k B T c κ 0 is a right eigenvector of L 0 † with corresponding eigenvalue −2µκ 0 . The resulting kinetic coefficients are functionals of the protocols γ α (t). Note that, besides the general reciprocity relation (42), these coefficients also satisfy the special symmetry relation (43), since the factorization condition (44) is fulfilled in the example discussed here.

B. Optimization
The optimal protocol γ * w (t, η) for the strength of the harmonic trap for a given time dependence of temperature γ q (t) maximizes the power output at fixed efficiency η. It is determined by the variational condition where the power P and the heat flux J q are regarded as functionals of γ w (t) and γ q (t). As we show in appendix D, this constrained optimization problem has the general solution for 0 ≤ t ≤ T , where we used the abbreviations (56), and γ 0 denotes an arbitrary constant. Using this protocol, the maximum power can be extracted per operation cycle at efficiency η.

C. Comparison with General Bound
In order to compare this result with the general bound (64), we evaluate the normalization constant defined in (61) and rewrite (76) as whereP 0 is the standard power introduced in (64). The dimensionless factor quantifies how close the maximum power found in the optimization comes to the general bound (64). Since, 0 ≤ γ q (t) ≤ 1, it is reached only forγ q → 0. This limit, however, requires γ q (t) → 0 for any t and thus, inevitably, leads to vanishing absolute power.
For an illustration of this issue, we chose γ q (t) as the step function (47) such that the system is alternately in contact with a hot and a cold bath, respectively during the time intervals T 1 and T − T 1 . We then find ψ = 1 − T 1 T and P 0 ∼ T 1 T . Thus, as T 1 is decreased, ψ comes arbitrarily close to 1 and P 0 decays linearly to zero. This example shows that our bound is asymptotically tight.
A particular advantage of our approach is that it allows to treat situations with a continuously varying temperature of the environment on equal footing with the scenario proposed in [14], which involves instantaneous switchings between a hot and a cold reservoir. In order to illustrate this feature, we consider the specific choice for the protocol γ q (t), which, in the linear response regime, is proportional to the temperature T (t). The function (80), which interpolates between a step function (d → 0) and a simple sine (d → ∞) [52], is plotted together with the corresponding optimal protocol γ * w (t, η, d) [51] in Fig. 4. We find that, for d = 0, this protocol shows two sudden jumps occurring simultaneously with the instantaneous changes of the bath temperature. Such discontinuities were shown to be typical for thermodynamically optimized finite-time protocols connecting two equilibrium states [53]. Since, here, we are concerned with periodic states generated by permanent driving rather than a transient process with equilibrium boundary conditions, it is, however, not surprising that both, the temperature and the trap strength protocol, become continuous as the shape parameter d is increased and, respectively, converge to two sines in phase in the limit d → ∞.
For the protocols (73) and (80), the parameter ψ becomes This function decays monotonically from 1 2 for d = 0 to 1 3 for d → ∞. Consequently, ψ can not reach its maximum 1 within the class of protocols (80). This limitation can be understood from the argument given below (79), since, for any d, we haveγ q = 1 2. The standard power P 0 is proportional to the function 2 + d − d(1 + d), which decays monotonically from 2 to 3 2 as d increases from 0 to infinity. Thus, P 0 exhibits the same qualitative dependence on the shape parameter d as the efficiency. We can therefore conclude that, at least within the model considered here, a steeply rising and falling temperature performs better than a smoothly changing one.

V. CONCLUDING PERSPECTIVES
In this work, we have demonstrated that nonequilibrium periodic states, which emerge naturally in periodically driven systems, can be endowed with the universal structure of irreversible thermodynamics. Moreover, by using a quite general stochastic approach, we have proven the generalized reciprocity relation (42) for the kinetic coefficients characterizing the linear response regime.
Our new framework is particularly useful for a systematic study of the performance of cyclic heat engines. Within the linear regime, bounding the power of these machines in such a way that the rather peculiar option of Carnot efficiency at finite power is ruled out requires the new relation (62). This constraint is beyond the laws of thermodynamics and the principle of microscopic reversibility and has been proven here for the first time on a general level. Remarkably, up to the normalization factor, an identical bound has been discovered only recently in a numerical analysis of a particular class of thermoelectric heat engines [39]. Whether this intriguing similarity suggests the existence of a so far undiscovered universal principle that applies to periodic as well as to steady states leading to a bound on power for any type of heat engine that operates in the linear regime remains an exciting topic for future research.
A promising starting point for investigations in this directions might be found in the Green-Kubo relations, which follow from first principles and provide general expressions for the conventional kinetic coefficients in terms of equilibrium correlation functions [43]. Using a Fokker-Planck approach, we have shown that an analogous representation for the kinetic coefficients exists in periodically driven systems. The quantities related by the relevant correlation functions are, however, well defined irrespective of specific dynamics governing the time evolution of the phase space density. It might therefore be possible to obtain Hamiltonian-based expressions also for the periodic kinetic coefficients introduced in this work. Finding a proper way to take the time dependence of temperature into account is arguably the major challenge here.
This problem also prevents an immediate extension of our formalism to the quantum realm. While the first part of our analysis, the identification of proper fluxes and affinities, carries over to quantum mechanics line by line, it is not clear at the moment whether the constraints on the kinetic coefficients obtained here classically can be likewise transferred or properly generalized. This topic appears all the more urgent in the light of recent developments [5,54,55] showing that the emerging field of quantum thermodynamics nowadays comes within the range of experiments. For an expression of the kinetic coefficients (19) depending only on equilibrium quantities the perturbations L X (t) showing up in the linear response solution (39) have to be eliminated. To this end, we invoke the property (37) of the full Fokker-Planck operator. Substituting (38) into (37), expanding the exponential in ∆H and ∆T and collecting linear order terms provides us with the relations where we used the definition (16). Up to corrections of order ∆ 2 , the periodic distribution p c (x, t) can thus be rewritten as In the second line, we replaced g α (x, t) with its equilibrium fluctuation defined in (26). This modification does not alter the right hand side of (A2), since L 0 p eq (x) = 0. However, it ensures that the function, on which the exponential operator acts in the third line, which was obtained by an integration by parts with respect to t, has no overlap with the nullspace of L 0 and thus the integral with infinite upper bound is well defined. Note that the upper boundary term vanishes, since the operator L 0 is nonpositive [41]. Inserting (A2) into (12) and (15) yields (A5) Herein, obviously,ġ α (x, t) can be replaced by its equilibrium fluctuation δġ α (x, t) in the first term. Since any constant lies in the left nullspace of L 0 and, hence, is orthogonal to the function e L 0 τ δġ α (x, t−τ )p eq (x), the same replacement can be carried out in the second term without leading to additional contributions such that (40) is, finally, obtained.

Appendix B: Reciprocity Relations
In this appendix, we establish the reciprocity relations (42) and (43). To this end, we first recall formula (40), which becomes using the definitions (25) and (41). By applying the detailed balance relation (32) and changing the integration variable t according to t → T − t, this expression can be transformed to where we introduced the shorthand notationL 0 ≡ L 0 (εx). Furthermore, to transfer the variable τ from the argument of δġ β in (B1) to the argument of δġ α in (B2), we used the identity which holds for any two T -periodic functions a(t) and b(t). Finally, we reverse the magnetic field as well as the driving protocols and apply the change of integration variables x → εx, whose Jacobian is 1. By exploiting the symmetry δg α (εx, t) = δg α (x, t), which follows from condition (31), and carrying out one integration by parts with respect to t in the first summand, we obtain thus completing the proof of the reciprocity relation (42).
We now turn to the special case where the function g w (x, t) can be separated in the form g w (x, t) = g w (x)γ w (t).
(B6) Plugging this expression into (B1) and invoking the definition g q (x) ≡ −H 0 (x) yields the expression which, by virtue of the detailed balance relation (32), equals Relation (43) can now be obtained by following the same steps as in the general case, that is, by applying the transformation x → εx, reversing the magnetic field and interchanging the arguments γ w (t) and γ q (t) of L αβ , using the symmetry δg α (x) = δg α (εx) and performing one integration by parts in the first term.

Appendix C: Positivity of A
The proof that the matrix A defined in (60) is positive semidefinite consists of two major steps. First, for arbitrary numbers y w , y q ∈ R, we consider the quadratic form Q 0 (y w , y q ) ≡ α,β=w,q L αβ y α y β , (C1) where we used the expression (40) for the kinetic coefficients, defined G(x, t) ≡ α=w,q y α δġ α (x, t) (C2) and applied the detailed balance relation (32) to obtain the second line. We recall the definitions (25) for the meaning of the angular brackets. The crucial ingredient for this first step consists of the identity Here, we used the relation which holds for any functions A(x), B(x) by virtue of the detailed balance condition (32), to obtain (C4) from (C3). Expression (C5) follows by applying an integration by parts with respect to τ and τ ′ , respectively, in the first and the second summand of (C4). Finally, the second contribution in (C5) vanishes after carrying out the tintegration, since the function G(x, t) is T -periodic in time.
Next, we note that (32) implies ⟨A L 0 † +L 0 † A⟩ (C8) = d n x A(x) L 0 p eq (x) + p eq (x)L 0 † A(x) = d n x [p eq (x)] Since the Hermitian part of this operator is negative semidefinite [41], it follows that (C8) is non-positive for any A(x). Hence, we can conclude that the quadratic form Q 0 (y w , y q ) is positive semidefinite, since it can be written in the form (C3). For the second step of the proof, we introduce the quadratic form Q(y w , y q , z) ≡ Q 0 (y w , y q ) + Q 1 (y w , y q , z) (C10) with Q 1 (y w , y q , z) ≡ N qq z 2 + 2z α=w,q L qα y α (C11) where y w , y q , z ∈ R and F (x, t) ≡ zδg q (x, t).
Note that, in (C11), we used (40) and (61) as well as the detailed balance condition (32). The expression (C11) can be rewritten as This assertion can be proven by expanding (C13), invoking (C7) as well as the identity which is implied by condition (36), and integrating by parts, first with respect to τ and then with respect to t, respectively, in the second and third term showing up in (C13). Finally, putting together (C3) and (C13) leads to Q(y w , y q , z) = − 1 The average showing up in this expression is of the form (C8) and thus must be non positive. Consequently, we have Q(y w , y q , z) ≥ 0 for any y w , y q , z. Moreover, since Q(y w , y q , z) = y t Ay (C16) with y ≡ (z, y w , y q ) t it follows that the matrix A must be positive semidefinite and thus the proof is completed.