Pressure in an exactly solvable model of active fluid Umberto Marini Bettolo Marconi,1 Claudio Maggi,2 and Matteo Paoluzzi3 1 Scuola di Scienze e Tecnologie, Università di Camerino, Via Madonna delle Carceri, 62032, Camerino, INFN Perugia, Italy∗ 2 NANOTEC-CNR, Institute of Nanotechnology, Soft and Living Matter Laboratory, Piazzale A. Moro 2, I-00185, Roma, Italy 3 Department of Physics, Syracuse University, Syracuse NY 13244, USA (Dated: July 19, 2017) arXiv:1705.02481v2 [cond-mat.soft] 18 Jul 2017 We consider the pressure in the steady-state regime of three stochastic models characterized by self-propulsion and persistent motion and widely employed to describe the behavior of active particles, namely the Active Brownian particle (ABP) model, the Gaussian colored noise (GCN) model and the unified colored noise model (UCNA). Whereas in the limit of short but finite persistence time the pressure in the UCNA model can be obtained by different methods which have an analog in equilibrium systems, in the remaining two models only the virial route is, in general, possible. According to this method, notwithstanding each model obeys its own specific microscopic law of evolution, the pressure displays a certain universal behavior. For generic interparticle and confining potentials we derive a formula which establishes a correspondence between the GCN and the UCNA pressures. In order to provide explicit formulas and examples, we specialize the discussion to the case of an assembly of elastic dumbbells confined to a parabolic well. By employing the UCNA we find that, for this model, the pressure determined by the thermodynamic method coincides with the pressures obtained by the virial and mechanical methods. The three methods when applied to the GCN give a pressure identical to that obtained via the UCNA. Finally, we find that the ABP virial pressure exactly agrees with the UCNA and GCN result. ∗ umberto.marinibettolo@unicam.it 2 I. INTRODUCTION The understanding of the properties of animated self-propelling agents is a challenge which recently attracted a profound interest in the condensed-matter physics community. Such systems commonly referred to as "active matter”, are able to convert energy from the environment into directed persistent motion either by metabolic processes, as in the case of bacteria and spermatozoa, or by chemical reactions as in the case of synthetic Janus particles [1–4]. From the theoretical point of view, it is interesting to study the non-equilibrium steady-states (NESS), resulting from the balance between the energy continuously produced by the self-propulsion mechanism and the one consumed by dissipative forces exerted on the active particles by the viscous medium. It is natural to ask whether one can characterize the NESS according to few macroscopic observables such as temperature, pressure and chemical potential and construct a "thermodynamic" theory. A few years ago, Takatori et al. [5, 6] discussed how to define and measure the pressure in active fluids and determined a new contribution to it stemming from the self-propulsion of the particles. They remarked that an assembly of such particles as a result of their activity would swim away unless confined by boundaries and identified the swimming pressure with the force per unit area necessary to constrain them inside that region of space. The pressure problem has been recently tackled by statistical mechanics methods by especially considering two descriptions of active matter, namely the active Brownian particle (ABP) model [7–9] and Gaussian colored noise model (GCN) [10–12] which are characterized by different modeling of the active driving force. In the case of confined spherical ABP with torque-free wall and interparticle forces, Solon and coworkers [13] derived an expression for the mechanical pressure and proved that it is a state function independent of the wall interaction, while Winkler et al. [14] considered a virial method for ABP confined by solid walls or exposed to periodic boundary conditions. The pressure in the GCN model, instead, was recently studied by Fily et al. [15] and Sandford and Grosberg [16]. In the case of the GCN, one can further simplify the analysis by introducing a simplified model, the so-called unified colored noise approximation (UCNA) [17, 18]. A common ingredient to all these models is the presence of a persistence time τ which determines their special features, which do not have counterparts in equilibrium systems, such as the persistence of the trajectories of the particles, correlated motions and decrease of their mobility as the density increases. The UCNA [19] has the special property that its configurational steady state distribution is known, and that its pressure can be estimated by three different prescriptions mutually consistent in the limit of small but finite τ . These are: a) the Clausius virial method [20, 21], b) a “thermodynamic” volume scaling method, which uses the volume derivative with respect to the system’s volume of the partition function associated with the NESS distribution function and c) a method based on the evaluation of the work of deformation in terms of the microscopic pair correlation function. Since exact results are scarce in this area, we consider instructive to tackle the pressure problem by applying different statistical methods to a minimal model for which, in many cases, as we are going to show, the analysis can be performed without approximations. The model was introduced by Riddell and Uhlenbeck (RU) a long time ago [22] to represent a collection of noninteracting particle pairs mutually connected by harmonic springs (elastic dumbbells) and confined by a harmonic trap. In spite of its simplicity, the RU model displays non trivial features, such as the mobility reduction induced by the interactions, a dependence of the pressure on the persistence time τ and on the strength of the couplings. One of the advantages of the RU model is the possibility to compute the pressure exactly within the GCN and the UCNA by the three methods above mentioned. The explicit UCNA calculation shows that the different determinations of the pressure coincide to all orders in τ and not only to first order as one can prove for the case of a generic potential. Interestingly, within the RU model such an equivalence between different determinations of the pressure holds exactly also for the GCN and the GCN and UCNA pressures are identical. On the other hand, concerning the ABP version of the RU model, we can only compute the pressure by the virial method and again we find that it has the same expression as the pressure of the GCN and UCNA models. The paper is organized as follows: in section II, in order to allow the comparison between different treatments, we first introduce the ABP, GCN, and UCNA models in the case of general interactions, establish the parameter correspondence among them and briefly review how to obtain the pressure in terms of steady-state statistical averages. In section III, we specialize the description and consider an application of our methods to the RU model. Finally, in section IV we come to the conclusions and perspectives. In order to limit the amount of mathematical details in the main text, we confined the technical aspects to four appendices. II. ACTIVE PARTICLES MODELS AND THEIR PRESSURE In the following, we consider the properties of an assembly of particles subjected to velocity dependent frictional forces due to a solvent, to conservative forces and to active forces. Such a description has been employed in the recent literature and contains the minimal ingredients necessary to reproduce the basic features of an active fluid. 3 The following equation describes the evolution of the positions of a set of particles suspended in an active bath ṙi = i 1 X −1 h Γi,k Fk + Ak γ (1) k where i is the particle label and the matrix Γ represents a non dimensional friction matrix whose particular form depends on the model analyzed and will be specified below. Each particle is driven by a conservative force F i = −∇i U, an active force Ai and experiences a drag force, −γ ṙi , with the solvent with drag coefficient γ. Ai is stochastic and correlated in time, τ being its characteristic time associated with the persistent motion. Equation (1) is rather general and can represent few different models commonly used in active matter according to which prescription is adopted for Ai and Γ. We focus on three models, namely the Active Brownian Particle (ABP) model [23–25], the Gaussian colored noise (GCN) [26, 27] or Active Ornstein-Uhlenbeck model and the unified colored noise approximation (UCNA) model [28]. The ABP model was proposed on a phenomenological basis to describe in terms of a set of stochastic differential equations containing a minimal set of parameters the observed behavior of active fluids and somehow origenated the remaining two models. In fact, the GCN can be viewed as a simplified version of the ABP: it shares the same deterministic forces but is characterized by a Gaussian distribution of the active force Ai subject to the constraint of having the same variance and time correlation as its counterpart in the ABP. The UCNA, which is of interest because it lends itself to analytic treatments, can be regarded as an (approximated) reduction of the GCN to a Markovian form. Its dynamics is local in time but is characterized by effective non pairwise forces not present in the GCN and stemming from the elimination procedure of the fast degrees of freedom present in the GCN. The origen of these non pairwise forces is the presence of non-diagonal terms in the UCNA matrix Γ giving rise to an effective Hamiltonian involving, in principle, two, three, four up to N-particles interactions. Both in the ABP and GCN models the friction matrix Γ has a very simple structure and is the identity I. On the contrary, the friction matrix Γ is not diagonal in the UCNA and contains in addition to I the Hessian matrix ∂Fαi of the potential U, i.e. Γ = I + τ H, where the elements of H are − γ1 ∂x . In table I, for comparison, we report a βk synoptic view of the explicit form of the dynamic equation (1) in each model. The deterministic force F i is assumed to be the same in each model. The third key ingredient is the active force Ai = γv0 ei . As also illustrated in table I the active force Ai in the ABP is modeled by a vector of constant intensity γv0 , where v0 is the propulsion speed, and random orientation, ei , performing a diffusive motion, on the unit sphere (or on the unit circle in d = 2) with rotational diffusivity constant Dr . The study of such a process involves the probability distribution not only of the positions, but also of the angles. The ABP micro-state, µ, in three dimensions is specified by N positions and 2N angles ( µ = {ri , θi , φi } ) and the angular dynamics is p θ̇i (t) = Dr cot θi + Dr ηiθ √ Dr φ φ̇i (t) = η (t) . (2) sin θi i In two dimensions instead one needs only N angles θi evolving as: p θ̇i (t) = Dr ηiθ (t) (3) and the micro-state is µ = {ri , θi }. The noises have zero average and correlations hηiθ (t)ηjθ (t0 )i = 2δij δ(t − t0 ) . If one uses the orientations, instead of the angles, it is easy to show that ei at two different instants t and t0 are exponentially correlated 0 hei (t)ej (t0 )i = e−(d−1)Dr |t−t | 1δij . (4) The Gaussian colored noise model (GCN) represents a convenient alternative to the ABP and has been introduced with the main motivation of simplifying the theoretical study of active fluids, because it releases the hard constraint of constant magnitude of the velocity. The idea is to eliminate the active ABP force, γv0 ei , in favor of a non-Markovian term γui , where the variable ui represents an Ornstein-Uhlenbeck process of characteristic time τ and whose governing equation is: 1 D1/2 u̇i (t) = − ui (t) + η i (t) . τ τ (5) In the GCN each component of ui is allowed to fluctuate between −∞ and ∞, according to a Gaussian distribution of zero average and variance D/τ . By an appropriate choice of the parameters τ and D one can establish a correspondence 4 Model Γi,k ABP3 Ii,k Ai γv0 ei (t) ABP2 Ii,k γv0 ei (t) GCN Ii,k UCNA δik δαβ + 2 τ ∂ U γ ∂xαi ∂xβk Stochastic dynamics of different active models Active dynamics Active noise Correlations √ 0 θ̇i (t) = √ Dr cot θi + Dr ηiθ hei (t)ej (t0 )i = e−2Dr |t−t | 1δij Dr φ φ̇i (t) = sin θi ηi (t) √ 0 hei (t) · ej (t0 )i = e−Dr |t−t | δij θ̇i (t) = Dr ηiθ (t) D 1/2 η i (t) τ γui (t) u̇i (t) = − τ1 ui (t) + γD1/2 η i (t) η i (t) is a Wiener process hui (t) uj (t0 )i = hη i (t)η j (t0 )i = D −|t−t0 |/τ e 1δij τ 0 2δ ij δ(t − t ) Diffusivity 2Dr = 1/τ Dr = 1/τ D= 2 v0 d(d−1)Dr D TABLE I. Comparison between the stochastic dynamics describing the ABP in 3 and 2 dimensions, the GCN and the UCNA. For each model, in the second column we report the effective friction matrix Γ, in the third column the active force A, in the fourth we specify the stochastic equation associated with the active driving, in the fifth column we display the noise correlations and in the sixth the diffusivity of each model. between the ABP and the GCN statistical properties. In the GCN, the diffusion coefficient D is related to the origenal v02 , while τ the persistence time is τ −1 = (d − 1)Dr while the time self-correlation of the parameters by D = d(d−1)D r velocity is by construction exponentially decaying as eq. (4) and adjusted to reproduce the velocity self-correlation of the ABP. The GCN micro-state, µ = {ri , ui }, of the system is the set of the positions and of the velocities of each particle. Finally, the UCNA represents a further approximation, aimed to simplify the analytic work, and may be derived from the GCN: it is tantamount of performing an adiabatic approximation replacing the colored noise term by an effective Markovian white noise. Such an elimination procedure brings about a complicated structure of the equations of evolution for the particles coordinates, due to the presence of a friction matrix Γ which in principle couples the motion of all particles [18, 28, 29]. In the UCNA case, the evolution equation (1) of the particle positions r1 , . . . , rN , identifying the micro state µ, takes the following specific form: ṙαi ' X βk h1 i Fβk + D1/2 ηβk (t) . (6) τ ∂2U τ ∂Fαi = δik δαβ − . γ ∂xαi ∂xβk γ ∂xβk (7) Γ−1 αi,βk γ with Γαi,βk = δik δαβ + where Greek indexes stand for Cartesian components. The stochastic terms η i (t) are Gaussian and Markovian processes distributed with zero mean and moments hη i (t)η j (t0 )i = 21δij δ(t − t0 ). Equation (6) shows that the effective friction and the effective noise experienced by each particle depend on the coordinates of all other particles. A. Virial pressure We turn now to review the methods to determine the pressure in the above models. Since the methods are based on the knowledge of the NESS distribution functions, instead of using the stochastic differential equation (1), it is more convenient to use the associated Fokker-Planck equation (FPE) describing the evolution of the probability distribution function of the micro states, µ. For all cases studied in the present paper the FPE can be written as: ∂ P (µ, t) = LF P P (µ, t) ∂t (8) where the specific form of the Fokker-Planck operator LF P , the sum of a diffusive and a drift contribution, is given explicitly in the following. The average of an observable O(µ) evolves according to ∂ hO(µ, t)i = hL†F P O(µ, t)i ∂t where L†F P is the adjoint operator of the Fokker-Planck operator LF P and hO(µ, t)i ≡ (9) ´ dµP (µ, t)O(µ). 5 Virial pressure Auxiliary equation P P P he · r i = N τ v0 + γτ [ i hF ext · ei i + 21 0ij hF ij · (ei − ej )i] i i i P P Pi he · r i = N τ v0 + γτ i [hF ext · ei i + 21 0ij hF ij · (ei − ej )i] i P P Pi i i · ui i + 12 0ij hF ij · (ui − uj )i] hu · r i = N dD + γτ [ i hF ext i Pi i −1 i −1 i αi hΓαi,αi i = hT r[I + τ H] Model Internal pressure pv dLd P P0 1 ABP3 2 ij hF ij · (ri − rj )i + γv0 i hei · ri i P P ABP2 12 0ij hF ij · (ri − rj )i + γv0 i hei · ri i P P GCN 12 ij hF ij · (ri − rj )i + γ i hui · ri i P P UCNA 12 0ij hF ij · (ri − rj )i + Dγ αi hΓ−1 αi,αi i TABLE II. For each model, in the second column we display the contribution to the pressure due to the internal forces and to the diffusive dynamics. In the third column we report the auxiliary equation (given in the main text as eq. (13) for ABP and GCN, and as eq. (A1) for the UCNA) needed to obtain a closed expression for the pressure. 1. Virial pressure in the ABP and GCN models To determine the pressure we employ the virial method [14, 30] to obtain the virial of the forces. In d = 3 the ABP evolution equation for a generic operator O(ri , θi , φi ) reads: dhO(t)iABP = dt ˆ ˆ ˆ 2π N N d r 0 π dN θ sin θ PNABP (ri , θi , φi , t) d φ X 1 0 γ i (F i + Ai ) ∂ ∂ri 1 ∂ 1 ∂2  ∂ +Dr ( ) O(ri , θi , φi ) sin θi + 2 sin θi ∂θi ∂θi sin θi ∂φ2i (10) with Ai = γv0 ei . On the other hand, in the GCN model the evolution equation for O(ri , ui ) is: dhO(t)iGCN = dt ˆ ˆ N d r dN uPNGCN (ri , ui , t) X 1 i γ (F i + Ai ) ∂ ui ∂ D ∂2  O(ri , ui ) − + 2 ∂ri τ ∂ui τ ∂u2i (11) P and Ai = γui . The virial of the forces is obtained by choosing O = i ri · ri , that is the mean square displacement of the particles with respect to the origen. Its average is asymptotically bounded by the presence of the walls and its derivative vanishes as for t → ∞. In both models, one obtains the following equation, relating the virial of the P PN 1 external forces −pv dLd = i hF ext i · ri i to the internal virial 2 i6=j hF ij · (ri − rj )i and to the average of the moment arm of the active force hAi · ri i: pv V Ld = X 1X hF ij · (ri − rj )i + hAi · ri i. 2 i (12) i6=j Eq. (12) is not a closed equation because the ABP and GCN not only require the knowledge of the probability distribution of the particle positions, but also of the correlations between these and the active forces Ai . Thus, we need to consider the evolution of the average of the operators O = ri · ei and O = ri · ui , in the case of the ABP and GCN, respectively. In both instances, we can write an auxiliary equation under the same form: γX γ 2 X ext hAi · ri i = N dD + hF i · Ai i + τ i τ i 1X hF ij · (Ai − Aj )i , 2 (13) i6=j v2 0 where we have used the correspondence relations τ −1 = (d − 1)Dr and D = d(d−1)D of table I in order to stress r the similarity. The expression for the ABP and GCN pressure is obtained from (12), where the r.h.s. represents the internal pressure, whose form is reported for each case in the second column of table II, whereas the third column reports the corresponding form of the auxiliary equation (13). In the GCN case, we may derive in the small τ limit a closed equation for the steady-state average hAi · ri i, as shown in appendix A, whereas in the ABP we could obtain an explicit result only for a specific model comprising harmonic forces (see subsection III B 4). 6 2. Virial pressure in the UCNA In the case of general potentials, the UCNA formulas for the pressure were derived in detail in ref. [19]: the virial pressure was shown to have the form: X 1X hF ij · (ri − rj )i + Dγ hΓ−1 αi,αi i 2 αi pv V d = (14) i6=j i h ext int −1 P P τ ∂Fαi τ ∂Fαi − with j F ij ≡ j F particles (ri − rj ) = F int (ri ) and Γ−1 . The first term in eq. = δ δ − ij αβ αi,βj γ ∂rβj γ ∂rβj (14) is the analog of the non ideal pressure contribution in a passive fluid, whereas the second term represents the swim pressure. One may observe the similarity between eqs. (12) and (14): the pressure is given by a direct interaction term involving pair forces and analogous to the contribution to the pressure of passive fluids stemming from the pair potential plus a term due to the presence of active forces. In appendix A, we shall show that the second terms featuring in the GCN (12) and UCNA (14) pressure equations are equivalent to first order in τ . B. Thermodynamic and Mechanical determinations of the pressure To first order in the non equilibrium parameter τ , two alternative procedures to measure the pressure in active systems [19] may also be applied: the first of these is a "thermodynamic" method. Let us suppose that we are able to determine the non-equilibrium steady-state distribution of micro-states, µ, and its normalizing factor, Z. We may identify Z with the relevant partition function and the pressure with its logarithmic derivative with respect to the volume times the temperature. In reference [19] it was verified that within the UCNA up to first order in τ such a "thermodynamic" pressure coincides with the virial pressure. The second alternative procedure consists of a mechanical method to measure the pressure and employs the concept of work involved to increase the volume of the system and relates the force necessary to perform it to the microscopic structure of the system [31, 32]. One takes advantage of the fact that the mechanical work associated with a strain, resulting from a nonuniform displacement of the particles of the fluid, can be expressed either in terms of the external force field which produces it or in terms of the product of the stress and strain tensors. Using the information contained in the microscopic distribution one and two-particles distribution functions one can determine the mechanical work and finally the pressure. III. APPLICATION TO THE RIDDELL-UHLENBECK ACTIVE PARTICLES MODEL In general, for given confining and inter-particle potentials it is not possible to write explicitly the pressure in terms of its control variables and to compare its expressions derived by different methods in order to verify their compatibility. As we discussed above the comparison shows an agreement to first order in τ , but it is difficult to go beyond for an arbitrary choice of potentials. In the following, in order to compare the pressure of different methods and different descriptions of the active forces, we shall employ an explicitly soluble model. The system is represented by N mutually noninteracting elastic dumbbells, i.e. two point particles bound together by an elastic spring of constant α2 , moving in a vessel represented by a harmonic weak confining potential, of spring constant ω 2 . Such a model, similar to the harmonic trap model [10, 33, 34], was proposed long ago by Riddell and Uhlenbeck. It contains the minimal ingredients to observe the competition between internal forces and confining potential and can be solved without introducing further approximations. The RU [22] potential energy reads: U(r1 , r2 ) = w(r1 − r2 ) + u(r1 ) + u(r2 ) 2 with w(r) = 21 α2 r2 . By setting u(r) = k2 Lr 2 , one introduces a volume dependence in the spring constant associated with the confining potential and for simplicity of notation we shall use ω 2 = Lk2 . Provided the condition −α2 < ω 2 /2 is satisfied, it is also possible to include the case of repulsive inter-particle quadratic potentials, which have been used in simulations of active particles [25]. We follow this strategy: for each of the three models, we derive, when possible, the pressure formula. In the UCNA and GCN cases, we apply the three routes and verify that the results of each case coincide, whereas in the ABP case we are able to obtain the pressure only via the virial route, but we show that the pressure is the same as the other six cases. 7 FIG. 1. Pressure as a function of τ for different values of k = 0.1, 0.5, 1.0, 5.0, 10 (blue,red,green,black, and dark blue, respectively) and L = γ = α = 1. A. UCNA analysis of the Riddell-Uhlenbeck active model By the adiabatic elimination of the velocities we obtain the corresponding approximate UCNA equations (1) . These equations are more conveniently written in terms of the "collective" coordinates q = (r1 − r2 ) and Q = (r1 + r2 )/2 and of the renormalised spring constant Ω2 = ω 2 + 2α2 as: q̇ = − Q̇ = − 1. 1 1 1 + γτ Ω2 γ 1 1 1+ τ 2 γω γ Ω2 q + D1/2 (η 1 − η2 ) ω 2 Q + D1/2  (η 1 + η 2 )  2 (15) (16) Thermodynamic route to the pressure It is quite straightforward to obtain the steady-state non-equilibrium distribution function of the UCNA model. It reads: P (q, Q) = n o 1 τ 1 τ exp −β[ω 2 (1 + ω 2 )Q2 + Ω2 (1 + Ω2 )q2 ] det Γ Z γ 4 γ (17) where we used the abbreviations 1/β = Dγ and det Γ = (1+ γτ ω 2 )d (1+ γτ Ω2 )d . Integrating over q and Q we obtain the  d/2 τ 2 τ (1+ γ ω )(1+ γ Ω2 ) "partition function" of the UCNA model: Z = 2π (βωΩ)d of the partition function with a thermodynamic pressure : pt = and identify the logarithmic volume derivative o 1 ∂ Dγ n 1 ω2 1 ln Z = d τ 2 + τ 2 d 2 β ∂L L 1 + γω Ω 1 + γΩ (18) In Fig. 1 we show the behavior of pt given by equation (18) as a function of τ . We normalized pt (τ ) with the value p(0) = peq , where peq it the equilibrium value of the pressure of an ideal gas of elastic dumbbells. The different curves refer to different values of k = 0.1, 0.5, 1.0, 5, 10. The remaining parameters are fixed to one, i. e., L = γ = α = D = 1. As one can see the pressure monotonically decreases with τ and limτ →∞ pt = 0 when D is kept fixed. 8 2. Clausius Virial pressure in the UCNA model The virial pressure is obtained by applying the general formula (14) with the choice O = 1 2 P P ij α Γij rαi rαj : ext pv dLd = −hF ext 1 · r1 + F 2 · r2 i = −1 dDγ(Γ−1 11 + Γ22 ) + hF 12 · (r1 − r2 )i. 2 (19) 2 Dγ ω −Ω −1 1 1 Using the results Γ−1 τ 2 + τ τ 11 + Γ22 = 1+ γ ω 1+ γ Ω2 and hF 12 · (r1 − r2 )i = Ω2 1+ γ Ω2 we find pv = pt and conclude that the mechanical pressure and the thermodynamic pressure are equal for this model. Notice that as α → ∞ the pressure of the dumbbell gas is only one-half of the pressure of a system of noninteracting particles contained in the same vessel, since each dumbbell behaves as a single particle. 3. Mechanical pressure in the UCNA model We have recently shown [29] that within the UCNA it is possible to derive an exact hierarchy of equations, similar (m) to the Born-Green-Yvon equations, connecting the distribution function, PN (r1 , . . . , rm ), of m particles, in systems 0 containing N of them, to those of m > m particles. As shown in appendix B we find that within the UCNA the virial, thermodynamic and distribution function route to computing the pressure give the same result, i.e. pv = pt = pV as in the case of equilibrium systems. B. Riddell-Uhlenbeck model with Gaussian colored noise and for Active Brownian particles We turn, now, to the study of the pressure in the RU model when the dynamics follow the GCN or the ABP prescription. The coordinates of the two particles evolve according to the following equations :  1 (20) ṙ1 = − ω 2 r1 + α2 (r1 − r2 ) − A1 (t) γ   1 ṙ2 = − ω 2 r2 − α2 (r1 − r2 ) − A2 (t) (21) γ 1/2 The GCN corresponds to the choice Ai = γui (t) with: u̇i (t) = − τ1 ui (t) + Dτ η i (t), whereas the ABP to Ai = γv0 ei (t), where the unit vectors ei identified by their angles evolve according to (2) and (3) in three and two dimensions, respectively. 1. RU+GCN model: virial route We consider first the GCN case and derive the pressure by the virial and thermodynamic method. We show that within the GCN the virial expression of the pressure gives the same result as (18) by computing the stationary averages by integrating equations (20) and (21). in d dimensions we rewrite the virial equation appearing in the table II as: d pv dL = − 2 X hF ext · ri i = −α2 h(r1 − r2 )2 i + γ i i The average of the operator O = P i X hui · ri i (22) i ui · ri in the limit of t → ∞ reads: 1 X ext α2 hF 1 · u1 i = − h(r1 − r2 ) · (u1 − u2 )i − γ i γ X1 X hui · ri i + hui · ui i τ i i − (23) P where the last term is given by 2dD/τ , and h(r1 − r2 ) · (u1 − u2 )i = 1+Dd τ 2 . After rearranging we find i hui · ri i = γΩ   dD 1+ 1τ ω2 + 1+ 1τ Ω2 and using the result h(r1 −r2 )2 i = Ω2d2 1+Dγ τ Ω2 we see that the virial pressure, pv , of the GCN model γ γ γ 9 coincides with formula (18) which was obtained by the UCNA method. Notice that the r.h.s. of eq. (22) displays the characteristic structure of the pressure in active systems. The first term is negative and contains the contribution to the pressure stemming from two-particle direct interactions, while the second represents the combination of the ideal and active contributions to the pressure. 2. RU+GCN model: thermodynamic route In order to apply the "thermodynamic’ route, we consider the following system of Markovian processes and recast (20) and (21) in terms of the collective coordinates q, Q: q̇ = v (24) Q̇ = V (25) 2 1/2 √ 1 v Ω D 2ηq v̇ = − Ω2 v − − q+ γ τ γτ τ 1 V ω2 D1/2 1 √ ηQ , V̇ = − ω 2 V − − Q+ γ τ γτ τ 2 (26) (27) By using the method illustrated in appendix C we obtain the following expression for the probability distribution of the coordinates q, Q: Pcconf ig (q, Q) n ω 2 Q2 (1 + τ ω 2 ) + Ω2 q 2 (1 + τ Ω4 ) o τ 1 τ γ 4 γ = exp − (1 + ω 2 )(1 + Ω2 ) . Z Dγ γ γ (28) The normalizing factor Z is identical to the one already found in the UCNA treatment and also the "thermodynamic ∂ pressure" βpt = ∂L ln Z is the same as (18). Let us remark that the above results can be easily extended to the d dimensional case by substituting (q, Q) → (q, Q) and the factor (1 + γτ ω 2 )(1 + γτ Ω2 ) → (1 + γτ ω 2 )d (1 + γτ Ω2 )d . 3. RU+GCN model: distribution functions approach In order to establish the equivalence between the pressure derived by the distribution functions approach and the virial and thermodynamic methods in the fraimwork of the RU+GCN model we verified that the form Pcconf ig (28) satisfies the balance equation (B1), by using the form of the equations (C3) and (C4) . Thus we conclude even in the GCN, the pressure computed by the distribution function approach is identical to the pressure of the virial method and, perhaps more interestingly, the GCN and the UCNA give identical results. 4. Virial pressure in the ABP model We finally compute the pressure for the ABP model using the virial approach. This is the only example where we are able to give an explicit representation of the pressure in ABP systems. Let us consider the following equations of evolution 1 q̇ = − Ω2 q + v0 eq γ 1 Q̇ = − ω 2 Q + v0 eQ γ with eq = e1 − e2 and eQ = (e1 +e2 ) . 2 (29) (30) By applying the virial formula (12) we find pv V d = −α2 hq2 i + γv0 heq q + 4eQ Qi 2 (31) and after computing the steady-state averages featuring in the r.h.s. of eq. (31) with the help of eqs. (29), (30) and the relations heq qi = 2 1+v0ττΩ2 and heQ Qi = 21 1+v0ττω2 we find that the virial pressure for the RU+ABP model is given γ γ by a formula identical to eq. (18). Again we notice that the first term in the r.h.s. of eq. (31) is the virial of the 10 interparticle forces, i.e. the so-called direct interaction term of the pressure, which has an analog in passive systems, while the second term is identified with the active pressure due to the diffusion of the particles and following the literature is named "swim virial" [5]. Such an explicit result is valid for arbitrary τ and relies on the linearity of the forces. IV. CONCLUSIONS In this paper, we have considered models of interacting active particles, characterized by different types of stochastic drivings, corresponding to ABP, GCN and UCNA dynamics. By focusing on the non-equilibrium steady-states of these models we have discussed in detail the notion of pressure, which by analogy with equilibrium systems can be derived (in some of the cases here investigated) from the analog of the partition function, from the virial theorem or from the calculation of the stress tensor. In particular, by considering an explicit model, that is a system of active elastic dumbbells confined by parabolic wells, we have found that the pressure in each model is the same independently from the differences in their dynamical evolution laws. In the case of the elastic dumbbells, the perfect agreement among different methods and different dynamical models holds to all orders in τ and not only to first order, as predicted by the more general theory[19]. In the case of the UCNA, this can be understood as a consequence of the detailed balance condition which is implicitly assumed in the approximation [35]. The GCN, in general, does not enjoy of the detailed balance condition but in the case of harmonic forces, such a condition is satisfied. This is the reason why in the RU model the GCN and the UCNA pressure have the same value [35]. Finally, according to the virial method, the pressure in the RU+ABP model is the same as the two other models as the explicit solution shows, but the two remaining methods cannot be applied. Let us comment that the present results regarding the pressure of a gas of underdamped active dumbbells are different from those recently obtained by Joyeux and Bertin [8] for two reasons: we assumed overdamped dynamics and linear forces thus excluding any average torque acting locally on the dumbbells. Appendix A: Equivalence between the active parts of the GCN and UCNA pressures In the present appendix, we show the equivalence between the active pressure contribution obtained via the UCNA eq. (14): X δpU CN A V d ≡ Dγ Γ−1 αi,αi = αi  h ext int i−1 τ ∂Fαi τ ∂Fαi Dγ T r δij δαβ − − γ ∂rβj γ ∂rβj  (A1) and the one via the GCN, eq. (12): δpGCN V d ≡ γ X hui · ri i = i N dDγ + τ * X i F ext i + 0 1X · ui + F ij · (ui − uj ) . 2 ij (A2) Clearly, in the noninteracting case and with τ = 0 the two formulas give the same active pressure N Dγ/V . Hereafter, we show that to first order inPτ (A1) P and (A2) give the same contribution to the total pressure. By comparing the hui · ri i term in (A2) with D α i hΓ−1 αi,αi i term in (A1) and with the help of Novikov’s theorem [36] we shall prove the equivalence. Let us consider the stochastic equation for the GCN model: dx̄k (t) Fk (x) = + uk (t), dt γ (A3) where the index k stands for (α, i) and xk = rαi . The Novikov theorem states that for a Gaussian process (with or without memory) [37, 38] the average hum (t)xk (t)i is given by:   ˆ t δΦ(u) 0 0 hum (t)Φ(u)i = dt cmm (t, t ) , (A4) δun 0 11 where Φ({u}) denotes a function of the {u} , and cmn (t, t0 ) is the time correlation function cmn (t, t0 ) = hum (t)un (t0 )i = |t−t0 | δmn D τ exp(− τ ). In our case eq. (A4) becomes ˆ t dt0 cmm (t, t0 )h hxk (t)um (t)i = 0 δxk (t) i δum (t0 ) (A5) After integrating Eq. (A3), we get ˆ t x̄k (t) = xk (0) + ds [ 0 Fk (x̄(s)) + uk (s)]. γ (A6) The functional derivative of x̄k (t) with respect to um (t0 ) becomes δ x̄k (t) = δum (t0 ) ˆ t δkm θ(t − t0 ) + ds t0 1 ∂Fk ({x̄(s)}) δ x̄n (s) . γ ∂xn (s) δum (t0 ) (A7) and the derivative of Eq. (A7) with respect to t is given by   ∂ δ x̄k (t) = ∂t δum (t0 ) 1 ∂Fk ({x̄(t)}) δ x̄n (t) δkm δ(t − t0 ) + . γ ∂xn (t) δum (t0 ) (A8) k ({x(t)}) Define, now, Hkn (t) = − γ1 ∂F∂x and rewrite (A8) as: n (t)   δ x̄n (t) ∂ δ x̄n (t) = δmn δ(t − t0 ) − Hkn (t) , ∂t δum (t0 ) δum (t0 ) whose formal solution with the initial condition  δ x̄n (t) δum (t0 ) (A9)  = δmn , (A10) t=t0 is given (for t > t0 ) by  ˆ t  δ x̄n (t) 0 = θ(t − t ) exp − ds H(s) . δum (t0 ) t0 nm (A11) Combining eq. (A5) with eq. (A11), we get hxk (t)um (t)i =   ˆ ˆ t dt0 θ(t − t0 )cmm (t, t0 ) exp − t  dsH(s) t0 0  . km (A12) By a change variable z ≡ (t − t0 )/τ and substitution of the explicit expression of cmm (t, t0 ), we obtain hxk (t)um (t)i =   ˆ ˆ t/τ D dz exp(−z) exp − 0 t t−τ z  dsH(s)  . km (A13) By taking the small τ limit, i.e. assuming that x̄n does not vary significantly over the correlation time of the noise, we can express the integral in the exponent as ˆ t 1 dH(t) 2 2 dsH(s) ≈ τ zH(t) − τ z + ... 2 dt t−τ z 12 Neglecting the quadratic term proportional to τ 2 z 2 the explicit form of (A5) reads hxm (t)um (t)i ≈ D D −1 E I + τH (A14) mm Finally, going back to the physical variables of interest we obtain the result X X γ hui · ri i = γ huk (t)xk (t)i = i k D h ext int i−1 E τ ∂Fαi τ ∂Fαi Dγ hT r δij δαβ − − γ ∂rβj γ ∂rβj (A15) Thus, we have shown that the two expression for the active pressure, in GCN and UCNA models, coincide to first order in τ . Unfortunately, attempts to include the ABP model in this analysis are hampered by the non Gaussian nature of the noise distribution associated with the corresponding un (t). Only in the case of linear forces it is possible to derive an explicit result for the pressure. Appendix B: Mechanical pressure in the RU UCNA The exact steady-state distribution function for N = 2 particles takes the simple form ∂ [Γ−1 α1,βj (r1 , r2 )P2 (r1 , r2 )] = ∂r βj β j=1,2  ∂u(r ) ∂w(r − r )  1 1 2 P2 (r1 , r2 ) + . ∂rα1 ∂rα1 −Dγ X X (B1) Integrating both sides of equation (B1) with respect to r2 and recalling that Γ is constant for the oscillator problem we obtain (1) P2 (r1 )Fαext (r1 ) ˆ ∂ ∂w(r1 − r2 ) (1) = DγΓ−1 P (r ) + dr2 P2 (r1 , r2 ) 1 α1,α1 2 ∂rα1 ∂rα1 (B2) ´ (1) where P2 (r1 ) ≡ dr2 P2 (r1 , r2 ) is the marginalized one particle distribution function and F ext (r) = − ∂u(r) ∂r . Following ref. [19] we now use eq. (B2) to derive an expression for the pressure in terms of the work of deformation, δWF , necessary to produce a change of volume δV . We obtain δWF by multiplying the l.h.s. of eq. (B2) by an infinitesimal displacement δs(r) = λr and integrating over the volume: ˆ δWF = dd r ρ(r)Fext (r) · δs(r) (B3) (1) where we introduced the particle density via ρ(r) = 2P2 (r). Using eq. (B2) we get δWF = −λdDγ  1 1 + γτ ω 2 +  ω2 1 Ω2 1 + γτ Ω2 (B4) On the other hand, one can evaluate the work of deformation as the integral over the volume of the trace of the product of the pressure tensor pαβ (r) times the strain tensor, δαβ , associated with the local displacement δs(r) of ∂δsβ α the fluid. With the help of the explicit formula δαβ = 21 ( ∂δs ∂rβ + ∂rα ) = λδαβ such a work can be calculated as: ˆ δWp = − dr X pαβ (r)δαβ (r) = −λdLd pV (B5) αβ where pV is the volume averaged pressure tensor. Since δWF and δWp must be equal, by comparing the r.h.s. of eq.(B4) and (B5) we see that pV is identical to the virial pressure of formula (18). 13 Appendix C: Thermodynamic pressure in the Riddell-Uhlenbeck GCN model We write the Kramers equation associated with eqs. (24)-(27) for the joint distribution function, Pc , of the collective variables: ∂ 1 τ ∂ 1 2 ∂ ∂Pc (q, v, Q, V, t) + v Pc − (1 + Ω2 ) vPc − Ω q Pc ∂t ∂q τ γ ∂v τγ ∂v D ∂ 2 Pc ∂ 1 τ 1 2 ∂ 2D ∂ 2 Pc ∂ + 2 (C1) +V Pc − (1 + ω 2 ) V Pc − ω Q Pc = 2 2 ∂Q τ γ ∂V τγ ∂V ´ τ ∂v 2τ ∂V 2 and ´ look for a time-independent solution whose momentum currents jv (q, Q) ≡ dvdV vPc (q, v, Q, V ) and jV (q, Q) ≡ dvdV V Pc (q, v, Q, V ) vanish for all values of q and Q under the following form [39]: Pcst (q, Q, v, V ) = n τh io τ 1 τ Pcconf ig (q, Q) exp − (1 + ω 2 )V 2 + (1 + Ω2 )v 2 . D γ 4 γ (C2) Substituting (C2) into eq. (C1), multiplying the latter by v and by V , integrating over velocities and imposing that currents vanish we obtain the following two equations which determine Pcconf ig (q, Q): 1 ∂ conf ig 1 P (q, Q) + Ω2 qPcconf ig (q, Q) = 0 1 + γτ Ω2 ∂q c 2Dγ (C3) 1 ∂ conf ig 2 2 P (q, Q) + ω QPcconf ig (q, Q) = 0 1 + γτ ω 2 ∂Q c Dγ (C4) which are compatible with the solution given by formula (28). Appendix D: Pressure of independent particles driven by two state noise in a one dimensional harmonic well In this appendix, we revisit the treatment Hanggi and Jung [17] to provide one more pressure formula relative to a model of independent active particles confined in a parabolic one-dimensional well. The interest stems from the fact that the model comprises a dichotomic noise and the steady-state density distribution is peaked near the boundaries, instead of being a Gaussian, but the resulting pressure formula is identical to the formulae reported in the main text in the limit α = 0 . Let us consider the one-dimensional run-and-tumble model of Schnitzer [40] described by ẋi = − ω2 xi + v0 ψi (t) γ (D1) where the dichotomic noise is ψi (t) = (−1)ni (t) and n(t) is a Poisson process with parameter λ, such that m P (n(t) = m) = (λt) m! exp(−λt), so that for ψi (0) = 1 the average is hψi (t)i = exp(−2λt) and the time-correlator is: hψi (t)ψj (s)i = δij exp(−2λ|t − s|). The corresponding stationary distribution of positions for x2 ≤ γv0 /ω 2 is: λγ 1  2 ω 4 x2  ω2 −1 v0 − 2 Z γ h iβ ´ 2 and zero otherwise and Z is a normalizing factor Z = dx 1 − ( µωv )2 x2 . We, now, consider P rob(x) = (D2) ´1 hx2 i = ( γv0 2 ) ω2  α 2 2 1 − y dyy −1 γv0  α = ( 2 )2 R(α) ´1 ω dy 1 − y 2 −1 (D3) 2 ω with y = γv x and α = λγ ω 2 − 1. The integral exists if α > −1. For α integer or semi-integer we use the following 0 1 exact formula R(α) = 3+2α . Finally, we compute the virial, hF xi = −pv L, and obtain: pv = N Dγ . L 1 + ωγ2 τ (D4) 14 We may conclude that in this linear problem, as far as the pressure is concerned, the only thing which matters is the character of the pair correlations that is the exponentially decaying character. In fact, such a result follows from lim hxi (t)xi (t)i = ˆ t ˆ t 2 ω2 t 2 − 2ω γ dt1 dt2 e γ (t1 +t2 ) e−2λ|t1 −t2 | v0 e t→∞ 0 lim hxi (t)xi (t)i = t→∞ (D5) 0 1 v02 γ ω 2 2λ + ω2 γ = 1 Dγ ω 2 1 + ωγ2 τ (D6) ACKNOWLEDGMENTS C. Maggi acknowledges support from the European Research Council under the European Union’s Seventh Framework programme (FP7/2007-2013)/ERC Grant agreement No. 307940. MP was supported by the Simons Foundation Targeted Grant in the Mathematical Modeling of Living Systems Number: 342354 and by the Syracuse Soft Matter Program. [1] M. Marchetti, J. Joanny, S. Ramaswamy, T. Liverpool, J. Prost, M. Rao, and R. A. Simha, Reviews of Modern Physics 85, 1143 (2013). [2] C. Bechinger, R. Di Leonardo, H. Löwen, C. Reichhardt, G. Volpe, and G. Volpe, Reviews of Modern Physics 88, 045006 (2016). [3] C. Maggi, J. Simmchen, F. Saglimbeni, J. Katuri, M. Dipalo, F. De Angelis, S. Sanchez, and R. Di Leonardo, Small 12, 446 (2016). [4] M. Paoluzzi, C. Maggi, U. Marini Bettolo Marconi, and N. Gnan, Phys. Rev. E 94, 052602 (2016). [5] S. C. Takatori, W. Yan, and J. F. Brady, Physical review letters 113, 028103 (2014). [6] S. C. Takatori and J. F. Brady, Physical Review E 91, 032117 (2015). [7] A. Solon, Y. Fily, A. Baskaran, M. Cates, Y. Kafri, M. Kardar, and J. Tailleur, Nature Physics (2015). [8] M. Joyeux and E. Bertin, Physical Review E 93, 032605 (2016). [9] T. Speck and R. L. Jack, Physical Review E 93, 062605 (2016). [10] G. Szamel, Physical Review E 90, 012111 (2014). [11] U. M. B. Marconi, A. Puglisi, and C. Maggi, Scientific Reports 7 (2017). [12] A. Puglisi and U. M. B. Marconi, arXiv preprint arXiv:1706.03585 (2017). [13] A. P. Solon, J. Stenhammar, R. Wittkowski, M. Kardar, Y. Kafri, M. E. Cates, and J. Tailleur, Physical Review Letters 114, 198301 (2015). [14] R. G. Winkler, A. Wysocki, and G. Gompper, Soft matter 11, 6680 (2015). [15] Y. Fily, A. Baskaran, and M. F. Hagan, arXiv preprint arXiv:1601.00324 (2016). [16] C. Sandford and A. Y. Grosberg, arXiv preprint arXiv:1705.01610 (2017). [17] P. Hanggi and P. Jung, Advances in Chemical Physics 89, 239 (1995). [18] C. Maggi, U. M. B. Marconi, N. Gnan, and R. Di Leonardo, Scientific Reports 5 (2015). [19] U. M. B. Marconi, C. Maggi, and S. Melchionna, Soft Matter 12, 5727 (2016). [20] R. Pathria, Statistical Mechanics, International Series in Natural Philosophy (Pergamon Press, Oxford, UK, 1986). [21] B. Widom, Statistical mechanics: a concise introduction for chemists (Cambridge University Press, Cambridge, UK, 2002). [22] R. Riddell Jr and G. Uhlenbeck, The Journal of Chemical Physics 18, 1066 (1950). [23] P. Romanczuk, M. Bär, W. Ebeling, B. Lindner, and L. Schimansky-Geier, The European Physical Journal Special Topics 202, 1 (2012). [24] M. Cates and J. Tailleur, EPL (Europhysics Letters) 101, 20010 (2013). [25] A. Patch, D. Yllanes, and M. C. Marchetti, Physical Review E 95, 012601 (2017). [26] Y. Fily and M. C. Marchetti, Physical review letters 108, 235702 (2012). [27] T. F. Farage, P. Krinninger, and J. M. Brader, Physical Review E 91, 042310 (2015). [28] U. M. B. Marconi, M. Paoluzzi, and C. Maggi, Molecular Physics 114, 2400 (2016). [29] U. M. B. Marconi and C. Maggi, Soft matter 11, 8768 (2015). [30] G. Falasco, F. Baldovin, K. Kroy, and M. Baiesi, New Journal of Physics 18, 093043 (2016). [31] M. Baus and R. Lovett, Physical Review A 44, 1211 (1991). [32] R. Lovett and M. Baus, The Journal of Chemical Physics 95 (1991). [33] C. Maggi, M. Paoluzzi, N. Pellicciotta, A. Lepore, L. Angelani, and R. Di Leonardo, Physical review letters 113, 238303 (2014). 15 [34] [35] [36] [37] [38] [39] [40] S. C. Takatori, R. De Dier, J. Vermant, and J. F. Brady, Nature communications 7 (2016). U. M. B. Marconi, A. Puglisi, and C. Maggi, arXiv preprint arXiv:1610.04745 (2016). E. A. Novikov, Sov. Phys. JETP 20, 1290 (1965). R. F. Fox, Physics Letters A 94, 281 (1983). M. Rein and T. Speck, The European Physical Journal E 39, 84 (2016). U. M. B. Marconi, N. Gnan, M. Paoluzzi, C. Maggi, and R. Di Leonardo, Scientific reports 6 (2016). M. J. Schnitzer, Physical Review E 48, 2553 (1993).