Stochastic efficiency of thermodiffusion: an extended local equilibrium approach

TL;DR

This paper develops a stochastic model for coupled heat and matter transport at small scales, analyzing efficiency fluctuations and their distributions using fluctuating hydrodynamics, with implications for various thermodynamic systems.

Contribution

It introduces a mesoscopic stochastic model based on fluctuating hydrodynamics to study efficiency distributions in thermodiffusion, extending the local equilibrium approach.

Findings

Bimodal efficiency distributions can occur under certain conditions.

Efficiency fluctuations evolve distinctly in large-size and long-time limits.

Results are applicable to a broad class of systems with ratio-based efficiencies.

Abstract

The recently established connection between stochastic thermodynamics and fluctuating hydrodynamics is applied to a study of efficiencies in the coupled transport of heat and matter on a small scale. A stochastic model for a mesoscopic cell connected to two macroscopic reservoirs of heat and particles is developed, based on fluctuating hydrodynamics. Within this approach, the fluctuating separation and thermodynamic efficiencies are defined. The conditions required to observe bimodal distributions of these efficiencies are determined, and the evolution of these distributions is investigated in the large-size and in the long-time limits. The results obtained are not restricted to thermodiffusion and can be generalized to systems where efficiency is defined as a ratio of stochastic state variables or dissipation rates.

Equations

$${\rho }_{tot}{c}_v\frac{\partial T}{\partial t}$$

$$\frac{\partial \rho }{\partial t}$$

$${\text{J}}_q$$

$${\text{J}}_1$$

$${\text{J}}_q$$

$${\text{J}}_1$$

$$\kappa =\frac{L_{qq}}{T^2}$$

$$D_T=\frac{L_{q1}}{\rho T^2}$$

$$D=\frac{L_{11}}{T}{\mu }^{\prime }$$

$$\mu ={\mu }^{\circ }(T,p)+R_{s}T\ln \rho ,$$

$${\text{J}}_q$$

$${\text{J}}_1$$

$$\frac{dU}{dt}=-\sum _k{\int }_{S^k}{\text{J}}_q^k\cdot {\text{dS}}^k\equiv \sum _k{F}_q^k$$

$$\frac{dM}{dt}=-\sum _k{\int }_{S^k}{\text{J}}_1^k\cdot {\text{dS}}^k\equiv \sum _k{F}_1^k,$$

$$\nabla T=\frac{\Delta T}{l};\nabla \rho =\frac{\Delta \rho }{l}.$$

$$F_q^a$$

$$F_q^b$$

$$F_1^a$$

$$F_1^b$$

$${\kappa }^{\prime }=\frac{S}{l}\kappa ;;D_T^{\prime }=\frac{S}{l}{D}_T;D^{\prime }=\frac{S}{l}D.$$

$$F_q^k=-L_{qq}^{\prime }\frac{\left(T-T_k\right)}{T^2}-L_{q1}^{\prime }\frac{\mu -{\mu }_k}{T}$$

$$F_1^k=-L_{q1}^{\prime }\frac{\left(T-T_k\right)}{T^2}-L_{11}^{\prime }\frac{\mu -{\mu }_k}{T},$$

$$L_{qq}^{\prime }=\frac{S}{l}{L}_{qq};L_{q1}=\frac{S}{l}{L}_{q1};L_{11}^{\prime }=\frac{S}{l}{L}_{11}.$$

$$\frac{dT}{dt}$$

$$\frac{d\rho }{dt}$$

$$\gamma =\frac{S}{\Omega l}$$

$$\lambda =\frac{\kappa }{{\rho }_{tot}{c}_v}$$

$$\alpha =\frac{R_s}{{\rho }_{tot}{c}_v}.$$

$$\left(T_{s,1},{\rho }_{s,1}\right)=\left(\frac{T_a+T_b}{2},\frac{{\rho }_a+{\rho }_b}{2}\right)$$

$$T_{s,2-3}=\frac{\lambda \mp \sqrt{\lambda (\alpha ({\rho }_a+{\rho }_b)(2D-D_T(T_a+T_b))+\lambda )}}{\alpha D_T({\rho }_a+{\rho }_b)}$$

$${\rho }_{s,2-3}=\frac{D\left(\alpha ({\rho }_a+{\rho }_b)\left[2D-D_T(T_a+T_b)\right]+2\lambda \pm 2\sqrt{\lambda (\alpha ({\rho }_a+{\rho }_b)\left[2D-D_T(T_a+T_b)\right]+\lambda )}\right)}{\alpha {\left[2D-D_T(T_a+T_b)\right]}^2}.$$

$$T_b=\frac{\lambda +\alpha \left({\rho }_a+{\rho }_b\right)(2D-D_T{T}_a)}{D_T\alpha \left({\rho }_a+{\rho }_b\right)}$$

$$T_{s,2}$$

$${\rho }_{s,2}$$

$${\tilde{\text{J}}}_q$$

$${\tilde{\text{J}}}_1$$

$$\overline{{\xi }_q^i\left(t\right){\xi }_q^j\left(t^{\prime }\right)}$$

Full text

I Introduction

Thermodynamics was originally developed in the 19th century to understand conversion of energy in macroscopic systems such as pumps or steam engines. Over the years, along with progress made in the analysis and in the manipulation of micro- and nanoscales devices, the interest in extending thermodynamics to small systems has been growing. At such scales fluctuations of the state variables cannot be overlooked and need to be incorporated explicitly. The numerous efforts made in this direction have led to the emergence of a new field collectively referred to as stochastic thermodynamics Schlögl (1980); Jiu-li et al. (1984); Gaspard (2004); Sekimoto (2010); Seifert (2012); Van den Broeck (2013); Ciliberto et al. (2010).

One recurrent and important question in stochastic thermodynamics is how to quantify the dynamics of relevant thermodynamic quantities, starting from the mesoscopic-level rules of stochastic processes. Most approaches so far have been based on statistical definitions of the thermodynamic quantities (such as the Gibbs-Shannon expression for entropy) combined with evolution equations for the probability distributions of state variables. Recently, we proposed a novel framework which rests on an extension of classical non-equilibrium thermodynamics, and thus of the local equilibrium hypothesis, to derive stochastic differential equations describing the spatiotemporal evolution of thermodynamic quantities of interest De Decker et al. (2015). This formalism was already applied to chemical reactions De Decker et al. (2016) and Brownian motion Nicolis and De Decker (2017). In the present paper we use it to analyze the fluctuating thermodynamics of a model for thermodiffusion in small systems.

Thermodiffusion consists in a coupling of transport phenomena, whereby gradients of temperature can induce transport of matter and gradients of concentrations can lead to heat flow. Because of this coupling, it becomes for instance possible to transport mass from a diluted reservoir to a more concentrated one. A system where such phenomena occur can be seen as an engine due to its ability to convert heat in work done against a gradient. Depending on the aim of the device, different relevant measures of the efficiency of thermodiffusion can be defined. The two most common efficiencies are the separation efficiency, which is a ratio between the gradients of mass and temperature, and the thermodynamic efficiency, which quantifies the relative dissipation of energy needed to induce thermodiffusion.

Both these efficiencies share the common property of being ratios of thermodynamic quantities. Various ratios of this type were previously studied in theoretical and experimental works. Ratios of work and heat were for example observed to display bimodal distributions for quantum dots and for Brownian motors Proesmans and Van den Broeck (2015); Martínez et al. (2016). A theoretical analysis of the thermodynamic efficiency was carried out later in Polettini et al. (2015). Bimodality was there interpreted as a sign that the system can work in different distinct regimes, giving rise to two peaks in the distribution of efficiency.

In the present work, we use the extended local equilibrium hypothesis to study the statistical properties of efficiencies in a simple model for thermodiffusion in compartmentalized systems. We first present the model and discuss its properties in the macroscopic (mean-field) limit (Section II). In Section III we show how fluctuations can be incorporated thanks to the aforementioned extension of the local equilibrium assumption. Stochastic versions of the separation and thermodynamic efficiencies are then introduced and analyzed in Section IV. The distributions of these quantities are studied numerically, and analytically when possible. Special emphasis is put on the parametric conditions leading to the emergence of bimodal probability distributions. We also study in detail how the efficiencies behave in the thermodynamic limit, i.e., for long times and/or large volumes. Interestingly, the way the thermodynamic limit is taken can have a strong impact on the distributions and these effects are not the same for the separation and for the thermodynamic efficiencies. These results are summarized and put in perspective in Section V, where we also discuss possible extensions of our work.

II Deterministic Model

We consider in this work an ideal and non-reactive dilute binary solution, with no bulk velocity and no external force. In the framework of a spatially continuous description, the evolution equations for the local temperature $T$ of the solution and for the local mass density $\rho$ of the solute are given by the divergence of fluxes:

[TABLE]

with $c_{v}$ and $\rho_{tot}$ being the specific heat at constant volume and the mass density of the solution, respectively. In the linear regime of non-equilibrium thermodynamics, the heat flux $\textbf{J}_{q}$ and the flux of matter $\textbf{J}_{1}$ are both linear combinations of the thermodynamic forces involved:

[TABLE]

In these equations, $\mu$ is the chemical potential per unit mass and $\nabla_{T}$ stands for gradients taken at constant temperature. The coefficients before the thermodynamics forces are the Onsager coefficients and the off-diagonal terms obey the Onsager reciprocity relation $L_{q1}=L_{1q}$. The fluxes are usually expressed by relating the Onsager coefficients to experimentally accessible quantities de Groot and Mazur (2011),

[TABLE]

where

[TABLE]

is the thermal conductivity,

[TABLE]

is the thermal diffusion coefficient and

[TABLE]

is the diffusion coefficient. Note that in this formulation information is required on the partial derivative of the chemical potential of the solute, $\mu^{\prime}=\left(\frac{\partial\mu}{\partial\rho}\right)_{T,p,\rho_{solv}}$, where $p$ is the pressure and $\rho_{solv}$ is the mass density of the solvent. Since we consider an ideal and dilute solution, the chemical potential takes the form

[TABLE]

where $\mu^{\circ}$ is a standard chemical potential and $R_{s}$ is the specific gas constant (the gas constant divided by the molar mass of the solute). The expressions for the fluxes thus finally read:

[TABLE]

Our objective in this work is to analyze the properties of stochastic efficiencies for the simplest model where thermodiffusion can occur. We consider to this end a homogeneous (well-stirred) system of volume $\Omega$, which is connected to two homogeneous reservoirs from which heat and matter can flow, as depicted in Fig. 1. Each reservoir plays both the role of a chemostat and a thermostat.

This model could be representative of mesoscopic channels or pores, where thermodiffusion is already used to sort particles Lervik and Bresme (2014) or for DNA sequencing He et al. (2013). The evolution laws for the internal energy $U$ and for the mass $M$ of solute in the central subsystem are given by integrating the fluxes of heat and matter flowing across the surfaces separating the reservoirs and the system,

[TABLE]

where $k=(a,b)$ refers to the reservoir under consideration and $S^{k}$ is the surface separating the central system from a reservoir $k$. In what follows, we will consider that these surfaces all have the same area $S$. To evaluate these evolution equations explicitly, we suppose that all gradients can be written as differences between values in the subsystems involved divided by the length of the central system (see de Groot and Mazur (2011) Chapter XV),

[TABLE]

Inserting these expressions for the gradients in the fluxes (11) and (12) leads to the following expressions for the rates appearing in (13)-(14):

[TABLE]

Here, $T_{a}$ ($T_{b}$) and $\rho_{a}$ ($\rho_{b}$) are the constant temperature and density of reservoir $a$ ($b$). These equations include the rescaled transport coefficients

[TABLE]

Notice that in a way similar to what is found in the spatially continuous description, the rates can be expressed as linear combinations of (here discrete) thermodynamic forces de Groot and Mazur (2011):

[TABLE]

in which one can find rescaled version of the continuous Onsager coefficients

[TABLE]

Inserting the rates (16)-(19) in eqs. (13)-(14), and using the relations $dU=\Omega\,\rho_{tot}\,c_{v}\,dT$ and $dM=\Omega\,d\rho$, the evolution equations for the temperature and for the mass density of solute in the central system are finally obtained:

[TABLE]

We introduced in these equations the geometric factor

[TABLE]

which is the inverse of a specific area and the quantities

[TABLE]

and

[TABLE]

Because the binary mixture is ideal and dilute, the transport coefficients $\kappa$, $D$, $D_{T}$ and the density of the solution $\rho_{tot}$ do not depend on $\rho$. Furthermore, we will assume for simplicity that these parameters are independent of temperature. Non-linearities in eqs. (24)-(25) thus solely arise from switching from Onsager coefficients to transport coefficients, as well as from the derivative of the chemical potential (10). More complex situations could be encountered by considering non-constant transport coefficients, but we prefer to focus here on the simplest non-linear case.

Eqs. (24)-(25), which form our model, can be seen as the discrete analogues of the continuous equations (1) and (2). They form an autonomous 2-variable dynamical system which can admit up to three steady-state solutions. One of these solutions is

[TABLE]

and depends on external constraints only. Since it includes the equilibrium state as a special case, we will name it the thermodynamic branch. The two other solutions depend both on reservoir-related quantities and on kinetic parameters:

[TABLE]

Depending on the values of the different parameters involved, the system can thus present one or three physically relevant solutions.

A typical example of bifurcation diagram corresponding to this system is depicted in Fig. 2, where the temperature $T_{b}$ of reservoir $b$ is varied. At low $T_{b}$, the only stable solution is the thermodynamic branch. Increasing $T_{b}$, a transcritical bifurcation takes place at which the thermodynamic branch loses stability in favor of one of the two solutions (30)-(31). The chosen set of parameters is representative of typical aqueous solutions and we thus expect that this transition could be observed experimentally. Note however that the transcritical bifurcation only appears for moderate values of $T_{b}$ when the thermal conductivity is low: this would correspond to a system whose walls are made of a moderately insulating material. As $T_{b}$ is further increased, these two additional solutions merge at a limit point bifurcation whose coordinates are given by

[TABLE]

At this limit point, the degenerate steady-state values of $T$ and $\rho$ take simple forms:

[TABLE]

This limit point marks the end of the domain of validity of the model developed here. Beyond this point, the system does not admit any stable steady state and both the temperature and the density “explode”. The absence of physically meaningful states in this parametric region is more than probably due to the fact that for such large temperatures, the linear forms we used to estimate the fluxes are no longer relevant. A similar effect is also found for too large density gradients. We will thus restrict ourselves to cases where the limit point bifurcation has not been crossed.

III Stochastic model

We consider in this work mesoscopic systems where fluctuations are not negligible. In order to take them into account, stochastic differential equations for the mass of solute and for the internal energy of the system need to be derived. To do so, we will here rely on the recently introduced extended local equilibrium framework De Decker et al. (2015). In this approach, the local equilibrium hypothesis at the basis of classical non-equilibrium thermodynamics is supposed to hold even when the evolution equations for the state variables are stochastic. The statistical properties of the entropy production associated with simple diffusion processes De Decker et al. (2015), chemical reactions De Decker et al. (2016) and Brownian motion Nicolis and De Decker (2017) were already investigated in this framework.

Due to the presence of fluctuations, the heat and the mass fluxes contain an additional noise term and are fluctuating quantities,

[TABLE]

where the first terms in the right-hand sides are the systematic fluxes given by eqs. (11)-(12) and the noise terms are denoted by the vectors $\bm{\xi}$. The spectrum of fluctuations for the state variables can become intricate in non-equilibrium states, exhibiting for example long-range correlations. A large amount of literature has been devoted to elucidate the causes of such features. They were related to the presence of multiplicative noise, whose amplitude depends on the spatial position (inhomogeneous noise) Tremblay et al. (1981), on advection termsKirkpatrick et al. (1982); Ortiz de Zárate and Sengers (2004) or on state variables Pagonabarraga and Rubí (1994). Since the system introduced in this paper is homogeneous, well-stirred and with constant transport coefficients, we consider that Landau’s classical fluctuating hydrodynamics Lifshitz and Pitaevskii (1980); Ortiz de Zarate and Sengers (2006) holds. Consequently, the additive noise terms correspond to zero-mean Gaussian white noises whose amplitudes obey the fluctuation-dissipation relations Foch (1971); Ortiz de Zarate and Sengers (2006)

[TABLE]

in which $i$ and $j$ represent the different contributions in space and $k_{B}$ is Boltzmann’s constant. The lines over $\xi_{q}$, $\xi_{1}$ represent averaging over the noise. $T_{m}$ and $\rho_{m}$ are the temperature and the mass density evaluated along the macroscopic (deterministic) trajectories computed by solving the deterministic equations (24) - (25). The stochastic evolution equations for $U$ and $M$ are now given by

[TABLE]

The extended local equilibrium hypothesis has the consequence that the “equation of state” $dU=\Omega\,\rho_{tot}\,c_{v}\,dT$ holds even in the presence of fluctuations. Consequently, starting now from eqs. (40)-(41) the evolution equations of the intensive state variables read

[TABLE]

$R_{T}(t)$ and $R_{\rho}(t)$ are zero-mean Gaussian white noises obtained by combination of the fluctuating contributions to the fluxes:

[TABLE]

They are both zero-mean Gaussian variables, whose covariances read

[TABLE]

The noise terms being entirely characterized, the analytical distribution of deviations of state variables from their mean values, $\delta X_{i}=X_{i}-\mkern 1.5mu\overline{\mkern-1.5muX_{i}\mkern-1.5mu}\mkern 1.5mu$, can be computed in the weak-noise limit. To do so, we expand the evolution equations around the means:

[TABLE]

Keeping the first order terms and introducing the matrix $\beta$ with elements

[TABLE]

the evolution of deviations can be written in the form of coupled Ornstein-Uhlenbeck processes,

[TABLE]

which obey a bivariate Gaussian distribution:

[TABLE]

where $\bm{\delta X}$ is the vector containing the fluctuations of state variables, $\bm{\delta X}=\left(\delta T,\delta\rho\right)$, and $C$ is the equal-time covariance matrix of the state variables:

[TABLE]

The latter can be computed from the matrix of rates derivatives $\beta$ Risken and Frank (2013):

[TABLE]

where

[TABLE]

is the covariance matrix of the noise terms $R_{T},R_{\rho}$.

IV Efficiencies

The efficiency of a thermodiffusive process can be quantified in several different ways. The separation efficiency $\chi_{k}$ is the ratio between the gradient of density and the gradient of temperature with respect to a given reservoir $k$:

[TABLE]

Evaluating this quantity is especially relevant when thermodiffusion is used to concentrate or to deplete a solute with the help of temperature gradients. Alternatively, as mentioned before the system can also be seen as a sort of engine where a usually non-spontaneous transport process can arise thanks to coupling with another transport phenomenon. In this context, a relevant measure is the thermodynamic efficiency $\eta_{k}$,

[TABLE]

where $\sigma_{-}^{k}<0$ is the entropy production of the non-spontaneous process and $\sigma_{+}^{k}>0$ the entropy production of the spontaneous process, both related to the reservoir $k$. The thermodynamic efficiency quantifies the (relative) dissipative cost associated with the emergence of the desired transport phenomenon. In what follows, we will be interested in the study of the statistical properties of those two efficiencies in a non-equilibrium steady state in which the state variables $T$ and $\rho$ fluctuate around their stationary macroscopic values $T_{m}$ and $\rho_{m}$, as given by eq. (29) on the thermodynamic branch and eqs. (30)-(31) on the post-bifurcation branch. Under these conditions the probability of the state variables is expected to relax towards its stationary distribution. This relaxation is typically very rapid for the parameter values used in the manuscript. We thus rapidly reach a state where both the mean and the variance of $T$ and $\rho$ are stationary as well.

IV.1 Stochastic separation efficiency

Because we include fluctuations in our description of the system, the separation efficiency (57) is the ratio of two state variables obeying stochastic differential equations. In the weak-noise limit, the evolution equations (42)-(43) define two coupled Ornstein-Uhlenbeck processes and the variables are thus expected to be distributed according to Gaussian distributions. It has been shown in various cases that the probability distribution of similar ratios can exhibit bimodality Polettini et al. (2015); Rana et al. (2014). We investigate here whether the same phenomenon can be observed for the distribution of the separation efficiency. More precisely, we will be interested in the influence of external control parameters (such as the temperature or the density of the reservoirs), of the size of the system and of the duration of observation on the emergence of bimodality.

IV.1.1 Probability distribution of $\chi_{k}$

Inasmuch as $\chi_{k}$ can be seen as the ratio of two random Gaussian variables, its probability distribution $P(\chi_{k})$ is given by Polettini et al. (2015)

[TABLE]

where

[TABLE]

Bimodality is defined as the appearance of two peaks in the probability distribution $P(\chi_{k})$. To illustrate the conditions under which bimodality emerges, we consider a system whose parameters are the same as in Fig. 2. To make computations explicit, we also need to define the geometry of the system. We choose here a simplified cubic geometry, as depicted in Fig. 1, where $\Omega=S\,l$. This choice of parameters leads to a bifurcation in the system as depicted in Fig. 2. The temperature $T_{b}$ of reservoir $b$ plays the role of a control parameter.

Fig. 3 plots the ratio of the two peaks’ amplitudes for the distributions $P(\chi_{a})$ and $P(\chi_{b})$ given by Eq. (59), as a function of the control parameter $T_{b}$. Note that when the distribution is found to be unimodal, this ratio is set to zero. We observe (Fig. 3(a)) that $P(\chi_{a})$ is bimodal only in the vicinity of the deterministic thermal equilibrium condition, $T_{a}=T_{b}=T_{m}=295$ K. $P(\chi_{b})$ is also bimodal under these conditions, but moreover becomes so again after the transcritical bifurcation point has been crossed, more precisely when the mean temperature in the system and the temperature of reservoir $b$ become equal, $T_{b}=\mkern 1.5mu\overline{\mkern-1.5muT\mkern-1.5mu}\mkern 1.5mu=375.49$ K (see Fig. 3(b)). More generally speaking, additional simulations reveal that bimodality is observed for a given separation efficiency $\chi_{k}$ whenever the system is in thermal equilibrium with the corresponding reservoir. In view of the definition of efficiency (Eq. (57)) bimodality thus also corresponds to cases where the macroscopic separation efficiency becomes singular as it changes sign. Interestingly, bimodality of the distribution is observed while the system resides in a unique steady state, as shown in Fig. 2. Unique steady states can thus give rise to bimodality in separation ratios.

Numerical simulations of the full stochastic differential equations (42) and (43) were carried out with the Euler-Maruyama algorithm Gardiner (2009). Histograms over $N_{real}$ trajectories of a duration $t$ with a time step $dt$ were produced from the simulations. They reproduce accurately the predicted emergence and disappearance of bimodality around these points, as shown in Fig. 4. It should be noted that the numerical histograms are often very well reproduced by the analytical expression (59).

IV.1.2 Fluctuations and efficiency

A question that naturally follows is how fluctuations affect the overall efficiency of the thermodiffusive process, in particular for cases where bimodal distributions are found. A first difficulty that arises in this context is how to choose a statistical quantity associated to the fluctuating efficiency that can be compared to the macroscopic efficiency. Indeed, the probability distribution (59) has non-converging moments, which excludes the use of the mean of $\chi_{k}$ as a relevant quantity. This lack of convergence also appears with numerical simulations of the full stochastic differential equations. We will thus rather analyze how the most probable efficiency $\chi_{k}^{\star}$ compares with the macroscopic value.

Fig. 5 shows a comparison between the macroscopic value (green dotted curve) computed from the and the most probable value (black curve) of $\chi_{b}$ computed from the analytical expression (59). We focus here on parametric regions where bimodality is found, since in other cases the most probable value is always very close to the macroscopic one. As the system goes through the points corresponding to (partial) thermal equilibria, the macroscopic efficiency diverges and switches from positive to negative values, as noted above. The most probable efficiency, on the other hand, does not present such divergence. Its absolute value is always smaller than that of its macroscopic counterpart and the transition from positive to negative values, although still discontinuous, is nevertheless finite.

This result suggests that in the presence of fluctuations, the thermodiffusive system switches from a regime of positive efficiency to a regime of negative efficiency in a smooth way. To confirm this, we decompose the possible values taken by the efficiency into 4 different subcases (here explained for $\chi_{b}$):

• •

Case 1: $T<T_{b}$ and $\rho<\rho_{b}$ ($\chi_{b}$ is positive).

• •

Case 2: $T<T_{b}$ and $\rho>\rho_{b}$ ($\chi_{b}$ is negative).

• •

Case 3: $T>T_{b}$ and $\rho<\rho_{b}$ ($\chi_{b}$ is negative).

• •

Case 4: $T>T_{b}$ and $\rho>\rho_{b}$ ($\chi_{b}$ is positive).

This decomposition highlights the fact that, for example, a positive separation efficiency can hide two different physical realities: one where density and temperature gradients are both positive, and another one where these gradients are both negative. The probability $C_{i}$ for the system to be in any of the 4 above subcases can be calculated explicitly for different $T_{b}$ with the analytical expression for the probability distribution $P(\chi_{k})$, eq. (59) (see the Appendix for more details). The result is shown in Fig. 6. It can be seen for low $T_{b}$ that the small system smoothly changes from case 4 ($\chi_{b}>0$) to case 2 ($\chi_{b}<0$) in the vicinity of the first thermal equilibrium, where both states become equally probable. In a sense, one could say that fluctuations allow the system to anticipate the appearance of a different regime of separation. As $T_{b}$ is further increased, the system then goes to the case 1 and finally to case 3 in a similar fashion. It should be noted that the extent of the domain over which these smooth transitions take place is not the same for all transitions and is moreover size-dependent, since it tends to zero as $\Omega$ increases.

These results show that there is an intricate connection between bimodality in $\chi_{k}$ and transitions between the aforementioned regimes. The thermodiffusive system switches from a separation regime to another either when $\mkern 1.5mu\overline{\mkern-1.5mu\rho\mkern-1.5mu}\mkern 1.5mu=\rho_{k}$, or when $\mkern 1.5mu\overline{\mkern-1.5muT\mkern-1.5mu}\mkern 1.5mu=T_{k}$. Close to these transitions, the amplitude of fluctuations is sufficient to make the system work in both regimes at the same time. However, this constant switching will not always lead to the appearance of two peaks in the probability distribution of separation efficiency. Comparing Figs. 6 and 3(b), we observe that bimodality in $P(\chi_{k})$ is found for some of the transitions, but not all of them. Indeed, bimodality appears in $P(\chi_{b})$ only for transitions where the denominator of the macroscopic separation efficiency vanishes: for instance, the transition point at $T_{b}=343.03$ K, where the numerator $\mkern 1.5mu\overline{\mkern-1.5mu\rho\mkern-1.5mu}\mkern 1.5mu-\rho_{b}$ is zero (violet to orange curve), does not produce a bimodal distribution. To emphasize the importance of the denominator, we can check the presence of bimodality in the inverse of the separation ratio,

[TABLE]

Ratios of the peaks in $P(\chi_{b}^{-1})$ are plotted in Fig. 3, showing bimodality around the transition point which was absent from $P(\chi_{b})$. Moreover, the broad interval of $T_{b}$ where bimodality is found corresponds to the interval where the transition from case 1 to case 2 takes place. These differences between the distribution of the separation efficiency and its inverse were also observed in the distributions $P(\chi_{a})$ and $P(\chi_{a}^{-1})$. In order to detect conditions under which different separation regimes can coexist in a small thermodiffusive system, both the distribution of $\chi_{k}$ and of its inverse should thus be considered.

IV.1.3 Emergence of the thermodynamic limit

We now turn to the question of the thermodynamic limit of the system under consideration. The thermodynamic limit is usually defined either as the behavior of a given system in the limit of very long times, or as its behavior for very large systems.

In the limit of large systems, the separation efficiency tends to the expected macroscopic value as can be seen in Fig. 7, where the analytical distribution (59) is plotted for different values of $\Omega$. We considered here a situation where $P(\chi_{b})$ is initially bimodal. As $\Omega$ increases, one of the peaks vanishes while the other peak, which is initially centered close to zero, gradually moves so that the maximum of the probability distribution converges towards the macroscopic value of the separation efficiency. The evolution of $P(\chi_{b})$ displays however a different behavior in the long-time limit, as shown in Fig. 8: bimodality is conserved and the shape of the distribution is not impacted by the duration of trajectory. Two different separation efficiencies will be both highly probable, instead of one as in the large size limit.

The reason behind this difference lies in the way the system size and time affect the variance of the stochastic variables appearing in the separation efficiency. On one hand, and in view of the fluctuation-dissipation relations used in the stochastic approach developed in Section III, the variance of the state variables scales as the inverse of the system size. As $\Omega$ becomes larger, it becomes more and more difficult for the system to fluctuate between the different separation regimes. The system thus spends most of its time in a state whose efficiency remains close to its macroscopic value. On the other hand, since the state variables are in the weak noise limit coupled Ornstein-Uhlenbeck processes, their variance rapidly (exponentially) converges to a stationary value that does not depend on the duration of the trajectory. The probability to switch from a regime to another is thus not affected by the time of sampling and the bimodal shape of the probability distribution is conserved.

IV.2 Stochastic thermodynamic efficiency

Since the thermodynamic efficiency is defined as a ratio of entropy productions, it is important to first discuss how we define stochastic entropy production.

IV.2.1 Stochastic entropy production

In classical non-equilibrium thermodynamics, the total entropy production per unit volume (or in short, the dissipation) of a system can be expressed as the sum of the product of fluxes $F_{i}$ and associated thermodynamic forces $X_{i}$:

[TABLE]

In the case hereby considered, $\sigma$ can also be expressed as the sum of contributions related to heat transfer and to matter transport for each reservoir $k$:

[TABLE]

The different contributions read de Groot and Mazur (2011)

[TABLE]

where $F^{k}_{q}$ and $F_{1}^{k}$ are given by eqs. (16)-(19). In the framework of the extended local equilibrium approach, the contributions to entropy production have the same structure, with the difference that the fluxes now contain a fluctuating contribution,

[TABLE]

Consequently, each contribution to $\sigma$ can be split into a systematic part $\left(\sigma_{i}^{k}\right)_{sys}$

[TABLE]

and a part arising from the effect of fluctuations $\left(\sigma_{i}^{k}\right)_{fluct}$

[TABLE]

At the macroscopic scale, the thermodynamic efficiency $\eta_{k}$ associated to thermodiffusion with respect to a given reservoir $k$ is defined as the ratio of instantaneous values of the different contributions to entropy productions (see eq. (58)). In the stochastic approach that we adopt, these contributions are proportional to Gaussian white noises, which could result in the appearance of singularities. To avoid such situations, we will consider in what follows the time-averaged entropy productions

[TABLE]

Notice that the noise present in the different contributions $\sigma^{k}_{i}$ to the stochastic entropy production is not additive but multiplicative. This makes the derivation of general analytical expressions for $\sigma_{i}^{k}$ or for its time-averaged counterpart $\Sigma_{i}^{k}$ a rather difficult task. However, $\Sigma_{i}^{k}$ can be shown to obey a Gaussian distribution at the steady state and in the weak-noise limit. Indeed, developing the dissipation $\sigma$ around its deterministic non-equilibrium stationary state $\overline{\sigma}$, we obtain to the leading order in $\Omega^{-1}$:

[TABLE]

The two first terms are integrals of Ornstein-Uhlenbeck processes and define zero-mean normally distributed variables in the long-time limit. The last term is the time average of a Wiener process, which is also normally distributed and has a mean of zero. The sum of all these contributions will thus also be a Gaussian variable (see De Decker et al. (2016) for more details).

As an illustration of this, the histograms of the different contributions to

[TABLE]

are depicted in Fig. 9 for a given value of $T_{b}$. Despite the small size of the system considered in this example, each contribution defines an approximately Gaussian distribution centered around a well-defined mean, as can also be seen in Fig. 10. Some contributions are centered around negative values ($\Sigma_{1}^{a}$), others around positive values ($\Sigma_{q}^{a}$) and one of them remains close to zero ($\Sigma_{q}^{b}$).

IV.2.2 Probability distribution of the thermodynamic efficiency

We consider in this Section the statistical properties of the stochastic thermodynamic efficiency, which we thus need to define in more detail. As previously stated, the macroscopic thermodynamic efficiency is usually defined as the ratio of the dissipation of the non-spontaneous process over the dissipation of the spontaneous one. This means that the very definition of $\eta_{k}$ can change as the parameters are varied, since the system can switch from a regime to another (for example, from a heat pump regime to a matter pump regime). To avoid ambiguities in defining thermodynamic efficiency for each situation, we choose to denote by $\eta_{k}$ the ratio of mass transport dissipation over heat transport dissipation:

[TABLE]

As mentioned above, we rely on time-averaged dissipations along a stochastic trajectory of duration $t$:

[TABLE]

to avoid singularities. Note that we will be interested in the statistical properties of both $\eta_{k}$ and $\eta_{k}^{-1}$

The fact that $\Sigma^{k}_{1}$ and $\Sigma^{k}_{q}$ are distributed in an almost Gaussian fashion suggests that $\eta_{k}$ and/or $\eta_{k}^{-1}$, just like the separation efficiency $\chi_{k}$, could obey non-trivial probability distributions including bimodality. However, the lack of analytical results prevents using the same method as in Section IV.1 to detect the emergence of such bimodal distributions. We will instead rest to do so on the “dip test”, which was originally proposed in Hartigan and Hartigan (1985).

The dip is a statistical indicator measuring the departure of an empirical distribution from unimodality. The computed dip is the maximum difference between the empirical cumulative distribution function and the unimodal cumulative distribution function which minimizes that difference. When the distribution is unimodal, the dip vanishes for exponentially decreasing probability densities and tends to a small positive constant otherwise. When multimodality is present in the probability density, the dip increases significantly. Fig. 11 displays the computation of the dip in $P(\eta_{b})$ and $P\left(\eta_{b}^{-1}\right)$ for 1000 realizations of stochastic trajectories at different values of $T_{b}$. Notice that the weak fluctuations of the dip result from the finite size of the sampling. As previously observed in the analysis of the separation efficiency, multimodality appears for well-defined values of the control parameter. Multimodality appears in $P(\eta_{b})$ at $T_{b}=294.95$ K, $T_{b}=312.00$ K and $T_{b}=375.5$ K, and for $P\left(\eta_{b}^{-1}\right)$ at $T_{b}=312.85$ K and $T_{b}=343.15$ K. Computations of the histograms at those points confirm that these multimodalities indeed correspond to bimodalities.

In a way similar to what was observed for $\chi_{k}$, the emergence of bimodality can be correlated to situations where the mean of the denominator of $\eta_{k}$ (or $\eta_{k}^{-1}$) vanishes. At the macroscopic scale, situations where the denominator of either $\eta_{k}$ or its inverse vanish would correspond to changes in the spontaneity of one of the processes, since the entropy production associated to it changes sign for these conditions. For example, the bimodality appearing in $P(\eta_{b})$ at $T_{b}=294.95$ K and $T_{b}=375.5$ K corresponds to cases where $\mkern 1.5mu\overline{\mkern-1.5muT\mkern-1.5mu}\mkern 1.5mu=T_{b}$ so that to the leading order in $\Omega^{-1}$, $\overline{\Sigma_{q}^{b}}=0$ since the mean thermodynamic force is zero (see eq. (69)). At these points, both $\chi_{k}$ and $\eta_{k}$ thus obey bimodal distributions. However, while for the separation efficiency a vanishing mean of the denominator can only be achieved for $\mkern 1.5mu\overline{\mkern-1.5muT\mkern-1.5mu}\mkern 1.5mu=T_{k}$, more possibilities exist for cancelling the mean of the entropy production appearing in the thermodynamic efficiency. Indeed, following eq. (69), the mean entropy production due to heat dissipation also vanishes (again, to the leading order in $\Omega^{-1}$) whenever the mean thermodynamic flux becomes zero. Predicting the parametric conditions under which this happens is not trivial, but the cancelling of the fluxes can be observed numerically. For $\eta_{b}$, this leads to the appearance of an additional region of bimodality around $T_{b}=312.00$ K. A similar situation is observed for the distribution of $\eta_{b}^{-1}$: the peak at $T_{b}=343.15$ K corresponds to a vanishing mean force, and the peak at $T_{b}=312.85$ K corresponds to a vanishing mean flux. To emphasize the importance of the value on the denominator, note that for a same value of control parameters, the distribution of an efficiency will be bimodal but its inverse will not, even though the physical situation is identical.

IV.2.3 Size and sampling time effects on the distribution of thermodynamic efficiency

Despite the absence of analytical expressions for $P(\eta_{b})$ and $P(\eta_{b}^{-1})$, the large size and long time limits can be investigated through numerical simulations. Both distributions were computed at temperatures $T_{b}=375.6$ K and $T_{b}=343.31$ K, where bimodality was observed for short times and small systems. At $T_{b}=375.6$ K, the originally bimodal distribution $P(\eta_{b})$ converges in the limit of large systems to a unimodal distribution whose most probable value is the macroscopic value (see Fig. 12). Moreover, the amplitude of fluctuations around this mean decreases as the system becomes larger. The originally unimodal distribution $P(\eta_{b}^{-1})$ follows the same tendency, with a decreasing variability and a convergence to the macroscopic value. Similar conclusions can be made for $T_{b}=343.31$ K, where $P(\eta_{b}^{-1})$ was observed to be bimodal. Even though the most probable value first peaked on the right of the distribution in small systems, the peak converges to the macroscopic value in macroscopic systems. The situation is here analogous to what was observed for the effect of the system size on the distribution of the separation efficiencies.

The long-time limit shows interesting departures from the distributions taken at large size limit. In both cases bimodality, when present, tends to disappear in favor of a unimodal distribution (see Fig. 13). The amplitude of fluctuations also decreases, but contrary to the large system limit the most probable value for large times does not converge to the macroscopic value but to an entirely different one, having sometimes an opposite sign. This shift can be attributed to the non-linearities present in the stochastic entropy production. In the limit of large systems high order moments such as the average of $T^{k}$ or $\rho^{k}$ tend to factorize as $\overline{T}^{k}$ and $\overline{\rho}^{k}$, respectively, since the mean deviations $\delta T^{k}$ and $\delta\rho^{k}$ scale as the inverse of the system size. This factorization is not expected to hold for small systems, even in the limit of long times.

Remember that for the separation efficiency, bimodality is preserved even for long sampling times. This is no more the case for the thermodynamic efficiency, for which the bimodality is only transient. This different behavior can be rationalized on the basis of the time dependence of the variances of the variables entering the ratios that define these two efficiencies. Thermodynamic efficiencies are not ratios of Ornstein-Uhlenbeck processes, as is the case for the separation efficiency, but ratios of integrals of Ornstein-Uhlenbeck and Wiener processes (at least, in the weak noise limit). The main difference between these two situations lies in the temporal dependence of the variances. A stationary Ornstein-Uhlenbeck process is characterized by a constant variance, while that of its integral and of the integral of a Wiener process are linear in time. So in the long time limit the mean of the dissipations will remain constant but their standard deviations (STD) will decrease as the square root of $t$, since the integrated dissipations are divided by $t$. The integral does not affect, however, the size dependence of the different moments of the dissipations. The mean dissipation is an intensive quantity and the dominant contributions to the variance (as computed from Eq. (74)) scale as $\Omega^{-1}$. Consequently the coefficient of variation (CV) of the integrated dissipations,

[TABLE]

decreases with time and system size, as previously mentioned. In long time limit, this means that the fluctuations around the mean dissipations will be sufficiently small to prevent the system from switching from one regime to another. The distribution will then become unimodal. The time needed to observe a unimodal distribution depends on the amplitude of fluctuations and thus varies with the size of the system.

V Conclusions

In this work, we relied on a stochastic extension of the local equilibrium hypothesis to study the efficiency of thermodiffusion in small compartmentalized systems. Starting from stochastic differential equations for the mass density and for the temperature of the system, two different measures of the efficiency of thermodiffusion were analyzed: the separation efficiency, which is a ratio between state variables and the thermodynamic efficiency, which is a ratio between time-averaged entropy productions. In most cases, the probability distributions associated with these quantities are Gaussian and centered around their macroscopic value. However, both of them were seen to obey bimodal distributions for specific values of the control parameters whenever small systems and short times are considered. These bimodalities were shown to correspond to situations where the system constantly switches between two different regimes, including for example transitions between a heat pump and a matter pump regimes.

This phenomenon could be traced back to situations where the mean of the denominators appearing in these separations vanishes. In the case of separation efficiency, bimodality appears when the mean temperature of the central cell equals that of one of the reservoirs, i.e., when a partial thermal equilibrium is established. In the case of thermodynamic efficiency, bimodality coincides with situations where the mean dissipation associated to one of the thermodynamic fluxes vanishes. This can happen either because the thermodynamic driving forces becomes zero or because the flux itself vanishes.

Another interesting property emerging from numerical simulations is the fact that the limit of large systems and the limit of long times are not equivalent. For both efficiencies, the distribution becomes unimodal and peaked around the macroscopic value in the large size limit. In the limit of long times, bimodality is conserved for the separation efficiency while it disappears gradually for the thermodynamic efficiency. This can be explained by the fact that the variance of the average entropy production decreases with time, while that of the state variables does not. Identical macroscopic constraints and stationary distributions of state variables can thus lead to different dynamics for the distribution of an efficiency, depending on how it is defined.

While the present work focuses on thermodiffusion, a similar approach could be used to derive stochastic models for other systems where dissipative phenomena are coupled, such as thermoelectric devices or active membrane transport. Other aspects of coupled transport phenomena like the trade-off between power and efficiency could be studied with our extension of the local equilibrium hypothesis. It would be particularly interesting to find relations between this trade-off and the various transport coefficients and reservoir parameters, a task for which our formalism is well suited. Beyond the application to coupled transport phenomena, the results about the behavior of stochastic ratios are expected to be quite general because they are solely based on the local equilibrium assumption and should thus apply to a broader range of problems. As long as the state variables describing a system can be described by Ornstein-Uhlenbeck processes, ratios of these quantities should be expected to obey bimodal distributions whenever the mean of the denominator vanishes. This could be observed, for example, when computing the yield of a process admitting several pathways. Similarly, under these conditions the extended local equilibrium hypothesis has the consequence that entropy production can always be expanded as a sum of integrals of Ornstein-Uhlenbeck and Wiener processes. Ratios of dissipations, including for example the ratio of mechanical power and heat flux used in Carnot cycles or the ratio of electrical current and heat flux in fuel cells, will thus present features that are similar to those observed here.

VI Acknowledgements

J.-F. D. acknowledges the F.R.S.-FNRS (Fonds National de la Recherche Scientifique) for financial support through the grant of a FRIA fellowship.

Appendix: Probabilities of separation categories

In this section the computation of the probability for each separation case, $C_{i}$, will be discussed explicitly. Suppose that we look for the probability distribution $P(y)$ of a generic ratio $y=\frac{x_{1}}{x_{2}}$, having only the joint probability distribution $P_{j}(x_{1},x_{2})$. After a simple change of variables, the probability distribution of the ratio can be split in two parts following the sign of the denominator $x_{2}$ :

[TABLE]

Using the same reasoning as in Polettini et al. (2015), an analytical expression for those two contributions can be computed if the joint probability $P_{j}(x_{1},x_{2})$ is a multivariate Gaussian distribution.

Considering the particular case of the separation efficiency $\chi_{k}$, the joint probability distribution of $P(T-T_{k},\rho-\rho_{k})$ will be given by the same expression as the one derived in Section III. Indeed, even if the probability distribution were computed for the fluctuations of state variables around their means, a simple linear variable shift enables us to find the same distribution with the same covariance matrix.

The distributions following the sign of the discrete gradient $T-T_{k}$ can be then written:

[TABLE]

The probability to observe a given separation category is finally obtained by computing the proper integral on the chosen domain of the efficiency

[TABLE]

Those probabilities were computed for each value of the parameter $T_{b}$, as depicted in Fig. 6.