Effective equilibrium states in the colorednoise model for active matter I. Pairwise forces in the Fox and unified colored noise approximations
Abstract
The equations of motion of active systems can be modeled in terms of OrnsteinUhlenbeck processes (OUPs) with appropriate correlators. For further theoretical studies, these should be approximated to yield a Markovian picture for the dynamics and a simplified steadystate condition. We perform a comparative study of the Unified Colored Noise Approximation (UCNA) and the approximation scheme by Fox recently employed within this context. We review the approximations necessary to define effective interaction potentials in the lowdensity limit and study the conditions for which these represent the behavior observed in twobody simulations for the OUPs model and Active Brownian particles. The demonstrated limitations of the theory for potentials with a negative slope or curvature can be qualitatively corrected by a new empirical modification. In general, we find that in the presence of translational white noise the Fox approach is more accurate. Finally, we examine an alternative way to define a forcebalance condition in the limit of small activity.
Contents
 1 Introduction
 2 Theory
 3 Effectivepotential approximation (EPA)
 4 Lowactivity approximation
 5 Conclusions
 A Simulation details
 B Effective manybody force in the Fox approximation
 C Integration of the first equality in Eq. (12)
 D Integration of the second equality in Eq. (12)
1 Introduction
Active Brownian particles (ABPs) provide a simple, minimal model system to study the collective behavior of active matter. Manybody Brownian dynamics simulations of these systems have provided considerable insight into a range of interesting nonequilibrium phenomena, such as the accumulation of particles at boundaries Elgeti (); Kantsler (); Xiao (); Yang (); Ni () and motilityinduced phase separation CatesReview (). Much of the phenomenology of ABPs can be captured using coarsegrained, hydrodynamic theories CatesReview (); Cates_PNAS (); fily (); SpeckLoewen (); gonnella2015motility (), which do not contain all information about the interparticle correlations. Some progress has recently been made in the linear response regime, which allows to decouple the equations of motions of the onebody density and polarization vector SharmaBrader ().
Due to the inherent difficulty of dealing simultaneously with both the translational and orientational degrees of freedom in active systems, attempts to develop a firstprinciples theory have largely focused on a simpler, related model, in which the particle dynamics are represented by a set of coupled OrnsteinUhlenbeck processes (OUPs). Within this model an exponentially correlated noise term, with a given correlation time, serves as proxy for the persistent trajectories of ABPs (connections between the two models were explored in Ref. faragebrader2015, ). While the removal of orientational degrees of freedom does indeed simplify the problem, it comes at the cost that one has to deal with the nonMarkovian dynamics of the translational coordinates. Fortunately, there exist various approximation methods grigolini0 (); grigolini (); fox1 (); fox2 (); ucna1 (); ucna2 (); hasegawa () which enable the OUP model to be represented using an effective Markovian, and therefore tractable, dynamics. Two different approaches to doing this, (i) the Unified Colored Noise Approximation (UCNA) of Hänggi et al. ucna1 (); ucna2 (), based on adiabatic elimination on the level of the Langevin equations, and (ii) the Fox approximation fox1 (); fox2 (), for which an approximate FokkerPlanck equation is developed, have recently been adopted in the context of developing simple theoretical tools to describe active particles faragebrader2015 (); wittmannbrader2016 (); maggi2015sr (); marconi2015 (); marconi2016mp (); marconi2016sr (); marconi2016 (); sharma2016 ().
When applied to active matter, both the UCNA and Fox approximations are referred to as ‘effective equilibrium’ approaches. The Markovian character of the dynamics implies that they obey a FokkerPlanck equation from which an effective probability distribution can in principle be obtained. Indeed, the possibility of mimicking the behavior of nonequilibrium ABPs using an equilibrium system of passive particles, interacting via effective interactions, was suggested by several researchers (see, e.g., Ref. SchwarzLinek, ) who observed that the phase separation induced by activity in systems of repulsive ABPs closely resembles that familiar from passive systems with an attractive interaction. Despite its appeal, several years were required before this observation could be turned into something more concrete. By starting from the simpler OUP model it became possible, via application of the UCNA maggi2015sr (); marconi2015 (); marconi2016mp (); marconi2016sr (); marconi2016 () and Fox faragebrader2015 (); wittmannbrader2016 (); sharma2016 () approximations, to put the notion of an effective equilibrium description on a firmer footing.
In this paper, we will compare and contrast the two different approaches to effective equilibrium. We will highlight the main approximations involved and assess the validity of the effectivepotential approximation (EPA), which has been employed in previous work to investigate activityinduced modifications of the microstructure and motilityinduced phase separation faragebrader2015 (). This analysis both clarifies the nature of the approximations involved and suggests ways in which the description can be improved.
The paper is laid out as follows: In Sec. 2 we first specify the model under consideration and describe the UCNA and Fox approaches to obtaining an effective equilibrium picture highlighting similarities and differences between them. In Sec. 3 we describe in detail the EPA, where the emphasis is placed on the UCNA due to its simpler structure. The resulting approximate effective potentials are compared to computer simulations using a standard softrepulsive and a nonconvex (Gaussian core) potential. In Sec. 4 we consider an alternative approach to obtain pairwise forces, i.e., the lowactivity limit, and make contact to the EPA. Finally, we conclude in Sec. 5.
2 Theory
In this section, we introduce the common starting point of both the UCNA and Fox approach. Since particles driven by Gaussian colored noise originally were not intended as a model for an active system, the choice of parameters in the literature may depend on the dimensionality and on whether contact to ABPs is made faragebrader2015 () or not maggi2015sr (). We will also clarify some notational issues.
2.1 Colorednoise model
We consider the coupled stochastic (Langevin) differential equations
(1) 
of particles. The motion of each particle at position is determined by conservative and stochastic forces and . The friction coefficient is related to the translational Brownian diffusivity and is the inverse temperature. We assume that the total interaction force can be written as the gradient of a pairwise additive manybody potential
(2) 
consisting of the onebody external fields and the interparticle potentials .
The vector represents the translational Brownian diffusion by a Gaussian (white) noise of zero mean and with the unit matrix . Here and in the following the dyadic product of two vectors with components results in a matrix. Any contraction as in a scalar product or a matrixvector product will be explicitly indicated by a “”. Hereafter, we shortly refer to the variable as (thermal) noise. The OUPs defined by
(3) 
with describe a fluctuating propulsion velocity as a nonGaussian (colored) noise of zero mean and
(4) 
Here we introduced the active time scale at which the orientation randomizes and the active diffusion coefficient , where is the average squared selfpropulsion velocity and the spatial dimension.
The colorednoise model for active particles contains two parameters describing the magnitude and persistence of the self propulsion. We now aim to clarify some notational differences in the literature. The persistence time can be explicitly related to the equations of motion of run and tumble particles cates2008 () or ABPs cates2013 (); faragebrader2015 (). The above definitions correspond to the latter case with , where is the rotational Brownian diffusion coefficient. Another common choice marconi2016mp () amounts to consider and . In the following, we use , where is the typical diameter of a particle and the dimensionless persistence time has been introduced in Refs. faragebrader2015 (); wittmannbrader2016 ().
Since the dimensionless active diffusivity implicitly depends on the persistence time, it constitutes the most general measure for the activity (together with the persistence length ). In order to connect to a system of ABPs faragebrader2015 (); wittmannbrader2016 (), it is convenient to consider instead of a dimensionless velocity , i.e., the Peclét number. In the literature, some other definitions of a Peclét number are used, which we will not consider here. One peculiar property of a system of active OUPs is that, even at vanishing selfpropulsion velocity , or , there is a contribution of to Eq. (1) arising from a finite reorientation time sharma2016 (). One thus does not recover the equation of motion of a passive (Brownian) particle, as in the ABPs model. In the longtime limit, however, the contribution to the dynamics becomes irrelevant and the same steady state is described as for a passive Brownian particle, see appendix A for more details. A proper passive system can only be recovered from Eq. (1) only in the limit , in the sense that the velocity correlation in Eq. (4) reduces to a white noise. A Brownian system is then represented by when the thermalnoise variable is removed, or, trivially, by setting which amounts to neglect the contribution of the OUPs.
2.2 Effective equilibrium approach
The most important step towards a theoretical study of the OUPs model is to derive from the nonMarkovian stochastic process (1) an equation of motion for the particle probability distribution . In this section, we will discuss the differences between the multidimensional generalizations of the UCNA ucna1 (); ucna2 () and the Fox fox1 (); fox2 () approaches to effective equilibrium and expound the surprising similarities between these two approximations in the (currentfree) steady state.
As a central quantity emerging in both cases, we define the friction tensor with the components
(5) 
resulting in the Hessian of and the diagonal block
(6) 
not to be confused with for particle. In the following, we briefly denote by the th block component of the inverse tensor .
The UCNA maggi2015sr (); marconi2015 () amounts to explicitly inserting the OUPs (3) into the overdamped limit of the time derivative of (1), resulting in the modified Langevin equation . It is now straightforward to obtain for this (approximate) Markovian system driven by white noise the Smoluchowski equation with the probability current (the superscript (u) denotes that the UCNA has been used)
(7) 
Note that the UCNA remains valid as long as the friction tensor (5) is positive definite.
The Fox approximation scheme applied to (1), on the other hand, only makes use of the correlator (4) of the OUPs, which, in turn, may also be interpreted as the correlator of corresponding to a coarsegrained equation of motion representing ABPs with a constant velocity in the direction of their instantaneous orientation that is subject to Brownian rotational diffusion faragebrader2015 (). This method directly yields the approximate Smoluchowski equation (superscript (f)) with sharma2016 (); SpeckCRIT ()
(8) 
where the regime of validity is the same as for UCNA. The major difference between Eq. (7) and (8) only affects the effective description of the dynamics as a result of the additional factor arising on the level of the Langevin equation within the UCNA. Note that in the original generalization of the Fox result faragebrader2015 () the tensor from Eq. (5) was incorrectly obtained as , which we will later identify as the (diagonal) Laplacian approximation. It will turn out that this (or another) approximation is necessary to obtain physical expressions for the effective interaction potentials.
2.3 Two versions of the steadystate condition
In contrast to the dynamical problem, the (currentfree) steadystate conditions
(9) 
for the stationary distribution can be cast in a coherent form, defining the effective diffusion tensor , such that only the components
(10)  
(11) 
differ between the UCNA (u) and Fox (f) results.
Multiplying Eq. (9) with and summing over repeated indices, the steadystate condition takes the more instructive (approximate) form marconi2015 ()
(12) 
introducing the effective force . The term is an approximation for , which becomes exact in the UCNA marconi2015 (). For the Fox approach, we argue in appendix B that this is still true in some important special cases, such that Eq. (12) is accurate enough for our purpose. For high particle numbers the contribution of the offdiagonal elements to becomes increasingly irrelevant marconi2015 (), which amounts to setting . Assuming this diagonal form, the determinant in Eq. (12) can be replaced according to as we have before approximating the expression in the last step (compare appendix B).
Putting aside the dynamical behavior, described in Sec. 2.2, the only difference between the UCNA or Fox approximation is manifest in the definitions, (10) and (11), of . Using UCNA the active diffusivity only appears as part of a prefactor in (10), so that the friction matrix , representing a correction due to activity, contributes to the steadystate result even in the case , that is when and . Hence, the logical parameter suggested by the UCNA to tune the activity is , with the passive system () restored only in the limit . This appears to be an artifact of the pathological contribution of to the displacement of the OUPs, whereas the connection to the experimentally more relevant system of ABPs is lost. For the latter it appears more natural to tune at constant . In the derivation of the Fox result (11), on the other hand, the explicit time evolution of the OUPs in Eq. (3) is irrelevant, suggesting a better approximate representation of ABPs faragebrader2015 (). This reflects that we recover the (same) passive system for either or (in the presence of noise).
Ignoring the noise contribution for , the UCNA and Fox approximations practically describe the equivalent effective steady states. The major advantage of this approximation, or the UCNA result in general, is that the inverse is pairwise additive, even if . Then the effective manybody potential defined as can be written in a closed form maggi2015sr (); marconi2015 (), admitting the explicit solution of Eq. (9). Due to the more nested form of Eq. (11) the Fox approximation does in general not admit an analytic result. As is not pairwise additive in either approach, some further approximations will become necessary to construct a predictive theory, which we discuss in the following sections.
3 Effectivepotential approximation (EPA)
Regarding the possible applications using standard methods of equilibrium liquidstate theory a desirable strategy is to approximate in Eq. (12) in terms of pair potentials. This approach allows to describe the phase behavior of ABPs approximated as particles propelled by a set of coupled OUPs, which has been discussed in detail faragebrader2015 (); wittmannbrader2016 () for passive softrepulsive and LennardJones interactions in three dimensions. However, it can be criticized that (I.i) a system which obeys detailed balance is used to represent the interactions in an active system, (I.ii) the validity criteria of the underlying theory might be violated so that further approximations are required and (I.iii) higherorder particle interactions are neglected, which are believed to be important for the phase separation in an active system. In the following, we define the effective pair interaction and motivate different approximations, which we compare to computer simulations of two ABPs and two particles propelled by OUPs. It is our objective to comment on the aforementioned points and illustrate the qualitative differences between the Fox and UCNA. For the sake of simplicity, we will restrict the presentation of technical aspects to the UCNA results.
3.1 Calculation of effective potentials
To identify an effective pair potential , we consider interacting particles, i.e., the lowdensity limit of Eq. (12). Ignoring the external forces for now by setting , it is easy to verify that
(13) 
and analog equation for , where we used and . Keeping in mind that we seek to employ this effective potential to approximately represent the interaction of many particles, it appears undesirable that an equal statistical weight is put to both the diagonal and the offdiagonal components of the diffusion tensor. As an alternative we propose the effective potential
(14) 
with , obtained for a diagonal form of with . For completeness we find in the oneparticle limit a quite similar formula
(15) 
for the effective external field , since for we have . Note that a quite different expression for an effective external potential can be derived starting from the equations of motion for ABPs wittmannbrader2016 (); pototsky2012 ().
Integration of the above equalities yields the desired formulas for the effective potentials depending only on the bare potential or and the activity parameters and faragebrader2015 (). Alternatively, we could have directly defined marconi2015 (); marconi2016mp () and from the manybody potential identified in the solution of (9), which is, however, inconvenient when the Fox approach is used. Assuming a bare potential obeying , the integrated form of (13) reads
(16) 
where and
(17) 
are the Eigenvalues of . We can further identify in the first term of Eq. (13) as the Eigenvector of corresponding to the Eigenvalue . Note that in (16) we could equally introduce an effective energy scale to absorb the factor marconi2015 (); marconi2016mp (). We refrain to do so as this interpretation would not be consistent with the the way enters within the Fox approach.
3.2 Limitations and possible corrections
Studying Eq. (16) more carefully, we notice that the effective potentials do not always behave in a physical way. This is because, in violation of the validity condition of both the UCNA and the Fox approximation, the diffusion tensor is not positive definite for a large number of relevant potentials. In general, we easily see that the logarithm will diverge whenever one of the Eigenvalues vanishes. Given a positive and convex bare potential , the Eigenvalue is strictly positive, which also means that the effective attraction solely arises from the term including the logaritm. However, as we have in this case, the Eigenvalue will vanish at a certain value of and we require a further approximation to remedy the unphysical behavior of in dimensions. At a highly nonconvex or negative region of the bare potential, the same problem occurs for . Interestingly, if we only require knowledge of an effective potential on a finite interval where the eigenvalues are positive, its overall unphysical behavior is irrelevant sharma2016 ().
First note that there is a broader range of admissible bare potentials when the diagonal approximation, Eq. (14), of the effective pair potential is used, or if we are interested in the onebody external field, Eq. (15). This can be understood from the explicit formula for , which we obtain from Eq. (16) by rescaling all terms proportional to with a factor . In the following, we propose different ways to generally rid the effective potential of possible artifacts of vanishing Eigenvalues in the last term of Eq. (13). A correction of the first term is not necessary and also has no noticeable effect.
Let us first assume that is convex, i.e., it represents a softrepulsive interaction. Then a sufficient criterion (due to the presence of the term in Eq. (5), some other potentials are allowed that are only slightly negative and slightly nonconvex) for the matrix to have strictly positive Eigenvalues would be that it depends on an elliptic differential operator rather than . Therefore, a convenient approximation is to redefine Eq. (5) by an elliptical operator, the simplest example of which is the Laplacian . Upon substituting
(18) 
the effective potential becomes
(19) 
where the additional term compared to (16) cannot be integrated in general. This Laplacian approximation has been successfully employed (together with the Fox and diagonal approximation) in explicit calculations faragebrader2015 (); wittmannbrader2016 (). In dimensions both differential operators reduce to the second derivative and Eq. (19) is equal to Eq. (16), which provides a good account of active particles interacting with a softrepulsive potential marconi2016mp ().
An alternative way is to empirically rectify the explicit formula for in Eq. (16). Most intuitively, one can expand the argument of the logarithm up to the first order in . In fact, this small approximation is quite similar to the Laplacian approximation (19) (and completely equivalent in one dimension), but we do not recover the additional term involving the integral. Performing the small approximation of the full expression (16) appears too crude, as an expansion of the logarithm does not converge for . The resulting effective potential will thus become totally uncontrolled for short separations of highlyrepulsive particles.
A more elaborate correction that may be applied also to highly nonconvex potentials is the inverse approximation, an empirical strategy maintaining the leading order in while not disregarding higherorder terms. This is achieved by substituting in Eq. (16) , where
(20) 
The major advantage of this approximation is that it yields quite similar results to the full potential whenever the validity condition is only slightly violated and the effective potential does not diverge if the bare potential is finite. The empirical motivation behind this correction is that constitutes the two leading terms of the “resummed” Taylor series of in the case . As described in appendix B, the most convenient implementation of the inverse approximation for the Fox result is to identify the expressions for in and use Eq. (20).
3.3 Comparison to computer simulations of two active particles
In Sec. 3.1 we introduced different strategies to define a suitable effective interaction potential in the effectiveequilibrium approximation for the colorednoise model. Now we illustrate under which conditions an approximate treatment according to Sec. 3.2 becomes necessary and compare the theoretical results to computer simulations. The easiest way to determine an effective potential numerically is to set up a twoparticle simulation, measure the radial distribution function and calculate . By doing so, we make the same approximation (I.iii) as in theory to ignore the manyparticle character of the interaction. However, the simulations for ABPs and OUPs, detailed in appendix A, take into account the orientation dependence and the nonMarkovian character of the dynamics, respectively.
3.3.1 The role of approximations, dimensionality and thermal noise
We first discuss some general observations in the UCNA for a softrepulsive system with the bare potential . The behavior of the Fox result is qualitatively similar. As expected, the full expression for the effective potential in Eq. (16) is impractical as it diverges at a certain distance , determined by the condition , which is when the first Eigenvalue within the logarithm vanishes, whereas is always positive. As suggested by Fig. 1a, this behavior is most problematic at larger values of , where should be rather negative, as it is the case in dimensions. We further see in Fig. 1a that this effect becomes more severe with increasing dimension. Both the inverse and Laplacian approximations successfully cure this unphysical divergence, which we see in Fig. 1b. As employing the diagonal form of the effective potential simply amounts to a rescaling of , we observe in Fig. 1b that it results in a smaller effective diameter of the repulsive part but a flatter potential well. Accordingly, becomes smaller.
Since the definition of the active diffusivity depends on the dimension , not all parameters , and can be kept constant upon varying the dimensionality. For a constant reorientation time and propulsion velocity the effective attraction in Fig. 1c is stronger in lower spatial dimensions for both approximations considered, which coincide with the full expression in . This behavior appears sensible, as two particles have less possibilities to avoid each other upon collision, and is in qualitative agreement with computer simulations of active OUPs. Moreover, we understand that motilityinduced phase separation is harder to observe in higher dimensions stenhammar2014 (). Keeping constant instead of (which then decreases with decreasing dimension) the same trend is observed in Fig. 1d for the inverse approximation and computer simulations, whereas the result for the Laplacian approximation barely changes with dimensionality. Comparing the qualitative behavior in Figs. 1c and d, we recognize in all spatial dimensions that the numerical effective potential is of longer range than the theoretical predictions in any approximation. This observation confirms the criterion discussed in Ref. ucna2, that the UCNA is expected to become less accurate for larger separations where a typical length scale of the active motion, closely related to the effective diffusion tensor, Eq. (10), exceeds the spatial scale over which the force field varies.
As the involved approximations become cruder in higher spatial dimensions, the quantitative agreement with the simulation results in Figs. 1c and d becomes worse. For , a remarkable agreement between the UCNA and simulation results for the radial distribution of two particles has been reported in Ref. maggi2015sr, , where in Eq. (1) was set to zero. Doing so also in our simulations, we observe in Fig. 2a that the effective potential deepens and its repulsive barrier becomes steeper. This curve is in excellent agreement with the theoretical result for zero noise, obtained by both the UCNA and Fox approach upon dropping the first term in Eq. (10) and Eq. (11), respectively. The full UCNA result is only slightly different in the repulsive regime. Intriguingly, we also recognize in Fig. 2a that the simulation data including the noise term are excellently represented by the Fox approach. We can understand these observations by recapitulating the idea behind the two approximations. The UCNA amounts to manipulate Eq. (1) by calculating the second derivative of in order to eliminate the variable in Eq. (3). The original discussion of the accuracy of this approximation does not account for the second stochastic variable . In contrast, the Fox approach is only dedicated to determine the approximate contribution of the colored noise to the effective probability current in Eq. (8), which is independent of other terms in Eq. (1). Therefore, the Fox theory has a broader range of applicability and should be accurate in both the presence and the absence of thermal noise.
Finally, we note that the excellent agreement between theory and computer simulations in one dimension implies that the diagonal approximation is not justified when it comes to describing a twobody system. It is, however, interesting to consider a single particle in an external field of the same form as the interparticle potential considered above. In agreement with the theoretical prediction, the computer simulations in Fig. 2b show nice agreement between the twobody system and a onebody system with the double value of the persistence time. We further observe that the theoretical result for one body is even closer to the simulation data than for two bodies.
3.3.2 Softrepulsive Brownian system in three dimensions
We also performed computer simulations of ABPs (described in appendix A) for . For a finite active diffusivity, Fig. 2c reveals that the numerical effective potentials for the two considered models with and without noise are nearly identical over the full range of separations. As for , the effective potential of active OUPs (and ABPs) in the absence of thermal noise has a deeper well and a larger repulsive diameter. Quantitatively, this difference is much more pronounced in three dimensions. In the following, we restrict ourselves to systems with thermal noise and compare in Figs. 2d and 3 to the predictions of the theory. As in Figs. 2c and d, our simulations of ABPs and OUPs, shown in the first column of Fig. 3, are in nice agreement for all sets of considered parameters. This is quite surprising since, on the manyparticle level, ABPs and OUPs have different steady states solonEPJST (); szamel2014 (). On the basis of our data for the simplistic twobody system we could rather conclude that OUPs subjected to thermal noise are an excellent model for ABPs at moderate activity faragebrader2015 (); wittmannbrader2016 ().
We see in Fig. 2d that the depth of the attractive well of all theoretical versions of the effective potential in dimensions is significantly overestimated when compared to the simulations of both ABPs and OUPs. The inverse approximation appears to provide the best guess of the point at which the effective potential changes its sign. To facilitate the further qualitative comparison we chose the axes in Fig. 3 according to the deviation from the simulations (first column), i.e., by a factor of 10 for the Laplacian approximation (second column) and 5 for the inverse approximation (third column). The chosen approximations exhibit a behavior similar to the full theoretical results of Eq. (16) in the physical region for , shown in the last column. The differences arising from using the exact effective force in Eq. (12) are discussed in appendix B.
For the considered softrepulsive bare potential, we observe in Fig. 3 some notable quantitative differences between the UCNA and Fox results, even at relatively high . The effective diameter of the repulsive part is generally smaller than in the UCNA, whereas the overall attraction is weaker in the Fox approach. This becomes most apparent in the Laplacian approximation. The first row of Fig. 3 contains the effective potentials evolving for a constant persistence time when the active diffusivity (or the Peclét number ) is increased. All approaches accordingly predict an increased effective attraction and the minimum of the potential is shifted to smaller separations faragebrader2015 (); wittmannbrader2016 (). The Fox results exhibit a stronger variation with , which is also more consistent with the numerical data. Similarly, the effective potentials in the second row deepen with increasing at constant , where the location of the minimum is nearly unaffected. The most interesting behavior is observed in the third column at constant . Again, all approaches agree that the minimum is shifted to larger separations with increasing , but the effective attraction predicted by the simulations is nearly constant, as simultaneously the magnitude of the selfpropulsion is decreased. This observation is not consistent with the UCNA results.
Based on the presented simple comparison, our conclusion is that the best choice for the theoretical effective potential is the Fox approach in the inverse approximation. Upon further increasing the activity (not shown) the quantitative discrepancy of the Laplacian approximation becomes even more pronounced. Regarding Fig. 1b one might get the impression that the additional assumption of the diagonal form of the effective potential also results in a slightly better (quantitative) agreement with the simulations. However, we stress that it is not clear in how far the effective pair potentials can accurately describe the manybody situation, as there are no higherorder interactions present in a twobody simulation.
3.3.3 Difficulties for nonconvex potentials
The most compelling argument in favor of the inverse approximation arises from considering nonconvex bare potentials, a case in which the Laplacian approximation becomes useless above a certain value of . This issue has already been discussed for a LennardJones potential SpeckCRIT () and we will consider here the potential of an active Gaussiancore fluid. Although this model has not received much attention in theories for ABPs, it is quite appealing from a theoretical perspective. Most prominently, this bare potential is a known exceptional case in which a simple meanfield theory is particularly accurate GaussianCore1 (); GaussianCore2 (), which might, in a way, also hold for the effective potential of the active system.
The effective potential of an active Gaussiancore fluid is discussed in Fig. 4 within the UCNA. As the absolute value of both the curvature and slope of this model potential is bounded, the Fox results (not shown) are quite similar, even for the moderate values of considered here. When the persistence time is sufficiently small, the effective potential, Eq. (16), does not diverge and the different approximations behave in a quite similar way. Interestingly, we observe in Fig. 4a that the divergence of the Laplacian approximation sets in at an even smaller value of than for the full potential. The latter diverges at two points, each related to one of the two Eigenvalues, if .
In Fig. 4b, c and d we discuss the only suitable form of the effective potential, i.e., the inverse approximation. Intriguingly, the predicted behavior depends on in which way, i.e., by means of which parameter, the activity is modified. Increasing the active diffusivity (or the Peclét number) at a constant value of results in a less repulsive core. Upon increasing , however, the height of the maximum increases and an attractive well develops at larger separations. Quite counterintuitively, we observe that in this case the effective interaction becomes more repulsive than in the passive case, also for the inverse approximation. Our computer simulations (also carried out for values of much larger than shown in Fig. 4) confirm that this is an artifact of the theory, related to the negative curvature of the bare potential. At constant the evolution of the theoretical results agrees qualitatively with the simulations. The simulation data are, however, not very sensitive to changes in the persistence time. At constant the theory predicts the correct trend upon increasing , whereas this is not the case at constant Peclét number.
To argue about the validity of the EPA, we consider two classes of bare interactions. Firstly, softrepulsive and convex potentials lead to quite accurate results in one dimension, but require an empirical correction in higher dimensions. Secondly, the understanding of the behavior of particles interacting with a bare potential which has a negative curvature remains one of the most urgent open problems in our theoretical framework. At the moment, the only way to obtain a workable theory in this case is by employing the inverse correction introduced in Eq. (20). Further numerical and theoretical analysis will be needed to fully clarify this issue. Finally, we note that potentials with attractive parts do not a priori constitute a problem for the theory, but usually have regions in which they are non convex. For a discussion of such problems see also Refs. SpeckCRIT, ; sharma2016, .
4 Lowactivity approximation
A second strategy to simplify the steadystate condition is to perform an expansion in the activity parameter . At linear order, the effective diffusion tensor becomes pairwise additive. In this lowactivity approximation, a YBGlike hierarchy can be obtained by successively integrating Eq. (9) over coordinates marconi2015 (), which allows defining a mechanical pressure and interfacial tension marconi2016 (). Moreover, for the active system evolving according to Eqs. (1) and (4) it has been demonstrated that there exists a regime for small values of , where the principle of detailed balance is respected fodor2016 (). This suggests that, at leading order in this parameter, the approximations resulting in Eq. (9) are perfectly justified.
Knowing, however, that Eq. (9) contains the same information as Eq. (12), which depends logarithmically on the parameter , we should clarify whether (II.i) the lowactivity expansion converges, (II.ii) it is sufficient to only consider the leading order and (II.iii) one can obtain similar results when employing the EPA. To do so, we demonstrate how the first member () of the YBG hierarchy can be rederived from Eq. (12) and discuss the consequences of approximating in terms of pair interactions. Again, we only discuss the UCNA, where, without any further approximation, the inverse diffusion tensor is found to be pairwise additive.
4.1 Alternative derivation of a local forcebalance condition
By integrating a multidimensional vector equation (label ) over coordinates we understand multiplying each side by , followed by summation over all particles and integation over all spatial coordinates. We further define the average of a vector , where denotes the full canonical ensemble average and is the average of the density operator . Approximating now the inverse mobility matrix in Eq. (10) as and integrating Eq. (9) over coordinates, we find the first member
(21) 
of a YBGlike hierarchy marconi2015 (); marconi2016 () for the active system, where, explicitly
(22) 
for any (nontrivial) differential operator acting on . In the derivation of (21) it turns out that the offdiagonal components of the mobility tensor do not contribute at first order in marconi2016 (). Hence, we might as well have assumed the diagonal form at linear order in beforehand.
In order to connect to the EPA, we derive a YBGlike hierarchy from Eq. (12). Assuming the diagonal form , the integration over coordinates of the first equality is carried out i n appendix C. Making use of the equilibrium version of the YBG hierarchy and expanding the expression up to first order in the result is
(23) 
introducing the averaged inverse diffusion tensor (compare Eq. (10))
(24) 
Multiplying Eq. (23) with it is easy to verify in appendix C that at first order in it becomes equivalent to (21) up to a term proportional to the expression in the second line, which we consider as a higherorder contribution.
In order to derive Eq. (23) in the Fox approach, an additional approximation is required, as the inverse of from Eq. (11) is not proportional to . We would thus need to redefine in Eq. (24) according to Eq. (11) where takes the role of . Regarding in general the presented alternative derivation of Eq. (21), its validity appears to be in question. This is because to derive the intermediate result in Eq. (23) it is necessary to expand a logarithmic term, the Taylor series of which only has a finite radius of convergence. We further explicitly assumed the diagonal form of the mobility tensor to avoid further terms that are not present in the original result. Employing in the next step the EPA will shed more light on these issues.
4.2 Lowactivity approximation of effective potentials
Having established a connection between (9) and (12) also at linear order in , we now turn to the case in which the second equality in (12) does not hold. This is when we assume along the lines of (2) but within the EPA using the results derived in Sec. 3.1. As detailed in appendix D, the obvious result is that all correlation functions between more than two particles vanish in the approximate integrated version
(25) 
of (12). Ignoring the interparticle interactions the approximation involving only becomes exact. This situation is the same as discussed in Ref. marconi2015, .
Considering the interacting system, we multiply Eq. (25) with as done previously for (23). According to appendix D, we can only approximately reproduce Eq. (21) by doing so. This reflects both the limitations of the EPA and an inconsistency between Eq. (9) and Eq. (12) when being subject to the same type of approximation, as we discuss in the following. We observe that (III.i) the coupling between external and internal interactions is ignored by Eq. (25) (III.ii) spurious threebody correlations appear on the lefthandside of Eq. (21) (III.iii) the second term on the righthand side of Eq. (21) is recovered but involves a seemingly unjustified expansion and (III.iv) if we do not explicitly assume a diagonal diffusion tensor, the last term in Eq. (21) changes by a factor two. As we are mainly interested in bulk systems, the first point is only briefly commented on in appendix D.
The term including the bare interaction force in Eq. (12) depends on the position of three bodies. Hence, the pairwise approximation, which amounts to setting
(26) 
should not be too crude. Moreover, we have discussed in Sec. 3.2 that such a contribution to the effective force is usually purely repulsive and thus plays only a minor role in characterizing a possible phase transition. However, we show in appendix D that the definition (24) of the averaged diffusion tensor is not fully compatible with the EPA, resulting in point (III.ii). This is in contrast to the clean derivation of Eq. (21), where Eq. (9) is recovered from Eq. (12) by multiplication with the manybody effective diffusion tensor before integrating over positions.
The last term in (12), although considered here for a diagonal diffusion tensor, constitutes a full body quantity. Recall from the discussion in Sec. 3 that an approximation as a pairwise quantity might be quite poor and an expansion of the logarithm does not converge. However, we demonstrate in appendix D that successively employing the EPA and expanding for small according to
(27) 
eventually results in full consistency with the respective term in Eq. (21), stated as point (III.iii). This suggests that the expansion to first order in ”implies” making the EPA when Eq. (12) is our starting point. Despite the aforementioned crudity of this expansion, we argue that Eq. (21) is valid, as its clean derivation from Eq. (9) does not require dealing with a logarithmic term. The last step in Eq. (27) is required to recover Eq. (21) without inducing undesired higherorder terms in , as, similar to point (III.ii), the integrated version is incompatible with the chosen . However, we stress that this approximation should certainly be avoided when calculating the fluid structure.
Finally, we demonstrate in appendix D that the offdiagonal elements of the diffusion tensor entering in Eq. (12) contribute to Eq. (25). Hence, the present approach would be even more inconsistent with Eq. (21) if we did not assume the diagonal form, as noted in point (III.iv). We also note that the same problem occurs for the according generalization of Eq. (23). In principle we could define in this case an additional averaged diffusion tensor , corresponding to the offdiagonal elements, which could counteract this inconsistency. Such a calculation would, however, not be useful when an effective pair potential is employed.
5 Conclusions
In this paper we studied different ways to define an effective pair interaction potential between active particles. Our numerical investigation reveals that a twoparticle system of ABPs and active OUPs exhibits a quite similar behavior. These results serve as a benchmark to test the approximations involved in recent effective equilibrium approaches, which have been reviewed and compared in detail. For spatial dimensions higher than one we introduced an empirical way to rid the theoretical result of possible divergences, which also appears to yield the best agreement with the simulation data, although the effective attraction is still significantly overestimated. Regarding the quite accurate onedimensional results and the qualitative features of the effective potentials in three dimensions, the Fox approximation is superior to the UCNA when the translational Brownian noise cannot be neglected. In the absence of noise both approximation schemes admit the same steadystate solution.
Further analysis is needed to better understand the role of the neglected manybody interactions in both the twobody simulations and the theory, which are thought to be imperative for a quantitative description of active systems SpeckCRIT (). The presented theoretical approach follows two major approximate steps to define the effective pair potential. First, we map the equation of motion (1) onto a deterministic FokkerPlanck equation (effective equilibrium picture) and then we define pair forces from the twoparticle limit assuming a vanishing probability current. It could well be that the mapping in the first step breaks down parts of the manybody nature of the interactions in the active system, such that the effective attraction in the manybody system becomes accessible already on the level of pair interactions. As a logical next step, it seems worthwhile to study the effective potential extracted from computer simulations of a manyparticle system, in order to clarify in how far the strong attraction of the effective potential needs to be seen as the result of a fortuitous cancellation of errors.
The lowactivity limit in the effective equilibrium picture also results in pairwise forces. Under this assumption, we revealed some minor inconsistencies between the two equivalent steadystate conditions in Eqs. (9) and (12), although the latter contains a logarithmic term. Relatedly, it was recognized in Ref. marconi2016, that different routes to define the active pressure only coincide at lowest order in the activity parameter . We suspect that further differences will occur at higher orders in and when employing further approximations, such as the EPA. We conclude that the route to follow should be carefully chosen for each problem, together with the underlying approximations.
The obvious purpose of both the lowactivity approximation and the EPA is to allow for an analytically tractable theory. It appears that the condition given by Eq. (12) supported by effective pair potentials is most convenient for accessing structural properties marconi2016mp (); faragebrader2015 (); wittmannbrader2016 (), whereas the lowactivity expansion of Eq. (9) provides a direct way to define mechanical properties marconi2016 (). Moreover, our analysis suggests that the thermodynamic results obtained from Eq. (12) can be rescaled in order to obtain a workable definition of mechanical active pressure and surface tension. Arguably, the most simplistic scaling factor would be the diffusivity of an ideal gas, which can be absorbed into an effective temperature marconi2015 (); marconi2016 (); wittmannbrader2016 (). A more general approach will be detailed in the second paper of this series.
Acknowledgements
R. Wittmann, A. Scacchi and J. M. Brader acknowledge funding provided by the Swiss National Science Foundation.
C. Maggi acknowledges support from the European Research Council under the European Unionâs Seventh Framework programme
(FP7/20072013)/ERC Grant agreement
No. 307940.
Appendix A Simulation details
We performed Brownian Dynamics simulations of a system composed of two particles of unit diameter interacting through a softrepulsive potential or a Gaussian softcore potential. The potential is truncated at a distance of . In the simulations of active OUPs, evolving according to Eq. (1), each particle is subjected to Gaussian thermal noise and a nonGaussian (colored) noise. The latter yields two distinct contributions to the displacement of each particle: one drift term, proportional to the reorientation time and one Gaussian process, proportional to sharma2016 (). For a vanishing active diffusivity , the drift term decays exponentially in time and is therefore irrelevant in the longtime limit. The integration time step is fixed to where is the time scale of translational diffusion. The total run time of the simulation is . For every , we calculate the distance between the two particles. The paircorrelation function is obtained in a standard way from the distance distribution. We have verified that, for the case of and finite , the obtained paircorrelation function is independent of , although the shorttime displacement is not.
We also performed Brownian Dynamics simulation of ABPs, for which the colorednoise variable in Eq. (1) is replaced with the vector describing a constant velocity of the selfpropulsion in the direction of the instantaneous orientation. The equation for the time evolution for the orientation vector of each particle is evaluated as an Ito integral, where is a white noise describing rotational diffusion. The integration time step is fixed to and the total run time is .
Appendix B Effective manybody force in the Fox approximation
In this appendix we discuss the accuracy of Eq. (12) of the main text in the Fox approximation, i.e., choosing the effective diffusion tensor from Eq. (11). The accurate definition of the effective force is
(28) 
since the conversion
(29) 
is only correct in the following cases:

for a system with no thermal noise. As stated in the main text, in this case the Fox and UCNA results are equivalent. Making use of the symmetry relation (with Greek indices labeling components and particles) and Jacobi’s formula the identity in Eq. (29) can be explicitly verified marconi2015 ().

for a passive system, since .

at leading order in the activity parameter , where and the same arguments as under point 1. can be used.

for particles in an effectively onedimensional symmetry, i.e., if there exists a coordinate frame in which the nontrivial contributions to reduce to an at most a tensor with identical diagonal elements. This can be easily shown by an explicit calculation

in dimensions

for a planar interaction potential

in the Laplacian approximation (18)

As a simple counterexample to the cases listed above, we note that Eq. (29) does not generally hold for and , since , whereas under point 4. we have . As for the approximate formulas (13) and (14), we find that the difference between the effective forces, Eq. (28), for particles with and without diagonal approximation is only a factor 2 in front of each factor .
Comparing the requirements for points 1. and 2. we can say that Eq. (12) is correct for both and , indicating that it should be a good approximation over all ranges of the parameter . Moreover, the assumption of a small persistence time is required in the derivation of the effective equilibrium approach faragebrader2015 (); maggi2015sr (); marconi2015 (). Considering point 3., this means that the approximation in Eq. (29) is consistent with the underlying theory. Indeed, Figs. 5a and 5b show that the approximation best for either small or large and small , respectively. In general, the difference is not significant compared to other approximations shown in Fig. 3 of the main text.
If the validity criterion for the Eigenvalues, given by Eq. (17), of is violated, neither side of Eq. (29) results in physical effective potentials. Therefore, the most important benefit of the approximate form on the righthand side is that it allows to employ the inverse approximation, as described in Sec. 3.2: the Eigenvalues of can be written as
(30) 
so that we can substitute according to Eq. (20) of the main text. The substitution of or a more general manipulation of is inconvenient since the effective potential would still diverge for although the bare potential does not. Therefore, the third column of Fig. 3 contains the optimal implementation of the the inverse approximation for the Fox approach. Also recall that, according to point 4.(c), the results in the Laplacian approximation are the same for both expressions in Eq. (29).
Appendix C Integration of the first equality in Eq. (12)
The derivation of Eq. (21) by integrating Eq. (9) over coordinates is quite similar to that of the YBG hierarchy in a passive system. The first member
(31) 
is recovered from (21) when setting . The second equality reflects the interpretation of the term on the lefthand side as the gradient of a chemical potential , which is constant in equilibrium. The second member reads
(32) 
and is related via the second equality to the first member. With the help of these exact equilibrium sum rules we will now derive Eq. (23) by integrating Eq. (12) over coordinates. Our presentation follows closely the derivation of a dynamical density functional theory including a tensorial diffusivity rexloewenDDFT (), whereas we only consider the steadystate condition.