General scaling in bidirectional flows of self-avoiding agents

The analysis of the classical radial distribution function of a system provides a possible procedure for uncovering interaction rules between individuals out of collective movement patterns. A formal extension of this approach has revealed recently the existence of a universal scaling in the collective spatial patterns of pedestrians, characterized by an effective potential of interaction \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V(\tau )$$\end{document}V(τ) conveniently defined in the space of the times-to-collision \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\boldsymbol{\tau }}$$\end{document}τ between the individuals. Here we significantly extend and clarify this idea by exploring numerically the emergence of that scaling for different scenarios. In particular, we compare the results of bidirectional flows when completely different rules of self-avoidance between individuals are assumed (from physical-like repulsive potentials to standard heuristic rules commonly used to reproduce pedestrians dynamics). We prove that all the situations lead to a common scaling in the t-space both in the disordered phase (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V(\tau ) \sim {\tau }^{-2}$$\end{document}V(τ)~τ−2) and in the lane-formation regime (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V(\tau ) \sim {\tau }^{-1}$$\end{document}V(τ)~τ−1), independent of the nature of the interactions considered. Our results thus suggest that these scalings cannot be interpreted as a proxy for how interactions between pedestrians actually occur, but they rather represent a common feature for bidirectional flows of self-avoiding agents.

The theoretical description of pedestrian flows represents a field of the greatest importance due to its direct impact on issues related to urban planning, monitoring of public spaces or optimization of evacuation protocols, to name a few [1][2][3][4] . In the last years, facilitated by an improvement in simulation capacities and in the availability of experimental data, this topic has attracted physicists for its interest as a multiagent system driven by non-trivial rules of ordering, alignment and self-avoidance, among other. Thus, a significant effort has been put both in (i) understanding these interaction rules in order to recover the patterns experimentally observed 5,6 and (ii) identifying the minimal or toy models which are able to capture the essentials of such phenomena 7,8 .
Lane or trail formation in crowds, in particular, has been extensively studied as a manifestation of self-organization in pedestrian flows, both theoretically [9][10][11][12][13] and through controlled experiments 14,15 . While in the case of humans one could attribute lane formation to intelligent and efficient decision-making based on visual information and subsequent prospection, such patterns have been shown to arise even in extremely simple models of agent interaction. The Vicsek model for swarming dynamics 16,17 is probably the best known, but many variations and generalizations based on lattice-gas 18,19 or social force 20,21 models have been developed. Additionally, the formation and sustainance of bidirectional (trail) flows represent also an intriguing situation of interest in behavioral biology too, since only a few social species in the animal kingdom are able to exhibit such behavior in natural conditions 22 , ants being the most renowned case 23,24 .
Despite all this multidisciplinary interest, the heterogeneity of models used nowadays to generate such flows/ dynamics sometimes goes against the possibility of finding general and far-reaching conclusions. So, existing models/works can sometimes provide different or even contradictory conclusions 15 . Works aimed at providing unified frameworks and/or at revealing the common properties of these approaches are then convenient to promote understanding and theoretical research within the field 25 .
Within this context, a valuable insight has been recently provided by Karamaouzas et al. 26 . By analyzing 1500 trajectories of pedestrians in outdoor environments they found consistent evidence for an effective potential of interaction τ V( ) between pairs, which was found to depend only on the time-to-collision τ between the individuals and not on the interparticle distance as in classical fluid systems. This is in line with recent approaches that introduce prospection of future outcomes as the main driving force for intelligent agents, as those based on causal entropic forces 27,28 . The analysis carried out in 26 2 , at least in the significant range of scales where self-avoidance between the pedestrians is relevant. So, a mechanical (Langevin-like) model based on this interaction potential could adequately reproduce the dynamics of pedestrians in real scenarios.
Intrigued by these results, we use here numerical simulations to explore the convenience of a description based on the τ-space for different bidirectional systems based on radically different rules of self-avoidance. In general, bidirectional flows emerge from the tension between pair interactions and the existence of two subpopulations that have preference for moving in opposite directions; when the density of individuals ρ is large enough the freedom of movement gets reduced and interactions will be then largely governed by self-avoidance rules.

Agent based model
We consider in the following a multiagent system where the dynamics of the (identical disk-like) agents are governed by a force Here, the first term accounts for the preference of each individual to maintain its preferred velocity v pr ( ) (denoting v as the actual velocity) with a certain intensity that we call the stubbornness, ξ. In order to generate bidirectional fluxes we assume two different subpopulations with the same number of agents each, whose preferred velocities are the same in modulus but have opposite directions. In particular, to avoid spuriors effects in the simulations we sample for each agent a stochastic preferred speed from a Gaussian distribution with mean = .  29,30 . Anyway, we have numerically checked that the results reported below are independent of the specific values chosen.
On the other hand, F sa ( ) stands for the pair (self-avoiding) interaction between agents. We here compare the results for three rules/mechanisms based on completely different grounds. The first one consists of a classical pair repulsive interaction in the radial direction, where r is the distance between pairs of individuals), with > k 0. The second one corresponds to the effective potential empirically obtained in 26 , this is, a repulsion in the time-to-collision (ttc) space, These two mechanisms then are reminiscent of typical interaction potentials from statistical mechanics, though in the second case we are considering that the relevant space for interactions is that of τ, so assuming that our intelligent agents can somehow prospect future collisions and adapt its behavior to the outcome of such prospections. So, the space of interactions is defined only for those individuals for which τ is finite, this is, for the set of pairs that will collide at some future time provided its present velocity is mantained (see Fig. 1). Finally, as a third case, F sa heu ( ) , we consider a nonphysical (heuristic) rule which has been found to reproduce most features (e.g stop-and-go waves, turbulent dynamics,...) of collective behavior in pedestrians 31 . This rule consists of recomputing continually the direction of motion in order to maximize at each step the distance that the agent could travel without colliding with other agents (see the Methods Section for the finer details for the implementation of the self-avoiding rules).
Numerical implementation of our multiagent system identifies, as expected, the existence of a phase transition from a disordered state to lane formation as a function of the values of ρ and ξ. To characterize this transition we use the order parameter cos( ) , (2) f θ = where θ is the angle between the actual and the preferred velocities (this is, v and v pr ( ) ) and the average is carried out over all the agents in the system. So, f → 0 corresponds to the disordered state in which individuals cannot Figure 1. Representation of the differences between the interactions based on the distance r (left) and those in the τ-space (right). The arrows represent the relative velocity of the agents respect to the orange agent at the origin. The individuals filled with solid blue are the only ones contributing to interactions with the orange one (it is, in the second case only those for which a finite and positive τ can be defined).
follow its preferred direction and spend their time avoiding collisions in all directions, while 1 f → represents the case where the agents are able to follow its desired direction of motion by adopting a collective pattern with alternate lanes in one direction and the other. Figure 2 confirms that the three mechanisms of self-avoidance exhibit qualitatively a very similar behavior for the order parameter f at the stationary state. Low stubbornness and high densities promote the disordered phase, while for high stubbornness and low densities the system self-organizes into a new phase where lanes are formed (we provide in the Supplementary Information short videos showing the dynamics observed in each of these two phases for the three self-avoidance mechanisms mentioned above). In the disordered state the difference between the actual direction of motion and the preferred one is relatively homogeneous in π (0, ). We can plot the probability distribution θ p( ) of angle θ to visualize this (Fig. 3, left). For the case of lanes, on the contrary, most of the individuals move in their desired direction and then the probability distribution becomes clearly peaked at θ = 0. Additionally, we observe how the heuristic mechanism exhibits a larger probability for large deviations in the lane state than the other interactions; this is due to the intrinsic properties of the algorithm, which allows larger reorientations provided they satisfy the maximization of the traveled distance, as explained above.

Results and Discussion
To understand in greater detail the properties of these two phases, we next study the spatial distribution of the agents in stationary conditions through the radial distribution function g r ( ) from the simulations carried out using the three self-avoiding mechanisms mentioned above. As in classical fluids, g r ( ) here compares the density of interacting agents at a distance r with the density obtained for a non-interacting system, with . The corresponding results are presented in Fig. 4. Despite some qualitative differences found due to the different nature of the self-avoiding mechanisms, we observe that the results are relatively consistent with those  www.nature.com/scientificreports www.nature.com/scientificreports/ from the classical Ornstein-Zernike (OZ) approximation for fluid systems 32 , which predicts an asymptotic decay ( ) 1) when the system is far from the critical region where the phase transition occurs. The oscillatory behavior observed there is also characteristic from similar statistical analysis on fluids 33 . All the points of the phase diagram that are far from the critical region satisfy approximately this scaling while those close to the critical region (separating disorder from lane formation) exhibit a much slower asymptotic decay. In consequence, g r ( ) cannot be easily used to detect or identify the particular state (disordered or lanes) in the system.
Going further, we reproduce the analysis in 26 by showing how g r ( ) gets modified if the data is split into three parts according to the relative speed between pairs of individuals i and j, relative speeds exhibit very different behaviors in all cases (Fig. 5, left column). However, we stress that, at least for the repulsive potential F rep sa for which interactions are not velocity-dependent, the function g r ( ) should be also independent of the actual velocity of the particles; this is confirmed in our simulations (see Suplemmentary Information file) for consistency.
The authors in 26 concluded that the differences observed in g r ( ) for different values fo v r , reflects that such function is not a very appropriate descriptor for capturing the effective interactions within the crowd or, stated in different words, the collective statistics of the system do not apparently yield a common behavior in the physical r space of the distances between individuals.
To explore here this idea we introduce a new magnitude * g r ( ), defined as the radial structure function but only for pairs of colliding agents, which are those for which τ is finite at that time step (Fig. 1). The splitting of * g r ( ) into different values of the relative velocities still shows that the results are strongly dependent on v r (Fig. 5, middle column), albeit the differences get reduced for F sa ttc ( ) and F sa heu ( ) (since these two interaction rules only apply to particles which are about to collide). Instead, for the rule F sa rep ( ) , which applies to all pairs of particles, the results found are almost the same as those for g r ( ). So, again collective statistics seem to depart from such descriptor. Our intuition, however, gained from the results in 26 is that the dynamics of self-avoidance should rather translate into a robust behavior within the τ-space, as the events with low τ are the ones which must be avoided first. So, we finally introduce τ † g ( ), which is the equivalent to g r ( ) but on τ-space, i.e. the density of agents found at a time-to-collision τ divided by the density we would find at the same τ for the case of non-interacting agents. The corresponding results are shown in Fig. 5 (right column). The idea that interactions should occur in the τ-space is of course introduced by hand in our rule F sa ttc ( ) , and also implicitly in the rule F sa heu ( ) , so when we explore the dynamics in the τ -space then we observe that the collapse between the three curves (for low, intermediate and high relative speeds) is almost perfect. However, we unexpectedly find that the collapse between the curves is moderately improved for F sa rep ( ) too, though this self-avoiding rule has nothing to do with τ. This suggests the existence of an underlying phenomena enhancing the relevance (at least at the level of how collective structures emerge) of the τ-space whenever self-avoidance and bidirectionality effects drive the system dynamics. In the Supplementary Material we carry out an alternative analysis of the radial distribution functions in order to provide additional support for this idea. Note that our results in Fig. 5 do not necessarily mean that pair interaction occurs in the τ-space (which is not the case for our repulsive potential, actually) but that at a collective level this is the effective situation produced. www.nature.com/scientificreports www.nature.com/scientificreports/ Next step is to derive an effective potential of interaction between agents in the τ-space. For the classical theory of fluids, the reversible work theorem 34 in the r space links the radial distribution function g r ( ) with interaction energy between pairs in the form ∝ V r g r ( ) ln[ ( )]. Using an analogy with this classical result, the τ-space also admits an equivalent derivation. Such derivation is based on the idea that the system satisfies in the τ-space a Boltzmann-like statistics. This cannot be justified from classical statistical mechanics since the concept of thermal equilibrium is meaningless in our context, but one can still invoke the Maximum Entropy Principle, which has solid foundations from information theory, to justify it at a statistical level 35 \cite{ . So that, the corresponding expression τ τ ∝ † V g ( ) ln[ ( )] must be interpreted as a statistical relation describing average properties in the τ -space; note, however, that τ V( ) does not represent the real potential of interaction between particles, and so τ ∇ V( ) r is not to be interpreted as a physical force. The effective potential τ V( ) obtained from our model is presented in Fig. 6, which represents the main result of our work. Surprisingly, we find that the three self avoiding mechanisms collapse for intermediate times-to collision (which is the significant region where most of the pair-pair interactions occur) into a common power-law relationship , with γ ≈ 2 for the disordered state and γ ≈ 1 for the state with lanes (in Table 1 we show the results obtained from fitting the curves presented in Fig. 6). This common scaling is then apparently independent of the self-avoiding mechanism, and would be a direct consequence of the bidirectional nature of the flow considered. Note also that, contrary to what happened in the r space (Fig. 4), the τ † g ( ) and the corresponding τ V( ) show a different decay for the disordered phase and the case with lanes, so τ † g ( ) can be effectively used to identify these two states. www.nature.com/scientificreports www.nature.com/scientificreports/ According to Fig. 3 (left), in the disordered state collisions can be produced in any orientation and the events corresponding to large τ are supressed by the shielding of closer events. The decay exponent γ takes then a value of 2, in agreement with the power-law proposed for pedestrians in 26 . While this is to be expected in the time-to-collision interaction F sa ttc ( ) by definition, there is no apparent reason to justify why the same behavior emerges for F sa rep ( ) and F sa heu ( ) . On the other side, the effective potential in the lane state exhibits a completely different behavior. The homogeneity in θ is there completely broken (Fig. 3, right) due to the two preferred directions of movement and most of the collisions are produced in the frontiers between opposite lanes. There is not shielding effect in this case and, as a consequence, the system is driven by a slower interaction decay, with γ ≈ 1. We have not been able, however, to find an analytical justification for the specific values of γ emerging in each state; this remains an open question.
To conclude, from our findings in Fig. 6 we obtain that the same effective potential, if computed through Eq. ?? as done here, could emerge from a very wide range of interactions between the agents. This suggests that the result τ τ − V( ) 2 experimentally reported in 26 is not necessarily determining the actual rule of interaction (or self-avoidance) used by pedestrians, but it could rather be the manifestation of a common dynamics exhibited by a wide range of systems combining self-avoidance and bidirectionality. In particular, our results in Fig. 5 confirm that it is not possible to discern whether pedestrians use a time-to-colision potential (as in 26 ) or a heuristic rule of path maximization (as in 31 ) only from examining the shape of the distribution function τ g( ), but additional analysis would be required. Still, the scaling τ τ γ − V( ) will presumably work as a useful effective rule in bidirectional flows for different situations of interest (either pedestrian movement, ant organization or complex plasma 36 , to cite some known cases). Such effective rule could be of great utility in order to simulate bidirectional fluxes without caring too much about the fine details of the interactions, and so it can be used as a toy or reference approximation to computational or analytical approaches in the field.
Our work is not able to delimit what is the exact range of validity of the scaling reported in Fig. 5. Given the strong differences between the three self-avoiding mechanisms explored here, we suspect that this range can be quite large as long as the initial scheme in Eq. 1 holds. However, this does not guarantee the existence of an equivalent effective potential in situations where additional forces get introduced. Hence the extension of our results in this direction can represent an stimulating area of future research.

Appendix A: Methods
Self-avoidance mechanisms. As we explained before, we propose in our work three very different mechanism for the self-avoidance potential F sa ( ) between the agents (which are considered as disk-like particles with the same diameter D). Here we provide a short description of each mechanism. www.nature.com/scientificreports www.nature.com/scientificreports/ Repulsive potentials. We introduce a repulsive pair energy which is a function of the distance r between the agents, where > A 0, u is a unit vector in the direction joining the pair of interacting agents, and the distance r is measured in units of D, so = r 1 corresponds to the distance between two adjacent agents. The parameter k regulates the decay of the force, so implicitly it determines the range of scales where its effect is relevant, with the limit → ∞ k , reproducing hard-disk interactions. This approach represents a reference model against which we subsequently compare the performance of more sophisticated mechanisms of self-avoidance between pedestrians. In this case the effect of F rep sa ( ) is to pull all the individuals apart, independently if their are moving forward to a collision or not.
Time to collision. The second potential used corresponds to an interaction explicitly occurring in the time-to-collision (or τ) space, defined as the time needed for the pair of agents to get in contact provided they followed their actual direction of motion. This time can be explicitly defined in terms of the relative velocity v r and the relative position r between two given particles, The idea of using the τ-space for driving interactions is directly borrowed from 26 and is justified from the empirical results on pedestrians dynamics therein obtained. The rule reads then , and u defined again as in 3. The exponential term is used as a cutoff to impair the effect of outermost collisions, so introducing the idea that agents possess a characteristic radius of perception (τ 0 , defined in the τ-space). The potential so defined only applies to agents moving forward to a collision, such that τ can be defined and is positive; this is, pairs for which a positive value of τ cannot be found are considered as noninteracting agents. As a result, the set of agents interacting with a given one is a dynamic object which is updated continuously throughout the simulation.
Heuristic rule. In order to consider very disparate mechanisms, we finally introduce a self-avoiding heuristic rule proposed for pedestrians dynamics in 31 . This algorithm relies on the idea of adapting the direction of motion by maximizing locally the accessible distance path. So, each agent samples its possible future trajectories by simulating internally (with a time horizon t m ) where it will reach by moving in a given direction (characterized by an angle α) for some fixed time, provided that the other agents are assumed to go on moving in the same direction they do have at present. After sampling for a range of values of α (up to maximum α max , to avoid sudden or extreme changes of direction) the agent will choose the one that maximizes the length covered. After the election, all the agents reorient synchronously and the internal simulation starts anew.
There is a second rule, which determines the walking speed modulus after the reorientation. This is introduced in order to maintain a certain time to collision between the agent and the obstacle in the chosen walking direction 31 . For this, we define a minimum time τ m such that times-to-collision are forced to stay always below τ m by reducing adequately the speed of the agents. That speed is then computed at practice as where d obs is the distance between the agent and the first obstacle in the desired direction α at that time step.
Implementation and technical details. In this section we provide additional details of how the simulations of our multiagent model have been carried out. For the first two interactions, the number of agents is fixed to = N 512, while for the heuristic rule the number is fixed to = N 128 due to the computational cost of its simulations (different time steps ∆t are also used in each case for the same reason, see below). The simulation time for the repulsive and time-to-collision mechanisms (it is, F rep sa ( ) and F ttc sa ( ) ) scales as ∝N 2 , as they are pair to pair interactions. Instead, the heuristic rule prospects into the future the different α paths. This algorithm implies a scaling time ∝ The agents are placed in a two-dimensional simulation box with density ρ using periodic boundary conditions. The units of length are re-scaled to σ, so = r 1 is equal to a diameter agent. The agent mass is settled as = m 1. The Verlet algorithm has been used to integrate the equations. The system is studied for different values of the density in the range ρ = .
. [0 05, 0 32], which is accomplished by fixing the number of individuals to a certain value N and changing the box size L, given ρ = N L / 2 ). The parameters used for the implementation of the self-avoidance mechanisms are as follows: www.nature.com/scientificreports www.nature.com/scientificreports/ • The repulsive interaction(3) is fixed to = k 4 (unless indicated otherwise), = .
k 1 5, τ = 10 0 and ∆ = . t 0 005 ttc according to 26 . The τ 0 value is defined in order to not affect the dynamics in the scaling region.
• The heuristic rule is fixed to τ = . To carry out the simulations 8 CPUs have been used, with a total simulation time around ~250 hours for each CPU.