Active-matter isomorphs in the size-polydisperse Ornstein-Uhlenbeck Lennard-Jones model

This paper studies size-polydisperse Lennard-Jones systems described by active Ornstein-Uhlenbeck particle dynamics. The focus is on the existence of isomorphs (curves of invariant structure and dynamics) in the model's three-dimensional phase diagram. Isomorphs are traced out from a single steady-state configuration by means of the configurational-temperature method. Good invariance of the reduced-unit radial distribution function and the mean-square displacement as a function of time is demonstrated for three uniform-distribution polydispersities, 12%, 23%, and 29%. Comparing to active-matter isomorphs generated by the analytical direct-isomorph-check method, the latter give somewhat poorer invariance of the structure, but better invariance of the dynamics. We conclude that both methods can be used to quickly get an overview of the phase diagram of polydisperse AOUP models involving a potential-energy function obeying the hidden-scale-invariance property required for isomorph theory to apply.


I. INTRODUCTION
Active matter involves particles that absorb energy from their environment and continuously perform motion dissipated into heat.This kind of motion, which in contrast to standard Newtonian or Brownian dynamics breaks time-reversal invariance [1,2], is relevant not only for describing biological systems ranging from bacteria to flocking birds [3][4][5][6][7][8][9][10], but also for microscopic artificial microswimmers and active Janus particles.
Many different approaches to the description of active matter exist, depending on whether point particles or particles with directional coordinates are considered and depending on the precise mechanism by which the particles autonomously perform mechanical work [5,6,8,9,11,12].Point-particle active-matter models include the Active Brownian Particle (ABP) [13,14] and Active Ornstein-Uhlenbeck Particle (AOUP) models; these models have been used to describe the motion, e.g., in active colloids [15].The AOUP model, which is simpler than the ABP model and has one less parameter, can be used to approximate ABP dynamics.Moreover, the AOUP model offers more possibilities to obtain theoretical predictions [16,17]; this is the model we choose to study in the present paper.Specifically, the AOUP model involves point particles subject to a colored-noise Langevin dynamics [16,[18][19][20].
In view of the variability of biological and other active systems, one cannot expect all particles to be identical.As a consequence, polydispersity has recently come into focus in connection with active-matter models [21][22][23][24].There is also currently great deal of interest in passive polydisperse systems coming from, in particular, their use in SWAPequilibrated supercooled liquids [25], in which context the question arises of how similar the dynamics of small and large particles are [26][27][28].Finally, it is worth mentioning that active matter at high density has recently been studied inspired by biological materials such as cells, both for monodisperse [24,29] and polydisperse cases [30], showing emerging collective phenomena with the spontaneous occurrence of spatial velocity correlations.
This paper studies the size-polydisperse AOUP Lennard-Jones (LJ) model.We recently demonstrated the existence of lines of approximately invariant structure and dynamics in the phase diagram of a binary LJ AOUP model; such lines are referred to as "active-matter isomorphs" [31][32][33].Inspired by the fact that the introduction of polydispersity into ordinary (passive) Newtonian models does not affect the existence of isomorphs [34], the present paper investigates whether the existence of isomorphs also survives the introduction of polydispersity into the AOUP model.This is worthwhile to investigate since the existence of isomorphs makes it possible to quickly establish an overview of the phase diagram because only a single point on each isomorph needs to be simulated.

II. THE AOUP EQUATION OF MOTION AND SIMULATION DETAILS
We consider a system of N particles in volume V and define the number density by ρ ≡ N/V .If the potential-energy function is denoted by U (R) in which R ≡ (r 1 , ..., r N ) is the vector of all particle coordinates, the AOUP equation of motion [16,[18][19][20] is Here µ is the mobility (velocity over force); the noise vector η(t) is colored according to an Ornstein-Uhlenbeck process, i.e., is a Gaussian stochastic process characterized by in which i and j are particle indices, α and β are xyz spatial indices, and D and τ are constants, respectively, of dimension length squared over time and time.We are interested in how the physics is affected when the density is changed, specifically in determining whether approximately invariant physics can be obtained by adjusting D and τ properly with density (µ is regarded as a material constant throughout).For the binary AOUP model this problem was studied in Ref. 32 that demonstrated how to change D and τ with density in order to achieve invariant structure and dynamics to a good approximation.The question is whether this is possible also for systems with large size polydispersity.In the AOUP model "reduced" quantities are defined by using l 0 = ρ −1/3 as the length unit and t 0 = τ as the time unit [32].Reduced quantities are marked by a tilde.When we speak about approximately invariant structure and dynamics, it refers to this particular state-point-dependent unit system.We studied a system of N = 5000 particles in three dimensions interacting by Lennard-Jones (LJ) pair potentials, which between particles i and j are given by v ij (r) = 4ε (r/σ ij ) −12 − (r/σ ij ) −6 with σ ij = (σ i + σ j )/2 (Lorentz-Berthelot mixing rule) and ε = 1.0.The particles sizes σ i are distributed according to a uniform distribution with unity average.As usual, the polydispersity δ is defined by δ 2 = ( σ 2 − σ 2 )/ σ 2 , which in our case reduces to δ 2 = σ 2 − 1.For a uniform distribution δ cannot exceed 1/ √ 3 ∼ = 58%.The three polydispersities studied below are δ ∼ = 11.5%,23.1%, and 28.9%, corresponding to the size ranges listed in Table I (for brevity these are henceforth reported as δ = 12%, 23%, and 29%).Note that the study entails substantially different particle sizes, with the ratio of largest to smallest particle volume equal to 27 in the 29% polydispersity case.All simulations used a shifted-force cutoff [35] of the ij particle interaction at the pair distance r = 2.5σ ij and the time step ∆t = ∆ t/(D ρ 2/3 ) in which ∆ t = 0.4 [32].The active-matter simulations were carried out on GPU cards using a home-made code, the MD simulations used RUMD [36].

III. STRUCTURE AND DYNAMICS ALONG AN ISOCHORE
Before discussing results for the variation of structure and dynamics along active-matter isomorphs, we briefly present analogous results along an isochore, i.e., for state points of the same density.This sets the stage by illustrating that structure and dynamics do vary significantly throughout the (ρ, D, τ ) AOUP phase diagram.Structure is studied by means of the average radial distribution function (RDF) denoted by g(r).In Fig. 1(a) RDFs are shown along the ρ = 0.85 isochore for the δ = 29% case, with values of D and τ taken from the δ = 29% DIC active-matter isomorph studied below.Figure 1(b) shows the same data in reduced coordinates, which in this case simply involves in a common scaling of the x-coordinate.The parameters used in the simulations are listed in insets of the figures (more decimals of these parameters are provided in the Appendix).
We find a substantial structure variation along the isochore.The same applies for the mean-square displacement (MSD) as a function of the time t, ∆r 2 (t) , which is plotted in a log-log plot in (c) LJ units and (d) reduced units.The short-time slope is two, reflecting the "ballistic" regime of the AOUP model, which is not present in ordinary Langevin dynamics [16,[18][19][20] because it results from short-time noise correlations resulting in an inertia-like persistence of the direction of motion.At long times the well-known diffusive behavior leading to unity slope is observed.We note that the dynamics varies significantly along the isochore, whether or not reported in reduced units.show the RDF as a function of r and of the reduced pair distance r, respectively (the curves are the same because r ∝ r along an isochore).We see a substantial variation in the structure, with the most pronounced structure found for the smallest values of the model parameter D (black curves).The MSD likewise shows no collapse along the isochore, whether plotted (c) as a function of the time t or (d) as a function of the reduced time t.The slowest motion is found for the smallest D (black curves).

IV. STRUCTURE AND DYNAMICS ALONG T conf -GENERATED ACTIVE-MATTER ISOMORPHS
Reference 32 used the configurational temperature T conf for determining how to change the AOUP model parameters D and τ with density in order to achieve (approximately) invariant reduced structure and dynamics.The assumption is that k B T conf is the relevant characteristic energy scale where T conf is defined by k B T conf ≡ (∇U ) 2 / ∇ 2 U [37 -39] in which k B is the Boltzmann constant, ∇ is the gradient operator in the 3N -dimensional configuration space, and the sharp brackets denote standard canonical-ensemble averages.In the thermodynamic limit the relative fluctuations of both the numerator and the denominator of T conf go to zero, which implies that it is enough to consider a single configuration R 0 using the expression The reasoning of Ref. 32 may be summarized as follows.Adopting e 0 = k B T conf as the energy unit supplementing the above introduced length and time units (l 0 = ρ −1/3 ; t 0 = τ ), we first note that the three quantities µt 0 e 0 /l 2 0 , Dt 0 /l 2 0 , and τ /t 0 are dimensionless.Assuming that these quantities cannot change with varying density if the structure and dynamics are invariant in reduced units, we conclude that µ ∝ l 2 0 /(t Since µ is assumed to be a material constant, this leads to τ ∝ ρ −2/3 /k B T conf and D ∝ k B T conf , i.e., to the following recipe for how D and τ changes with density in terms of their values D 0 and τ 0 at a reference state point of density ρ 0 : For a large system T conf (ρ 0 ) may be evaluated from a single (steady-state) configuration, Reference 32 demonstrated that this approximation introduces a negligible error for typical system sizes.In order to find T conf (ρ) one scales R 0 uniformly to the density ρ, i.e., substitutes R = (ρ 0 /ρ) 1/3 R 0 into the configurational temperature expression.This leads to We used these equations for generating three active-matter isomorphs starting in each case from the parameter values D = 1100 and τ = 10 at the reference densities ρ 0 =0.99, 0.91, and 0.85, respectively, for the polydispersities δ = 12%, 23%, and 29% (the reference densities were chosen to have the same virial, i.e., give the same contributions to the pressure coming from the interactions).
Results for the variation of the average RDF are given in Fig. 2. The left column reports the RDF for the three polydispersities as functions of the pair distance r, the right column shows the same data as functions of the reduced pair distance r.In the latter case we find a good, but not perfect, data collapse and conclude that the average structure is approximately invariant along the active-matter isomorphs.In view of the fact that the density varies by no less than a factor of two, this is not trivial.Figure 3 shows analogous data for the MSD plotted in the same way with the left column giving the MSD as a function of time and the right column giving the same data in reduced units.There is a good data collapse with, however, a somewhat faster motion at the higher densities.There is a good, but not perfect, collapse of the reduced MSD along the active-matter isomorphs.

V. COMPARING TO DIRECT-ISOMORPH-CHECK GENERATED ISOMORPHS
Above we demonstrated good invariance of the structure and dynamics along active-matter isomorphs generated by the T conf method [32].That method is easy to use and efficient because it requires just a single steady-state configuration at the reference state point in order to trace out the corresponding active-matter isomorph in the relevant phase diagram, in casu the (ρ, D, τ ) diagram of the AOUP model.An alternative method for tracing out active-matter isomorphs is the analytical "direct isomorph-check" (DIC) method, which in Appendix A of Ref. 32 was shown to result in somewhat better isomorph-invariance of the dynamics for the AOUP Kob-Andersen binary LJ model.
Consider first a standard passive Newtonian systems involving LJ pair interactions of any kind, i.e., singlecomponent, binary, or polydisperse systems, defined by some mixing rule.For such a system the analytical DIC recipe for tracing out a standard equilibrium isomorph [40,41] is Here h(ρ) is the following function of density [40,41] in which ρ 0 is the reference-state-point density and γ 0 is the density-scaling exponent at the reference state point.The latter quantity can be determined numerically by means of in which ∆W and ∆U are the deviations from the equilibrium values of virial and potential energy, respectively, and angular brackets denote N V T equilibrium averages [31,42].
The systemic temperature T s (R) is defined as the temperature of the corresponding thermal-equilibrium Newtonian system at the state point with the density of the configuration R and average potential energy equal to U (R) [43].In the thermodynamic limit of any system (passive or active, equilibrium or non-equilibrium) fluctuations in T s (R) go to zero, implying that one has at any time a well-defined systemic temperature T s .For, e.g., a driven passive or an active-matter system, a "systemic isomorph" is defined as a curve in the (ρ, T s ) phase diagram identical to an isomorph in the standard equilibrium Newtonian (ρ, T ) phase diagram [43].Thus in the analytical DIC, the systemic-temperature's variation with density is given by i.e., T s (ρ) ∝ h(ρ).The analytical DIC method for generating an active-matter isomorph of the AOUP model is arrived at by replacing the configurational temperature in Eq. ( 3) by the systemic temperature T s (this procedure is justified in Ref. 32).Via Eq. ( 8) this leads to Table II reports the systemic temperatures at the reference state points of the three polydispersities studied.As mentioned, the reference densities were chosen to have the same virial; we see that they also have almost the same systemic temperature.PolyLJ AOUP, Active -Matter Isomorph FIG. 4. Active-matter isomorphs for the polydispersities δ = 12%, 23%, and 29%, generated from the reference state points by two different methods, the T conf method of Sec.III (full curves) and the analytical direct-isomorph-check (DIC) method (dashed curves).The isomorphs are visibly different.
Fig. 4 shows the active-matter isomorph obtained from the T conf method (full curves) and the analytical DIC method (dashed curves), starting at the reference state point (ρ, D, τ ) = (ρ 0 , 1100, 10) in which the reference density is 0.99, 0.91, and 0.85, respectively, for δ = 12%, 23%, and 29%.The two methods for generating isomorphs result in visibly different curves; thus there is more than 50% difference in D and τ at the largest density in the 29% polydispersity case (green curves).How different are these active-matter isomorphs when it comes to average RDF and MSD data collapse?The RDF case is investigated in Fig. 5, which shows that the structure is somewhat more invariant along the T conf -generated active-matter isomorphs than along the DIC-generated isomorphs, albeit this is a minor effect because in both cases the structure is fairly invariant.The differences are most pronounced at higher polydispersity.
FIG. 5. Comparing the degree of structural invariance along T conf -generated and DIC-generated active-matter isomorphs.(a), (c), and (e) show the reduced average RDFs for polydispersity δ = 12%, 23%, and 29%, along the T conf -generated isomorphs, while (b), (d), and (f) show the corresponding reduced average RDFs along the DIC-generated isomorphs.There is a somewhat better data collapse along the T conf -generated isomorphs.
Figure 6 reports results for the MSD.Here we reach the opposite conclusion: the DIC method results in a somewhat better data collapse than the T conf method.The same conclusion was reached for the binary Kob-Andersen AOUP model in Ref. 32 (that did not investigate the average RDF).The invariance of the largest 20% particle RDF is much better than the small-particle RDF.For the MSD, however, both cases are isomorph invariant to a good approximation.
To illuminate why the structure is not always isomorph invariant, we studied for the 29% polydispersity simulation data the structure and dynamics of the 20% smallest and largest particles along the DIC isomorph (Fig. 7).The RDF is here defined by limiting the central particle to be either among the smallest or among the largest 20% and counting only neighboring particles of the same kind.In both cases, the peaks are narrower and higher than that of the average RDF (Fig. 5(f)), which reflects the limitation to similar-sized particles.We note that unlike in Fig. 5(f), in Fig. 7(a) the most pronounced structure is seen at the lowest density.Interestingly, the structure around the smallest particles is not DIC-isomorph invariant, while that around the largest particles is; a similar result applies for the T conf -generated isomorph (data not shown).In contrast to the findings for structure, the dynamics of both small and large particles is DIC-isomorph invariant to a good approximation, even though smallest particles move considerably faster than the largest ones.We conclude that the lack of perfect isomorph invariance of the overall RDF largely reflects the fact that the structure around the smallest particles is not isomorph invariant.

VII. SUMMARY AND OUTLOOK
We have shown that the uniform-distribution size-polydisperse LJ AOUP model has active-matter isomorphs for the polydispersities δ = 12%, 23%, and 29%.This demonstrates the robustness of the active-matter-isomorph concept, which for passive systems applies whenever the potential-energy function obeys the hidden-scale-invariance condition discussed in Refs.44 and 45.The existence of isomorphs means that the dimension of the polydisperse AOUP phase diagram is effectively reduced from three to two, since it implies that lines exist in the (ρ, D, τ ) phase diagram along which the reduced structure and dynamics are invariant to a good approximation.From a practical perspective, this fact makes it easy to quickly get an overview of the AOUP model's phase diagram.
Two methods have been studied for generating active-matter isomorphs, one based on the configurational temperature and one based on the systemic-temperature concept.We find that both methods work well despite the fact that they do not trace out identical active-matter isomorphs (Fig. 4).In practice, the latter method will be easier to use in the case of LJ active matter for which a simple expression is available for the function h(ρ) where the parameter γ 0 may be evaluated from a single passive-matter simulation (Eq.( 7)).
More work is needed to clarify how polydispersity relates to the existence of active-matter isomorphs in general.As regards the AOUP model, it would be interesting to investigate whether the introduction of energy polydispersity affects the existence of isomorphs.More generally, other models like the active Brownian particle model with a potential-energy function that obeys hidden scale invariance should be investigated in polydisperse versions in order to illuminate the robustness of the active-matter-isomorph concept.

Appendix
This Appendix provides more decimals than given in the figures of the parameters ρ, D, and τ used in the simulations.It follows from Eq. ( 4) and Eq. ( 9) that both methods of tracing out active-matter isomorphs result in Dτ ρ 2/3 =Const.Within the number of decimals given, this applies for the below reported parameters.S3.Values of ρ, D and τ for polydispersity δ = 23.09%along the T conf isomorph (Fig. 2 and 3).S7.Values of ρ, D and τ for polydispersity δ = 28.87%along the DIC isomorph (Fig. 5 and 6).

FIG. 1 .
FIG.1.Average radial distribution function (RDF) and mean-square displacement (MSD) for state points on the ρ = 0.85 isochore of the δ = 29% polydispersity LJ AOUP model (the D and τ values are those of the below studied δ = 29% activematter DIC isomorph).(a) and (b) show the RDF as a function of r and of the reduced pair distance r, respectively (the curves are the same because r ∝ r along an isochore).We see a substantial variation in the structure, with the most pronounced structure found for the smallest values of the model parameter D (black curves).The MSD likewise shows no collapse along the isochore, whether plotted (c) as a function of the time t or (d) as a function of the reduced time t.The slowest motion is found for the smallest D (black curves).

FIG. 2 .
FIG.2.Structure probed along T conf -generated active-matter isomorphs.(a), (c), and (e) show the average RDFs for polydispersity δ = 12%, 23%, and 29%, respectively, while (b), (d), and (f) show the same data as functions of the reduced pair distance r.In all three cases we see a good collapse of the reduced RDF along the active-matter isomorph.

FIG. 3 .
FIG.3.Dynamics probed along T conf -generated active-matter isomorphs.(a), (c), and (e) show the MSDs for polydispersity δ = 12%, 23%, and 29%, respectively, as functions of time, while (b), (d), and (f) show the same data in reduced units.There is a good, but not perfect, collapse of the reduced MSD along the active-matter isomorphs.

FIG. 6 .FIG. 7 .
FIG.6.Comparing the degree of invariance of the dynamics along the T conf -generated and DIC-generated active-matter isomorphs.(a), (c), and (e) show the reduced MSDs for polydispersity δ = 12%, 23%, and 29%, along the T conf -generated isomorphs, while (b), (d), and (f) show the corresponding reduced MSDs along the DIC-generated isomorphs.There is a better data collapse along the DIC-generated isomorphs.

TABLE I .
Values of the polydispersity δ, σ range, and ratio between largest and smallest particle sizes for the three cases of uniform polydispersity studied.