Thermodynamics of active field theories: Energetic cost of coupling to reservoirs

The hallmark of active matter is the autonomous directed motion of its microscopic constituents driven by consumption of energy resources. This leads to the emergence of large scale dynamics and structures without any equilibrium equivalent. Though active field theories offer a useful hydrodynamic description, it is unclear how to properly quantify the energetic cost of the dynamics from such a coarse-grained description. We provide a thermodynamically consistent framework to identify the energy exchanges between active systems and their surrounding thermostat at the hydrodynamic level. Based on linear irreversible thermodynamics, we determine how active fields couple with the underlying reservoirs at the basis of nonequilibrium driving. This leads to evaluating the rate of heat dissipated in the thermostat, as a measure of the cost to sustain the system away from equilibrium, which is related to the irreversibility of the active field dynamics. We demonstrate the applicability of our approach in two popular active field theories: (i) the dynamics of a conserved density field reproducing active phase separation, and (ii) the coupled dynamics of density and polarization describing motile deformable droplets. Combining numerical and analytical approaches, we provide spatial maps of dissipated heat, compare them with the irreversibility measure of the active field dynamics, and explore how the overall dissipated heat varies with the emerging order.

In recent years, a large number of works have focused on developing a thermodynamic approach to active matter. They are led by the search for generic observables to quantify, classify and predict the anomalous properties of active systems. For instance, the pressure and chemical potential allow one to distinguish systems depending on whether or not they obey equations of state [35][36][37][38], and each can be useful to predict phase diagrams [39,40]. Moreover, quantifying the irreversibility of the dynamics enables us to locate where and when activity mainly affects the system [30,41,42], and to explore the relation between irreversibility and phase transitions [43,44]. This has motivated several experiments to measure the dissipation associated with irreversibility in various sys-tems [45][46][47][48]. Furthermore, it has been shown that, for minimal active models, changing the dissipation by using a dynamical bias provides an alternative route to clustering and collective motion in active matter [49][50][51][52][53].
Progress in building the thermodynamics of active matter has been mainly achieved so far in particle-based descriptions. Indeed, such dynamics bear a natural mechanical interpretation of energy exchanges, with either an external operator or the surrounding thermostat, in terms of forces and displacements. The tools of stochastic thermodynamics, first introduced in thermal systems [54,55] and then extended to active ones [41,[56][57][58][59][60][61], offer a framework to quantify systematically work, heat, and entropy production from the microscopic dynamics. They also allow one to describe the consumption of chemical fuel, at the basis of self-propulsion, in a thermodynamically consistent manner [62][63][64]. In contrast, though some hydrodynamic approaches consider coupling with a momentum-conserving fluid [1,[65][66][67], nonequilibrium terms in many field dynamics do not rely on any explicit mechanical force [21,[29][30][31][32]. Then, a systematic definition of how active systems exchange energy with their environment at a coarse-grained level has been elusive: It remains to build the energetics of active field theories.
A major breakthrough of stochastic thermodynamics is to relate explicitly the irreversibility of the dynamics, as measured by the divergence of forward and time-reversed realizations, with the amount of energy dissipated in the thermostat [54,55]. This connection only holds for thermodynamically consistent dynamics, whose formulation is constrained so that the connection to the underlying thermostat is properly taken into account. Importantly, the active field theories postulated only from symmetry arguments do not generally satisfy these constraints a arXiv:2008.06735v2 [cond-mat.stat-mech] 19 Jun 2021 priori [21,[29][30][31][32]. Hence, it is unclear to which extent the measure of irreversibility in these models, often referred to as entropy production rate (EPR) and already evaluated in previous works [30,42,68], actually provides relevant information about energy dissipation.
Interestingly, linear irreversible thermodynamics (LIT) provides a definition of dissipation in terms of the thermodynamic forces and conjugated currents at a hydrodynamic level [69]. After identifying the relevant forces and currents for a given theory, the field dynamics are formulated by postulating linear relations between them. These theories were originally designed to capture the effect of external drives, such as temperature gradients or electric fields, yet extensions to systems with internal activity, such as active gels, have been shown to successfully reproduce the behavior of living materials [70][71][72]. Though some active theories do not follow linear forcecurrent relations a priori, it is tempting to draw analogies with LIT in order to systematically define dissipation in this broad class of dynamics beyond active gels. The challenge is then to embed active field theories with arbitrary nonequilibrium terms into the specific structure of LIT, thus enforcing a thermodynamically consistent framework. This is a non-trivial task with the benefit of drawing powerful generic results from thermodynamic considerations.
In what follows, we offer a framework to evaluate the energetic budget of active field theories. Starting from first principles, we demonstrate how to define the heat rate dissipated to the thermostat from the fluctuations of the active fields. Importantly, we show that the heat rate can be generically decomposed into a homogeneous background contribution, independent of active fields, and a contribution given in terms of the statistics of the active fields. This decomposition allows one to quantify how the structure and dynamics of active fields affect where heat is dissipated, thus opening the door to estimating and comparing the energetic cost associated with different emerging orders. To illustrate the relevance of our framework, we apply it to field theories which capture the emergence of a phase separation and/or a polar order. Overall, our results demonstrate the ability to estimate the rate of energy required to sustain a given active dynamics away from equilibrium.
Our approach relies on systematically constructing the dynamics of a set of underlying fields, which drive the system out of equilibrium, from that of the active field dynamics based on the force-current relations of LIT. Importantly, under non-restrictive assumptions (see Sec. II A), the evolution of active fields remains independent of that of the driving fields: The latter are hidden degrees of freedom that do not affect the emerging order. This leads to show that the heat rate, whose expression follows from the total EPR measuring the irreversibility of both the active and driving fields, can be evaluated from the fluctuations of active fields only. Importantly, the heat rate is distinct from the EPR quantifying the irreversibility of the active field dynamics alone, which we refer to as the explicit EPR in what follows.
We analyze in detail the heat rate in two popular models for active matter: (i) the dynamics of active phase separation, known as Active Model B [29,30], and (ii) the dynamics of polar motile droplets [73,74]. For model (i), we find that the heat-rate mainly varies at the interfaces between dense and dilute phases, where it reduces compared to its bulk value. We further evaluate the heatrate scaling with noise amplitude and driving parameter. Our results are compared with the explicit EPR, as an alternative measure of the deviation from equilibrium [30,42,68]. The analysis of model (ii) reveals that the rate of dissipated heat varies across the profile of polar droplets. We also report a hysteresis loop associated with the splitting and fusion of multiple droplets, and we discuss the scaling with noise amplitude and driving parameter.
The paper is organized as follows. First, we present in Sec. II how to embed generic active field theories within LIT and calculate the heat rate, which we then relate with the explicit EPR in Sec. III. In Sec. III A, we consider an application of our framework to dynamics that capture active phase separation in terms of a conserved scalar field [29,30]. To go beyond scalar field theories, we then study in Sec. III B the dynamics of motile droplets coupling density and polarization [73,74] as a prototypical model of deformable living cells [75][76][77]. Finally, we present our conclusions in Sec. IV.

II. THE ROLE OF UNDERLYING RESERVOIRS
We consider active dynamics of hydrodynamic fields which can be either obtained from explicit coarsegraining of microscopic dynamics, or written phenomenologically using symmetry arguments [1,3]. Our approach consists in introducing additional fields, associated for instance with chemical reactions which sustain the dynamics away from equilibrium, see Fig. 1. This amounts to identifying the nonequilibrium terms in the original dynamics as a coupling to chemical reservoirs following the framework of linear irreversible thermodynamics [69]. Below, we present in detail the procedure to enforce a thermodynamically consistent structure of the dynamics: First, for a conserved scalar field, and then for generalized field dynamics which couple a conserved scalar field and a polarization field. The key in providing a thermodynamic framework for active materials is to realize that such materials are typically a part of a larger system, which provides the drive needed to sustain nonequilibrium activity, as described in Fig. 1.

A. Coupling active and chemical fields
To introduce pedagogically our framework, we start by considering the dynamics of a conserved scalar field φ fuel reservoir, e.g. ! ! product reservoir, e.g. ! + ! active system: FIG. 1. Schematic representation of an active system (blue) put in contact with reservoirs of chemical fuel (red) and product (green) which set a constant, homogeneous chemical potential difference ∆µ in the active system. This is essentially a nonequilibrium grand-canonical ensemble for the active system (details in Appendix A). Within our framework, ∆µ embodies the driving parameter which controls the nonequilibrium terms in the dynamics (1-4) for the active density field φ and the rate of fuel consumptionṅ. The active system and the chemical reservoirs are surrounded by the thermostat (yellow) which maintains a fixed temperature T . The fluctuations of φ and n lead to dissipation of heat Q into the thermostat, which quantifies the energetic cost to maintain the whole system away from equilibrium.
representing the density of active components: where F is the free energy, λ is the mobility, ∆µ is the driving coefficient, and C is a vector-valued function of φ and its gradients. The noise term Λ is Gaussian with zero mean and correlations given by where T is the temperature of the surrounding heat bath. The term T ν is a generalization of the spurious drift that typically appears in ordinary stochastic differential equations with multiplicative noise. Its expression is determined by that of C, it depends on both time and space discretizations, and it obviously vanishes when fluctuations are neglected (T = 0). In Appendix B we generalize standard results for stochastic processes [78] to stochastic field theories and derive the expression for the spurious drift term. The dynamics (1) has been used extensively to reproduce the phase separation of active particles [26,27,[29][30][31][32]. In these works, the need for the additional term T ν was not addressed explicitly, mainly because previous studied were not concerned with thermodynamic consistency, and also since the noise Λ seems to be additive when considering only the fluctuations of φ. When we describe below the origin of ∆µ, by introducing additional field dynamics, it will become apparent that the noise Λ is in fact multiplicative due to its cross-correlation with the noise of the additional field, as described in Eq. (4). Our goal is to connect the emergent behavior of φ with the underlying consumption of energy resources. To this end, we describe explicitly the fluctuations of the degrees of freedom at the basis of nonequilibrium drive, referred to as chemical fields in what follows, though our framework extends to other types of drive. Inspired by recent works [68,79], we regard the driving coefficient ∆µ as the chemical potential difference between fuel and products of a chemical reaction, see Fig. 1, which applies, for instance, to the oxidation of hydrogen peroxide involved in the self-propulsion of Janus colloids [15][16][17]. This leads us to consider the dynamics of the chemical coordinate n, which is (half) the difference between the local number density of product molecules and that of the fuel molecules (see Appendix A). It is described as a field fluctuating in space and time, while ∆µ is kept constant and homogeneous.
We aim at proposing a systematic method to couple the active field φ and the underlying chemical field n. It relies on the fact that the active system is a part of a large nonequilibrium system that relaxes (slowly) towards equilibrium. With this assumption, the explicit dynamics of n can be deduced from linear irreversible thermodynamics (LIT) [34,69,72,80,81]. Identifying J and −∇(δF/δφ) as the current and the thermodynamic force associated with φ, respectively, LIT states that the currents {J,ṅ} can be written as a linear combination of the thermodynamic forces {−∇(δF/δφ), ∆µ}. It is clear from (1) that the factor coupling the current J and the force ∆µ is directly given by C. Accordingly, and because φ is even under time-reversal, Onsager reciprocity relations require that the coupling factor between the currenṫ n and the force −∇(δF/δφ) is also C [82], so that the dynamics of n follows aṡ where γ is the chemical mobility, which we take constant in what follows. As a result of this assumption, the equation for φ is autonomous and does not rely on knowing the fluctuations of the chemical field n. The noise term ξ is Gaussian with zero mean and correlations given by Note that, though LIT states linear relations between forces and currents, the coupling factor C need not be linear with respect to φ and its gradients. It is convenient to introduce the Onsager matrix L which gives the coupling between forces and currents in d + 1 dimensions. For d = 3, it is given by Then, the dynamics (1-4) can be expressed in a compact form as [83] J,ṅ = L − ∇ δF δφ , ∆µ + T ν, χ + Λ, ξ , where the noise correlations read and denotes transpose. The expression of {ν, χ} can be obtained from that of the Onsager matrix L following a systematic route, as detailed in Appendix B. In particular, it depends on the choice of how the gradient terms appearing in C are discretized in space, see Appendix B. In the specific examples considered below, a judicious choice of the discretization can be made such that the spurious drift vanishes. Moreover, one can show that ν = 0 for d = 1, and that {ν, χ} both vanish whenever C is a local function of φ independent of its gradients.
The dynamics (1-4) is thermodynamically consistent in the sense that it obeys detailed balance, and thus relaxes to an equilibrium state at temperature T , when ∆µ derives from a given chemical free energy F ch so that ∆µ = −δF ch /δn [69,72,80,81], see Appendix A. Equilibrium relaxation also requires that the Onsager matrix is positive semi-definite (det L ≥ 0). When considering the dynamics within the active system, ∆µ can be regarded as constant, see the nonequilibrium grandcanonical ensemble described in Appendix A. Then, the realizations of the active field φ are independent of that of the chemical field n, and the dynamics now operates away from equilibrium. Although the realizations of φ are independent of n, the presence of n determines the existence of the spurious drift term ν and thus affects the φ dynamics. Within this grand-canonical description,ṅ is the important field (rather than n) and it should be thought of as the local rate of chemical reactions.

B. Dissipation and irreversibility
The nonequilibrium drive ∆µ breaks time-reversal symmetry and leads to dissipation of energy in the form of heat Q from the system to the surrounding thermostat. Following stochastic thermodynamics [54,55,84], the heat along a trajectory can be evaluated from the irreversibility of the dynamics. It amounts to comparing the path probabilities of the forward and time-reversed dynamics, respectively denoted P and P R , which quantify the probability of observing a trajectory of the currents {J,ṅ} within a given time interval [0, t] 1 : The steady state heat rateQ is theṅ where the average is taken with respect to noise realization (or P {J,ṅ} t 0 ). In equilibrium, the dynamics are symmetric under time reversal with the same statistics for forward and backward trajectories (P = P R ), so that the system does not dissipate any heat (Q = 0). In the presence of nonequilibrium drive in steady state, timereversal symmetry is broken (P = P R ) which yields a constant rate of dissipation in steady state (Q > 0).
The irreversibility of the dynamics can also be evaluated at the level of active field alone: The irreversibility measure S, referred to as explicit entropy production rate in what follows, has been evaluated in various active dynamics, either particle-based [41,56,57,59,60,62] or field theories [30,42,68], to assess unambiguously the deviation from equilibrium. Our approach differs in that we not only account for the irreversibility of active fields, but also for that of underlying chemical degrees of freedom. In this extended phase space, provided that it accounts for all the relevant hydrodynamic fields, the irreversibility indeed measures the heat dissipated by the entire system. Following standard procedures [88][89][90], the dynamic action A which sets the path probability P ∼ e −A reads where´V is a spatial integral over the whole volume V of the system, and L −1 is the inverse of L. We regard the currents {J,ṅ} and forces {−∇(δF/δφ), ∆µ} as odd and even under time reversal, respectively. The action for the time-reversed dynamics A R is then deduced readily from (11) by flipping the sign of [J,ṅ]. Substituting P ∼ e −A and P R ∼ e −A R into (9), the dissipation rate follows from straightforward algebra as (see Appendix C) where lim t→∞ 1 t´t 0 · ≡ · t is the steady-state time average. For ergodic systems, the two averages are the same and one may be omitted. Hereafter, we use · to denote both averages. The expression (12) features the sum of products between thermodynamic forces and conjugate currents, analogously to the dissipation rate in LIT [34,69,72,80,81]: This confirms that we embed active field theories within a thermodynamically consistent framework. Note that the product is interpreted here and in what follows with Stratonovich convention.
As a result, the steady-state heat rateQ equals the rate of work injected by the nonequilibrium drive ∆µ to sustain the dynamics away from equilibrium: This is equivalent to the first law of thermodynamics, as expected when the path probabilities include all thermodynamically relevant fields. The expression (13) would actually be the same if insteadṅ was held constant and ∆µ allowed to fluctuate.
For an equilibrium dynamics where ∆µ derives from the chemical free energy F ch , (∆µ = −δF ch /δn), the heat rate rate vanishes in steady state (Q = −d F ch /dt = 0), as expected. Substituting the chemical dynamics (3) in (13), we deducė Hence, the heat rate can be separated into (i) a homogeneous contribution γV ∆µ 2 corresponding to a background term independent of the fluctuations of the active and chemical fields {φ, n}, and (ii) a contribution determined only by the fluctuations of the active field φ, namely independent of that of n. The existence of n however, is crucial in determining the form of the heat rate. This becomes clear below when we compare the heat rate with the explicit EPR, in which the dynamics of n are not accounted for, see Eq. (26). Note that fast-relaxing fields which are deliberately omitted in our hydrodynamic description can only contribute to heat rate through an additional background term. Interestingly, this homogeneous contribution is eliminated when considering the difference of heat rates at constant ∆µ, for instance by changing parameters of the free energy F: The heat-rate difference then depends only on how the fluctuations of the active field φ vary with such parameters.

C. Generalized field dynamics
To demonstrate that our framework is indeed relevant for a large class of active field theories, we now consider the coupled dynamics of a conserved scalar field φ and a polar field p: where λ Ω and ∆ Ω are respectively the mobility and the constant driving coefficient for Ω ∈ {φ, p}, and C Ω depend on {φ, p} and their gradients. The noise term Λ Ω is Gaussian with zero mean and correlations given by In what follows, we assume that ∆ φ and ∆ p are independent, so that each one of ν Ω is only determined by the corresponding C Ω . The dynamics (15-16) typically describe the coarse-grained dynamics of polar agents, ranging from vibrated grains [18,19] to bird flocks [20,21] and aligning bacteria [91,92]. In practice, the dissipation rate for systems featuring other types of order parameters, such as a nematic tensor [93][94][95][96] or a non-conserved scalar field [97], extends straightforwardly from the results detailed below for the specific dynamics (15)(16). Note that in all these examples both φ and p are structural order parameters, and are therefore even under time-reversal.
The spurious drift terms T ν Ω were not considered in previous work. In what follows, we address cases where the driving coefficients ∆ Ω are either odd or even under time reversal, and we assume that even (odd) driving represents a chemical potential difference ∆ Ω = ∆µ Ω (chemical current ∆ Ω =ṅ Ω /γ Ω ). We show in Appendix B that the expression for ν Ω in terms of C Ω depends on the choice for the parity of ∆ Ω . Besides, we put forward explicit cases where ν Ω vanishes for judicious choices of the spatial discretization of gradient terms in C Ω .
With the assumption that the fields {φ, p} are even under time-reversal, LIT enforces that the form of the chemical dynamics is identical for either choice ∆ Ω = ∆µ Ω or ∆ Ω =ṅ Ω /γ Ω [34,69,72,80,81]: where γ Ω is the chemical mobility, and ξ Ω is a zero-mean Gaussian noise with correlations The noises ξ Ω and Λ Ω are correlated only if the driving is even (∆ Ω = ∆µ Ω ), in which case The expression of χ Ω follows from that of C Ω , as detailed in Appendix B, and it differs according to whether ∆ Ω = ∆µ Ω or ∆ Ω =ṅ Ω /γ Ω . We stress that in both cases the realizations of the active fields {φ, p} are independent of the chemical dynamics (17). An alternative formulation of the dynamics, not considered explicitly here, consists in taking fluctuating ∆ Ω in (15) and setting the conjugated chemical degree of freedom constant. Within this formulation, the chemical dynamics affects directly the active field dynamics, but the results for the heat rate below will remain the same, namely they only depend on the driving mechanism and not on how it affects the active dynamics.
The steady-state heat rate is now defined bẏ and it can be obtained following a similar procedure as that in Sec. II B. It again differs from the explicit entropy production rate S given by Identifying the thermodynamic forces and their conjugated currents as The expression (22) is then valid for either ∆ Ω = ∆µ Ω or ∆ Ω =ṅ Ω /γ Ω . It extends to an arbitrary number of active fields, potentially including other types of order parameters such as nematic tensors, and it remains valid when each active field couples to several chemical fields, see Appendix C: For any of these cases, the heat rate actually follows directly from the dynamics ofṅ Ω . Substituting the chemical dynamics (17) in (22), when ∆ Ω is a force (∆ Ω = ∆µ Ω ), we geṫ In general, ∆ φ and ∆ p need not have the same parity, so that the heat rate can be a combination of the forms given in (23)(24).

III. APPLICATIONS TO ILLUSTRATIVE FIELD THEORIES
Before applying our generic theory to quantify the heat rate in specific models, we compare the heat rate (14) with a measure of deviation from equilibrium obtained in previous works [30,42,68]. Substituting in (14) the expression of ∇(δF/δφ) taken from the dynamics (1) yieldṡ (25) Because previous works did not account for the spuriousdrift terms, a proper comparison requires consideration of the special case in which {ν, χ} = {0, 0} and V C · Λ dr = 0. These expressions depend on the spatial discretization scheme. In Appendix B we provide a recipe for calculating these expressions and give examples in which they vanish. In such caseṡ and the explicit entropy production rate S, which was also computed in previous studies [30,42,68], reads Therefore, (26) provides a connection between the thermodynamic heat rateQ and the explicit entropy production rate S. From the semi-positivity of the Onsager matrix L, which ensures det L = λγ − C 2 ≥ 0, it then follows that T S is a lower bound toQ. The bound is saturated when J andṅ are proportional (det L = 0): In such a case, the fluctuations ofṅ are slaved to that of J, so that the irreversibility of the whole dynamics can be found from trajectories of J alone. For the generalized field dynamics that includes the dynamics of both φ and p, we again consider the case in which {ν Ω , χ Ω } = {0, 0} and´V C Ω · Λ Ω dr = 0. Then, substituting in (23) the expression of {∇(δF/δφ), δF/δp} taken from the dynamics (15) for ∆ Ω = ∆µ Ω we geṫ where This shows explicitly the difference between the heat ratė Q and the explicit entropy production rate S, similarly to (26). For ∆ Ω =ṅ Ω /γ Ω , we have insteaḋ in which case the heat rate differs from the explicit entropy production rate by a background term.
We are now in the position to apply our generic theory to two popular active field theories: (i) the dynamics of a conserved scalar field which reproduces active phase separation, and (ii) the coupled dynamics of a conserved scalar field and a non-conserved polar field that captures the behavior of motile deformable droplets.

A. Active phase separation
To illustrate how our framework can quantify the heat rate to sustain a phase separation away from equilibrium, we consider a popular active field theory for a conserved scalar field φ that is even under time-reversal, known as Active Model B [29,30]. Taking the coupling term as C = −∇(∇φ) 2 in (1) recovers the dynamical equation of Active Model B whenever the spurious drift term T ν vanishes. From symmetry arguments, this coupling term is the lowest order in gradients and in φ which cannot be integrated into a free energy [26,27,29,30]. A term of the form (∇φ)(∇ 2 φ) is potentially present at same order as ∇(∇φ) 2 [30][31][32], yet both terms are equivalent in one spatial dimension, as we consider below.
The spurious drift terms {T ν, T χ} appearing in dynamics (1-4) vanish when choosing a specific spatial discretization, as shown in Appendix B. Then, there is no need to actually modify the dynamical equation already used in [29,30] to embed Active Model B in a thermodynamically consistent framework. For a constant chemical potential difference ∆µ, the dynamics follows aṡ where we have set the mobility λ = 1, and {Λ, ξ} are zero-mean Gaussian white noises with correlations proportional to the temperature T , as given in (2) and (4). The free energy F captures a phase separation between dilute and dense regions: In what follows, most of our results are valid for a generic f , and the specific form (32) is used for explicit evaluation only.
The corresponding heat rate, as given in (14), readṡ The heat rate quantifies the irreversibility of the whole dynamics based on trajectories of the active current J and of the chemical currentṅ, see (9). The heat rate profileq(x) depends on the details of the dynamics via the parameters of the free-energy F, the driving coefficient ∆µ, and the temperature T which controls the amplitude of fluctuations. For strong fluctuations, namely high temperature T , we expect the heat rate to be uniformly dissipated in the system, with only a weak dependence on the details of the density profile. Conversely, in the regime of small T , the local heat rate should reveal the salient features of the density profile which require energy to be sustained.
To explore the connection between density profile and heat rate, we then rely on a small-noise treatment of the dynamics. Given that (33) is fully determined by the fluctuations of φ, independently of that of n, we focus on the dynamics of φ alone. Expanding the density field ) and substituting this ansatz in (31), the leading order equation yields the deterministic mean-field dynamics: Hence, φ 0 relaxes to a steady-state profile which can either be uniform or comprising phase-separated domains depending on free-energy parameters in (32), the global density (1/V )´V φ(x)dx and the driving parameter ∆µ [30,31]. At higher orders, φ 1 and φ 2 follow a set of coupled stochastic dynamics given bẏ where and Λ 0 is a zero-mean Gaussian white noise with corre- Owing to the linearity of D 1 in (35)(36), φ 1 has Gaussian fluctuations. Substituting the density ansatz in (33), we geṫ where The expressions (37)(38) give the leading orders of heat rate at small noise for an arbitrary f .
For a homogeneous profile (φ 0 (x) = cst), the leading and first orders of the small noise expansion (37-38) vanish (ε 0 = 0 and ε 1 = 0). Then, the non-trivial contribution to heat rateQ − γV ∆µ 2 scales like T 2 at small T , In the absence of phase separation, namely for a homogeneous average profile of density ( φ(x) = cst), the non-trivial contribution to heat rateQ − γV ∆µ 2 and the density field irreversibility T S scale like T 2 at small temperature T , and as ∆µ 2 at small driving parameter ∆µ, as shown respectively in (a,b) where solid lines are guide lines. Their difference also exhibits similar scalings in these regimes, as shown in (c,d), and it is in good agreement with our prediction ( and it also behaves like ∆µ 2 at small ∆µ, as confirmed by the numerical results in Figs. 2(a-b). Therefore, in the absence of any density pattern, one can makeQ−γV ∆µ 2 arbitrarily small by reducing either the amplitude of fluctuations T or the driving parameter ∆µ. In particular, at vanishing T , the uniform density profile is identical to that of Passive Model B, namely for ∆µ = 0, which explains whyQ − γV ∆µ 2 also vanishes.
Previous works quantified irreversibility based on trajectories of the active current J only [30,68], without considering that of the chemical currentṅ: as defined in (27). The irreversibility measure T S has similar scalings as that ofQ − γV ∆µ 2 at small T and small ∆µ, as shown in Figs. 2(a-b). Interestingly, while the heat rateQ converges to the finite value γV ∆µ 2 at zero T , the irreversibility measure T S vanishes in this limit: The former captures the consumption of underlying chemicals, whereas the latter only sees an effective equilibrium dynamics. Moreover, the difference between T S andQ − γV ∆µ 2 is ∆µ 2´V [∂ x (∂ x φ) 2 ] 2 dx according to (26). We compute analytically this difference in Ap- x pendix C to show that it also scales like T 2 and ∆µ 2 at small T and small ∆µ, respectively, as confirmed by our numerics in Figs. 2(c,d).
For a phase-separated profile (φ 0 (x) = cst), as shown in Fig. 3(a-b), the leading order ofQ − γV ∆µ 2 scales like T 0 , since now ε 0 = 0, and it reaches a finite value at T = 0. Hence, the heat rateQ is not only determined by the background term γV ∆µ 2 at zero temperature, it now also depends on the mean-field density profile. In contrast, T S scales like T and thus vanishes at T = 0, see Fig. 3(c), as already reported in [30]. The different scalings ofQ and T S in this regime reveal that the former is affected by the existence of two separated phases, whereas the latter does not allow one to distinguish the active phase separation from its passive counterpart. This clearly illustrates that the irreversibility shown by the active current J alone, when the underlying chemical fluxṅ is not monitored, cannot capture the full energetic cost of creating phase separation away from equilibrium. In other words, if one were to propose T S as a thermodynamically consistent measure of the full energetic cost, based on the explicit entropy production rate S which discards the driving field fluctuations, then a nonequilibrium phase separation could be sustained at zero cost, in contradiction with the basics of thermodynamics.
The heat profileq(x) given in (33) not only provides information about where heat is dissipated, it also quantifies how the average chemical current ṅ(x) varies in space. At small temperature, it is constant in the dense and dilute phases, where the density profile is flat, and it has a non-monotonic behavior across the interface, as shown in Fig. 3(a): The system dissipates less heat at the interface than in the bulk, and it does so by reducing locally the chemical current to accommodate for density gradients. Likewise, the profile σ(x) in (39) is flat in the bulk and varies strongly at the interface [30]. Yet, now both the bulk and interface contributions vanish at zero temperature, see Fig. 3(b), consistently with the fact that T S vanishes in this regime. Moreover, both T S andQ − γV ∆µ 2 scale as ∆µ 2 for small ∆µ, as shown in Fig. 3(d), similarly to the case of a homogeneous density profile: The energetic costQ and the irreversibility measure T S vanish at zero ∆µ, since Active Model B becomes Passive Model B in this regime.

B. Motile polar droplets
As a second example of how our formalism may be used, we now turn to study the coupled dynamics of a polar field p and a scalar field φ. In our case, they represent the local polarization and density of active components, respectively, for instance self-propelled particles with aligning interactions [18][19][20], and are even under time-reversal. Our aim is to capture the emergence of complex order beyond the case of a phase separation, already discussed in Sec. III A, by incorporating the possibility of observing a nonequilibrium polar order [33,34]. We focus on dynamics with one spatial dimension for simplicity. A minimal ingredient to allow for a nonequilibrium advection of the fields then consists in taking the coupling terms as C φ = φp and C p = −p∂ x p in (15). When the spurious drift terms {ν φ , ν p } vanish, one recovers the dynamical equations used to model actomyosin droplets in the absence of hydrodynamic flow, e.g. due to substrate friction, as detailed in [73] for instance. In general, such type of coupling terms appears naturally when coarse-graining the dynamics of aligning active agents [24,25], and they also follow from symmetry arguments [21].
To show that our framework is also applicable for odd driving, we now choose to treat constant chemical currents, namely ∆ φ =ṅ φ /γ φ and ∆ p =ṅ p /γ p . The spurious drift terms {ν φ , χ φ , ν p , χ p } in (15)(16)(17)(18)(19) all vanish when choosing an appropriate spatial discretization, as detailed in Appendix B. The dynamics are then given bẏ We have set the mobilities {λ φ , λ p , γ φ , γ p } all equal to 1, and {Λ φ , Λ p , ξ φ , ξ p } are zero-mean Gaussian white noises with correlations proportional to T , as given in (16) and (18)(19). Inspired by recent works [73,74], we take the free energy F which leads to the formation of motile and quiescent regions: where the coefficients {a,φ, A, K, κ} are all positive. At equilibrium (ṅ φ = 0 andṅ p = 0), the system undergoes a phase separation whenever the global density (1/V )´V φ(x)dx is positive and less than φ d = φ 1 + 1 + A/(2aφ 2 ) , yielding coexistence between the dilute isotropic phase {φ, p} = {0, 0} and the dense polar The associated heat rate (24) readṡ To explore how the heat rate behaves at small temperature T , we again expand the fields as φ = φ 0 + √ T φ 1 + O(T ) and p = p 0 + √ T p 1 + O(T ). The mean-field dynamics follows from (40-41) aṡ for Ω ∈ {φ, p}. The first correction to the mean-field profile readṡ in terms of In the absence of a droplet, namely for homogeneous average profiles of density ( φ(x) = cst) and polarization ( p(x) = cst), the nontrivial contribution to heat ratė Q − 2Vṅ 2 scales like T 2 at small temperature T , and asṅ 2 at small driving parameterṅ φ =ṅp ≡ṅ. The numerical measurements get closer to our analytical predictions (C18), shown respectively in markers and solid lines, when the lattice spacing ∆x decreases, as expected. Simulation details in Appendix D. Parameters: a = A = κ = K =φ = 1, T = 10 −3 , V = 64.
where Λ Ω,0 are zero-mean Gaussian white noises with correlations Λ Ω,0 (x, t)Λ Ω ,0 (x , t ) = 2δ ΩΩ δ(x−x )δ(t− t ). As for the expansion in Sec. III A, the active fields at first order {φ 1 , p 1 } have Gaussian statistics. The heat rate (43) can then be expanded in the form (37), where {ε 0 , ε 1 } now read As a result, (47) provides the leading orders of heat rate at small temperature for a given free-energy density f .
In the homogeneous state (φ 0 (x) = cst and p 0 (x) = cst), the mean-field contribution toQ − V (ṅ 2 φ +ṅ 2 p ) vanishes (ε 0 = 0), yet the first order correction provides a non-zero contribution (ε 1 = 0): The non-trivial contribution to heat rateQ − V (ṅ 2 φ +ṅ 2 p ) scales like T , in line with what was found in [68] and in contrast with the T 2 scaling for the conserved dynamics of the scalar field φ in Sec. III A. We compute analytically this contribution in terms of the dynamical parameters, as detailed in Appendix C. For simplicity, we choose the driving parameters in the dynamics of φ and p to be equal (ṅ φ =ṅ p ≡ṅ), in which caseQ − 2Vṅ 2 behaves aṡ n 2 : This scaling is confirmed by our numerical results in Fig. 4. Note that our analytical result (see Eq. (C18)) depends on the lattice spacing through an ultra-violet cutoff.
For a droplet state (φ 0 (x) = cst and p 0 (x) = cst), as shown in Fig. 5(a),Q − 2Vṅ 2 scales like T 0 since the leading order is now determined by ε 0 = 0. Analogously  5. (a-b) The average profiles of density φ(x) and polarization p(x) , reported in the co-moving frame of droplets moving towards x > 0, show that polarization is non-zero only within droplets and its profile has a front-tail asymmetry. For multiple droplets, each one of them moves at same velocity with a fixed separation distance. (c-d) The corresponding average profiles of heat rateq(x), given in (43), are negative at the interface and increase monotonically from tail to front of each droplet. (e) The nontrivial contribution to heat ratė Q − 2Vṅ 2 =´Vqdx increases with the driving parameteṙ n φ =ṅp ≡ṅ of the motile polar droplets. Above a critical value ofṅ, when the droplet splits into two droplets, we observe a discontinuity ofQ − 2Vṅ 2 associated with a hysteresis loop whose area increases with the speed at whichṅ varies. to the phase separated state in Sec. III A, such a scaling implies that the heat rateQ depends on the details of the density and polarization profiles even at vanishing temperature. Increasing the value ofṅ splits the droplet into several ones which move in the same direction with a fixed separating distance, as shown in Fig. 5(b). For one droplet, the heat profileq(x) in (43) is peaked with negative values at the droplet interface, and it increases continuously with positive values from tail to head, see Fig. 5(c). This behavior is qualitatively similar for two droplets, see Fig. 5(d). In contrast with the case of a purely scalar field theory in Sec. III A, the heat profile is now non-zero not only at the interface, but also in the dense phase: This stems from the density and polarization profiles of droplets being non-flat. Moreover, the fact thatq(x) can have both signs illustrates that the local heat rate can be either above or below the background dissipation 2ṅ 2 . Thereforeq(x) − 2ṅ 2 can potentially be negative locally, as long as the overall heat rateQ stays positive, which corresponds to extracting energy from the thermostat at specific locations.
Interestingly,Q−2Vṅ 2 as a function ofṅ displays a discontinuity when the number of droplets varies, see transition between one and two droplets reported in Fig. 5(e). This shows that the total heat rate is strongly affected by the transition between different patterns, hence it can potentially be regarded as a relevant observable to characterize transitions, in line with previous results in particlebased active dynamics [43,44]. Varyingṅ linearly in time, we observe a hysteretic behavior so that the area of the loop inQ − 2Vṅ 2 vsṅ space increases with the driving parameter velocity dṅ/dt. Moreover, Fig. 5(f) confirms thatQ − 2Vṅ 2 scales like T 0 at small noise, as predicted analytically.

IV. CONCLUSION
Building the thermodynamics of active matter is a major challenge of modern nonequilibrium statistical mechanics. By combining first-principles and phenomenological arguments, it aims at quantifying and predicting anomalous properties in terms of a few well-chosen observables [35-38, 43, 44, 50, 51]. Following this route, the irreversibility of active dynamics has recently attracted much attention, since it provides an unambiguous measure of the distance from equilibrium: It is quantified by the explicit entropy production rate (EPR) which compares forward and time-reversed realizations of the dynamics [30, 41, 42, 56-60, 62, 68].
At microscopic level, the particle-based EPR can be related to the amount of heat dissipated by the system, though this relation can be more intricate for active systems [60] than for thermal ones [54,55]. At the hydrodynamic level, the connection between heat and explicit EPR is generally lost, so that the physical motivation for evaluating the explicit EPR in active field theories is sometimes unclear. Indeed, the heat rate is proportional to the total EPR provided that the latter measures the irreversibility of all hydrodynamic fields [54,55]. In contrast, the explicit EPR, which focuses on the irreversibility of active fields alone and discards the fluctuations of underlying driving fields, captures only a partial contribution to the heat rate. In practice, evaluating the total heat rate is then a challenge of modeling properly the coupling between active and driving fields.
In this paper, we have shown that the heat rate can be decomposed into a background contribution, independent of the active field, and a non-trivial contribution that dictates how the emerging order affects the energy cost. Importantly, the latter can be deduced systematically from the active field dynamics alone, provided that the equations of motion for active and driving fields are thermodynamically consistent, namely that the connection to surrounding thermostat is properly taken into account [54,55]. To ensure such a connection, we have embedded active field theories within linear irreversible thermodynamics [69], inspired by previous works [68,[70][71][72]79]. It amounts to considering underlying degrees of freedom as the basis of the nonequilibrium drive of the dynamics. Yet, at variance with previous studies [68,[70][71][72]79], we now consider explicitly the fluctuations of these driving fields.
Thermodynamic consistency enforces some spurious drift terms in the dynamics that are proportional to noise amplitude. As such, our framework is distinct from the approach commonly followed to derive active field theories, based either on symmetry arguments [21,[29][30][31][32][33][34] or explicit coarse-graining of microscopic dynamics [24][25][26][27][28], since it enforces dynamical terms often neglected in these theories. In practice, while being conceptually important, the spurious drift terms can be made to vanish by judicious choices of spatial discretization. More generally, these terms do not affect the mean-field behavior of the system at vanishing noise, so that the emergent dynamics and structure are still consistent with the existing literature of active field theory in this regime.
Within our framework, the dynamics of the active fields can be read out independently of that of driving fields, so that the latter may be regarded as hidden degrees of freedom. Moreover, the spatial decomposition of heat rate can be evaluated in terms of active fields only. This supports the fact that the emerging dynamics and patterns of active fields alone provide direct access to spatial variations of heat rate, without need to measure the fluctuations of hidden degrees of freedom. Note that fast-relaxing variables are neglected by our hydrodynamic description, which potentially provide additional contributions to the heat rate. Yet, provided that there is indeed a clear time scale separation between the fluctuations of hydrodynamic fields, either active or driving fields, and that of other neglected variables, any additional contribution to heat rate only changes the constant background term: Thus, it does not affect the connection between active field patterns and spatial variations of heat rate.
To demonstrate the practical relevance of our approach, we have evaluated the heat rate in two popular active field theories: (i) the dynamics of a conserved scalar field which captures active phase separation [29,30], and (ii) the coupled dynamics of scalar and vector fields which describes the emergence of motile de-formable droplets [73,74]. Spatial decomposition has revealed that there is reduced heat rate at interfaces, and we have analyzed the leading order of heat rate at weak noise in relation with emerging patterns. For motile droplets, we have also shown that the heat rate undergoes a discontinuous transition when the droplets either split or merge.
Our work provides the relevant framework to study how heat rate relates to emerging patterns at hydrodynamic scale. The map of heat rate indicates which parts of the system mainly dissipate energy to sustain nonequilibrium fluctuations. At variance with the map deduced from explicit EPR, which was studied in [30], the heatrate profile allows one to decipher where the underlying degrees of freedom contribute to shape the emergent profile of active fields. Note that our framework relies on the assumption of linear deviation from equilibrium thermodynamics [69], which does not hold for all active systems. It would be interesting to consider chemical reactions beyond the linear assumption, for which a thermodynamically consistent framework has been proposed recently [98,99].
Among the active field theories encompassed by our framework, many of them describe living systems, for instance dense assemblies of cells [77,97] and swarms of bacteria [91,92]. While previous experimental works have already evaluated the dissipation of either isolated molecular motors [100,101], cilia and flagella [45], tracers in living cells [102,103], or in vitro cytoskeletal network [46], only little is known regarding where energy is dissipated in spatially extended living systems. Our work opens the door to establishing maps of heat rate in models of living matter, with a potential to relating high/low dissipation locations with specific biological functions. Moreover, our predictions for the overall heat rate could potentially be tested against experimental measurements of the energy dissipated by living systems, such as metabolic rates [104,105], using for instance calorimetric techniques [106].
From a broad perspective, our framework lays the groundwork to bridge the gap between the thermodynamics of microscopic and hydrodynamic active theories. The stochastic thermodynamics of particle-based active dynamics has already received much attention in the last few years [41,[56][57][58][59][60]62]. Using systematic coarsegraining procedures, some active field theories are derived from microscopic equations, yielding explicit correspondences between hydrodynamic kinetic coefficients and microscopic parameters [24][25][26][27]. Based on these theories, our work offers the opportunity to compare the predictions for the heat rate of particle-based dynamics and that of their hydrodynamic counterparts, as way to analyze critically the energetics of active models at different scales. Interestingly, a specific class of active models has considered explicitly the coupling between particle degrees of freedom and underlying chemical reactions, following the recipe of LIT [64,107]. It would be interesting to explore whether coarse-graining this microscopic model leads to our framework at the hydrodynamic level.
Moreover, one could examine the performances of work extraction for continuum models [108,109] and compare them with results obtained recently for particle-based engines [110][111][112]. Furthermore, changing heat rate by using dynamical bias, one could study dynamical phase transitions in hydrodynamic theories, and compare them with that reported in active particles [49][50][51][52][53]. One expects the collective states of particle-based models emerging at high/low heat rate to coincide, at least qualitatively, with instabilities of hydrodynamic models in the same regime. If not, one could potentially try and revise the hydrodynamic equations to find a better agreement with their microscopic counterparts. Following this route, our work not only opens the door to controlling heat rate in active field theories, it also potentially provides a way to constrain their formulation. This calls for deeper investigations and encourages further contributions to the thermodynamics of active matter. Consider a simple model of a conserved field dynamics as presented in Sec. II A and depicted in Fig. 1. This system can be described by the dynamics of three species: active particles (φ) that are only present in the active subsystem, fuel (n f ), and products of the fuel consumption by the active particles (n p ). The three corresponding continuity equations are: where J f,p = −D {f,p} · ∇µ {f,p} , and r is the rate of fuel consumption that is non-vanishing only within the active subsystem. Here, the chemical potentials of the fuel and products are defined as usual µ {f,p} = δF/δn {f,p} . Schematic of the active slab system. Small orange/green spheres are fuel/products and large hollow spheres are the active particles, which are confined to the surface and consume fuel at rate r. The products of the reaction diffuse out of the slab to the bulk system (reservoir of both fuel and products) while fuel molecules diffuse into the slab from the bulk to maintain constant chemical potential.
We continue by defining the chemical coordinates n = (n p − n f ) /2 and n t = (n p + n f ) /2 such that µ f,p = (δF/δn t ∓ δF/δn) /2, and the chemical potential difference is ∆µ = µ f − µ p = −δF/δn. When diffusion of fuel/products is fast enough compared to the rate of fuel consumption r and the dynamics of the active fields φ, µ {f,p} adjusts very fast compared with the active dynamics, so that it can be considered to be constant throughout the entire system. In such a caseṅ p = −ṅ f = r, n t = 0, andṅ within the active subsystem. Although the dynamics of the chemical coordinate n and the fuel/products are essentially the same, these fields are not equivalent. Specifically, the free-energy dependence on either fuel, products, or the chemical coordinates is generally different. Finally, we assume that the timescale in which a significant change in ∆µ occurs is very long compared to the timescales of interest. Then, the reservoirs of fuel/products can be regarded as having constant chemical potentials, and a constant chemical potential difference ∆µ is maintained throughout the small active subsystem (see also Fig. 1), which is the source of activity. This is what used in the main text Eq. (3). The construction described above essentially forms a non-equilibrium grand-canonical ensemble. Within this ensemble,ṅ must be thought of as being the rate of fuel consumption, while the connection to the free energy of the reservoirs is seemingly lost.
A concrete example for such an active subsystem, which is also a prototype, is a thin slab as depicted in Fig. 6. For instance, this would be an appropriate description for the experiment on light-activated selfpropelled colloids of Ref. [17]. In this geometry one may write D {f,p} = D {f,p} ⊥ê ⊥ê⊥ + D {f,p} I −ê ⊥ê⊥ whereê ⊥ refers to the direction perpendicular to the thin slab. Because the slab is thin, diffusion of particles in/out of the slab is much faster than within it, such that fuel/products do not flow within the slab, J f,p (J f,p ·ê ⊥ )ê ⊥ . Conservation of mass then dictates that J p ·ê ⊥ = −J f ·ê ⊥ (the active particles cannot leave the slab) so thaṫ When diffusion of fuel/products in/out of the slab is fast compared with the active dynamics we get Eq. (A2) within the slab. Note that in this example the diffusion of fuel/products within the slab does not need to be fast compared with the active dynamics. It is sufficient to have fast diffusion of fuel/products perpendicular to the slab. At times long enough that the fuel reservoir starts to become exhausted, one must consider the change of ∆µ due to fuel/products fluxes in/out of the active subsystem, as in Eq. (A3). On these timescales there should not be any steady-state heat production. This is evident from Eq. (13) after substituting ∆µ = −δF/δn,ṅ t = 0 and Eq. (A2), which givesQ = d F[n, n t ] /dt = 0.

Appendix B: Spurious drift terms
In this Appendix, we obtain the expression of the spurious drift terms {ν Ω , χ Ω } in (1-4) and (15)(16)(17)(18)(19). To this aim, we first derive the Fokker-Planck equations (FPEs) associated with the spatially-discretized dynamics. Then, we choose the spurious drift terms so that the Boltzmann distribution is the steady state solution of FPEs in the equilibrium regime. We focus on the onedimensional case for simplicity (d = 1), since the generalization to higher d is straightforward. To generalize the discussion to d dimensions, one only needs to use the d-dimensional version of the gradient matrix instead of the one-dimensional matrix used below.
The spatial discretization amounts to considering the variables {φ i (t), p i (t)}, where the indices i denote lattice coordinates, whose dynamics converge to that of {φ(x, t), p(x, t)} in the limit of small lattice constant ∆x, where x = i∆x. In particular, we introduce the gradient matrix A defined by A standard choice for A is given by A ij = (δ i,j−1 − δ i,j+1 )/(2∆x), though other spatial discretizations are possible. In what follows, we discuss the consequence of such a choice in the form of the spurious drift terms.

Conserved dynamics for scalar field
We first consider the conserved dynamics for a scalar field in (1)(2)(3)(4), and write them in the discretized form aṡ where ψ i = (δF/δφ)(x = i∆x), and the coupling term . . ) depends on φ and its gradients in general. The noise terms {Λ i , ξ i } are Gaussian with zero mean and correlations given by (B3) Given that the correlations between Λ i and ξ i depend on the variable φ i through the coupling term C i , one has to specify the temporal discretization scheme of (B2). In what follows, and in the main text, we choose Stratonovich convention, which allows one to use the standard rules of differential calculus [78]. This is particularly convenient when deriving the expression of the heat rateQ defined in (9).
The associated FPE for the probability density P ({φ i , n i }, t) can then be derived following standard methods as [78] where we have introduced the matrix M i defined by In the continuum limit of small ∆x, it follows using (B1) that (B4) converges to the standard functional FPE for the probability density P ([φ(x), n(x)], t) [22,113].
Importantly, by taking {ν i , χ i } as the stationary solution of (B4) is given by the Boltzmann distribution P s ∼ e −∆x F/T at equilibrium, namely when [ψ i , ∆µ i ] = [∂F/∂φ i , −∂F/∂n i ], as expected [83,85]. As a result, the expression of {ν i , χ i , L i } in (B3) and (B5) provide a systematic way to compute the spurious drift terms in terms of C i . When C i is independent of n i , as is assumed below, (B5) vanishes if C i only depends on φ i , namely when it is a local function of φ independent of its gradients. Besides, the extension of (B5) for d > 1 follows directly by substituting the d-dimensional version of the gradient matrix A.
When d = 1, the chain rule then leads to simplify (B5) as (B7) The matrix M i can be written as (B8) Substituting the expression of M i in (B7), we deduce that ν i always vanishes for any C i in d = 1, yet it can still potentially be non-zero in higher dimensions. Besides, we deduce the expression of χ i as (B9) To obtain Eq. (26) from Eq. (14), one has to evaluate i C i Λ i following standard stochastic calculus [78], which reads where we have used again that C i is independent of n i . From (B9-B10), it follows that the relation between the heat rateQ and the explicit entropy production rate S given in (26) holds whenever j A ij (∂C i /∂φ j ) = 0, which is the case considered in the main text. Let us focus on the specific coupling term C AMB = ∂ x (∂ x φ) 2 = 2(∂ x φ)∂ 2 xx φ corresponding to Active Model B, as considered in Sec. III A. This coupling term can be written using different discretization schemes, such as both of which converge to C AMB at small ∆x. A priori, one might expect the spurious drift terms to be independent of the discretization scheme, yet we now show that different discretizations yield different expressions for the spurious drift terms in general. For C (1) , we get where we have used A ij = −A ji . Taking A ij = (δ i,j−1 − δ i,j+1 )/(2∆x), we deduce [A 2 ] ik A ik = 0, so that (B12) is zero. Substituting (B12) in (B7), we conclude that there is no spurious drift associated with C (1) since both ν i and χ i vanish, yet this no longer holds when considering higher-order schemes for the gradient matrix A. For C (2) , we get where we used again A ij = −A ji . Given that A is anti-symmetric, any odd (even) power of A is antisymmetric (symmetric), so that [A 3 ] ii = 0 and [A 2 ] ii = 0. Then, (B13) is always non-zero for any form of the gradient matrix A. The examples in (B12-B13) illustrate that the choice of spatial discretization affects drastically the form of the spurious drift terms.
For the study of Active Model B presented in Sec. III A, we choose to discretize C AMB using C (1) . Since the corresponding spurious drift terms vanish, this choice allows us to embed Active Model B in a thermodynamically consistent framework without need to change the dynamical equations considered in [29,41].

Generalized field dynamics
We now consider the generalized dynamics for conserved and non-conserved fields in (15)(16)(17)(18)(19). When the driving parameter is a chemical potential difference (∆ Ω = ∆µ Ω ), the spurious drift terms follow from a straightforward extension of (B5) as The expression of {ν Ω , χ Ω , L Ω } in (B14-B15) can then be used to derive explicitly the spurious drift terms for given coupling terms C Ω . As discussed in Sec. B 1, the choice for spatial discretization of the gradient terms appearing in C Ω is crucial to determine the corresponding spurious drift terms: A judicious choice can potentially make {ν Ω , χ Ω } vanish. The case where the driving parameter represents a chemical current (∆ Ω =ṅ Ω /γ Ω ) deserves a more careful treatment, which we discuss now. For simplicity, we address the dynamics of a polar field p without any conserved scalar field φ, as given bẏ where h i = −(δF/δp)(x = i∆x), and {Λ i , ξ i } are Gaussian noises with zero mean and correlations Note that there is no longer any correlation between Λ i and ξ i in contrast with (B3). To derive the corresponding FPE for P ({p i , n i }, t), it is convenient to substitute the expression ofṅ i in the dynamics of p i , yieldinġ Here the noise term Γ i = Λ i + C i ξ i /γ is Gaussian with zero mean and correlations given by where the Onsager matrix K i reads Note that the dynamics can be written in a compact form, analogous to that given in (6) when the driving parameter is a chemical potential difference, as Thoughṅ is taken constant in (B21), we consider the case where both p i and n i are stochastic variables to obtain the expression of the spurious drift terms {ν i , χ i }. The FPE for P ({p i , n i }, t) then readṡ where we have introduced the matrix J i defined by J i J i = K i . Choosing the spurious drift terms enforces that the stationary solution of (B22) is P s ∼ e −∆x F/T when [h i , ∆µ i ] = [−∂F/∂h i , −∂F/∂n i ], as expected in equilibrium [83,85]. Then, (B19-B23) define {ν i , χ i } in terms of C i through J i when the driving parameter represents a chemical current (∆ p =ṅ/γ). This definition is in general different from that given in (B14-B15) when the driving parameter is a chemical potential (∆ p = ∆µ p ). For d = 1, the spurious drift terms {ν i , χ i } are proportional to ∂C i /∂p i . In particular, taking which both converge to C p = p∂ x p = (1/2)∂ x p 2 at small ∆x, as considered in Sec. III B, we get where we have used A ii = 0. It follows that {ν i , χ i } vanish for the coupling term C (2) p but not for C (1) p . As already noticed for the conserved dynamics of a scalar field, see (B12-B13), different spatial discretizations yield different spurious drift terms. For the motile polar droplets studied in Sec. III B, we take the coupling term to be discretized as C (2) p which leads to vanishing spurious drift terms. Note that also vanishes whenever ∂C p,i /∂p i = 0, in particular it does so for C (2) p . Moreover, in the case where the dynamics (B16) features an additional conserved scalar field with driving parameter proportional to chemical current, namely ∆ φ =ṅ φ /γ φ , the spurious drift terms follow straightforwardly by extending (B23). Indeed, since the FPE for P ({φ, p i , n φ,i , n p,i }, t) can be separated into two sectors, associated with derivatives given by either {∂/∂p i , ∂/∂n p,i } or { j A ij (∂/∂φ j ), ∂/∂n φ,i }, we get where (B28) The spurious drift terms vanish when C φ depends on φ only locally. In particular, this is the case for C φ = φp as considered in Sec. III B.
generic expression in terms of the driving parameter and its conjugated chemical field. Moreover, we derive explicitly the dependence ofQ on model parameters for active phase separation and motile polar droplets, to leading order in noise strength, as considered respectively in Secs. III A and III B.

Generalized field dynamics
We first consider a generalized dynamics for a conserved scalar field φ and a polar field p of the forṁ is Gaussian with zero mean and correlations given by and the spurious drift terms are given for d = 1 as where M is defined by MM = L. In contrast with (15)(16)(17)(18)(19), we now consider an arbitrary Onsager matrix L, with the only constraint that it should be positive semidefinite (det L ≥ 0).
Following [85,87], the path probability P ∼ e −A associated with (C1-C3) is defined by where, as a consequence of the Stratonovich discretization, the spurious drift terms do not appear in the expression (C4) [85]. Note that some terms which are even under time-reversal have not been written explicitly in (C4) since they are not relevant for deriving the heat rate. These terms could potentially be relevant if one or several of the order parameters were odd under time reversal. The time-reversed dynamic action A R follows from (C4) by changing the sign of [J,ṅ φ ,ṗ,ṅ p ]. From the definition in (9), the heat rate can be written aṡ (C6) where we have usedφ = −∇ · J. In steady state, using d F /dt = 0, we then deduce the final expression of heat rate given in (22).
To treat the case where the driving coefficients are odd (∆ Ω =ṅ Ω /γ Ω ), the first step is to substitute in the dynamics of {φ, p} the corresponding expressions forṅ Ω , as is done in Sec. B 2, for instance. Then, following the same procedure as in (C4-C6), it is straightforward to show that (22) also holds for odd driving coefficients.
The range of integration for the wavenumber integral in (C13) runs over [−π/∆x, π/∆x]. In practice, for the discretized version of Active Model B in (B2) taken with the gradient matrix A ij = (δ i,j−1 − δ i,j+1 )/(2∆x), the odd and even lattice sites decouple when the driving coefficient ∆µ is small. In this regime, the field dynamics effectively evolves with a lattice constant 2∆x, so that the appropriate range of integration is then [−π/(2∆x), π/(2∆x)]: This is the wavenumber integration range that we consider when reporting our analytic prediction (C13) in Fig. 2

Motile polar droplets
Finally, we compute the heat rate for the dynamics of motile polar droplets, as discussed in Sec. III B. The corresponding heat rate is given in (43). In the homogeneous state (φ 0 (x) = cst and p 0 (x) = cst), it can be expanded at small noise T , yielding at leading ordeṙ where we have eliminated some boundary terms.
The numerical scheme for Sec. III A is the same as in [30]. In particular, it corresponds to choosing the discretization C (1) i with A ij = (δ i,j−1 − δ i,j+1 )/(2∆x) in (B11), so that the spurious drift terms vanish. The time and spatial discretization constants are fixed to be ∆t = 0.01 and ∆x = 1. The numerical scheme for Sec. III B is as follows. We assume the fields {φ m i , p m i } live on-lattice with periodic boundary condition: whereas the current J m , Λ m p,i } are Gaussian random variables with zero mean and unit variance, independent for each i and m. Note that we have chosen discretization C (2) p,i from (B24) so that the spurious drift terms vanish.
The multiplying factor in the last term of (D2) and (D4) comes from the regularization of the delta function in (2). The conserved noise Λ m φ,i+ 1 2 lives off-lattice, which means that its values are specified only for halfinteger lattice sites, whereas the non-conserved noise Λ m p,i lives on-lattice. The numerical results above does not depend strongly on ∆t (we choose ∆t = 10 −3 to 10 −5 ), however it depends slightly on the spatial discretization ∆x, as shown in Fig. 7.
The average density and polarization profiles in Figs. 5(a,b) and the heat rate profiles in Figs. 5(c,d) are computed in the co-moving frame of the droplets, since the droplets are moving with constant velocity towards x > 0. First we compute the centre of mass of the droplet(s) R(t) =´xφ(x, t) dx ´φ (x, t) dx. The instantaneous density profile in the co-moving frame will then be φ(x, t) → φ(x − R(t), t). Fig. 5(a) is the longtime average of φ in the co-moving frame: φ(x) = lim t1→∞´t 1 0 φ(x−R(t), t) dt/t 1 . Similar procedure is performed for the average polarization p(x) (Fig. 5(b)) and heat rate profileq(x) (Figs. 5(c,d)).