Congestion Transition on Random Walks on Graphs

Abstract

The formation of congestion on an urban road network is a key issue for the development of sustainable mobility in future smart cities. In this work, we propose a reductionist approach by studying the stationary states of a simple transport model using a random process on a graph, where each node represents a location and the link weights give the transition rates to move from one node to another, representing the mobility demand. Each node has a maximum flow rate and a maximum load capacity, and we assume that the average incoming flow equals the outgoing flow. In the approximation of the single-step process, we are able to analytically characterize the traffic load distribution on the single nodes using a local maximum entropy principle. Our results explain how congested nodes emerge as the total traffic load increases, analogous to a percolation transition where the appearance of a congested node is an independent random event. However, using numerical simulations, we show that in the more realistic case of synchronous dynamics for the nodes, entropic forces introduce correlations among the node states and favor the clustering of empty and congested nodes. Our aim is to highlight the universal properties of congestion formation and, in particular, to understand the role of traffic load fluctuations as a possible precursor of congestion in a transport network.

1. Introduction

Modeling city mobility is a crucial aspect of planning future mobility infrastructures and implementing governance policies for sustainable mobility in future smart cities [1,2]. A comprehensive microscopic model to simulate urban mobility must delineate individual mobility demand [3,4], simulate the decision-making mechanisms of individuals to determine mobility strategies, and consider the effects of physical interactions on transport networks [5].

This formidable task requires an enormous amount of data to set up the model parameters and sophisticated methodologies to analyze the simulation results and detect the control parameters of the system. Modeling traffic dynamics has revealed some universal features, such as stop-and-go congestion when vehicle density surpasses a certain threshold [6], which do not depend on the specifics of vehicle interaction. However, the emergence of congestion in an urban transport system is a different phenomenon [7,8], as the dynamics at intersection points become more significant than the dynamics on the roads. Previous studies have proposed percolation theory as a key to understanding the emergence of congestion in a traffic network [9,10,11]. The very definition of congestion on a road network can be approached from different perspectives: individual and macroscopic. The individual perception is mainly influenced by changes in the average velocity of mobility paths, and the uncertainty of travel time can be fundamental in understanding individual behavior in urban mobility [12,13].

Conversely, the macroscopic approach considers the performance of the entire transport system by introducing a macroscopic fundamental diagram [14] and modeling the network dynamics. In this work, we propose a reductionist approach to model transport network dynamics and study congestion formation under stationary conditions by employing a random process on a graph, where each node represents a location, and the link weights (i.e., the transition rates to move from one node to another) are related to the statistical distribution of mobility paths on the network [15] and individual mobility demand. Our aim is to highlight the universal features of congestion transition and the role of traffic load fluctuations.

To this end, we introduce a reductionist model using a random walk on a graph based on two assumptions shared by queue models of transport networks: the existence of a maximum flow rate and a maximum capacity for a road. In the thermodynamic limit, where both the number of nodes and particles tend to infinity but with a fixed ratio, fluctuations remain finite, and an average field approach is not suitable to describe the macroscopic evolution of the system. To highlight the effects of traffic fluctuations, we explicitly study a transport network in a balanced condition where the average incoming and outgoing flows at each node are equal. This assumption reflects the existence of a Wardrop equilibrium [16,17] for urban traffic, where mobility paths distribute to optimize the use of the transport network.

The existence of a maximum flow rate and a finite capacity for each node implies that any displacement involves a maximum number of particles and is possible only if the number of particles in the destination nodes is smaller than the maximum capacity [18,19].

Our main assumption regarding the dynamics is that, in a stationary state, traffic fluctuations can be modeled by a Markov process [20,21]. We recall indeed that a coarse-grained description of a chaotic dynamical system justifies a stochastic approach. The Markov property, indicating short correlations in traffic fluctuations among connected nodes, is likely justified if the mobility demand is distributed across the urban fabric due to city complexity [2].

In the approximation of the single-step process, we can characterize the traffic load distribution of the nodes according to a local maximum entropy principle [22,23]. The application of the entropy concept to understand the statistical laws of urban mobility has been proposed in previous works [24,25,26] using a big-data approach. We extend this result to show how entropy principles allow the study of properties of non-equilibrium states near a stationary state of the transport network model [27,28,29].

We are able to study how the congestion transition can be detected from the fluctuation statistics of the node traffic loads, whose variance reaches a maximum when a peak at the congested nodes appears in the traffic load distribution, but there is no singularity in the thermodynamic limit. Using numerical simulations, we show that macroscopic congestion in the network, due to the emergence of a congested macroscopic cluster, can be accurately explained by a percolation transition, where the appearance of a congested node is an independent event. However, in the more realistic case of synchronous dynamics for the transport network, we demonstrate the appearance of entropic forces [30] that tend to cluster empty and congested nodes, increasing the size of the congested clusters even before the formation of macroscopic congestion. The existence of small congested clusters can introduce significant variance in the travel time distribution for individual paths, and our results suggest that this distribution can be used to characterize congestion formation before the percolation transition.

The paper is organized as follows: In Section 2, we discuss how to use nonlinear stochastic Markov systems as dynamical models of a transport network near a stationary state. In the Section 3, we illustrate the properties of random walks on graphs as models of transport networks and introduce an entropy-based approach to study traffic load fluctuations. In Section 4, we study congestion formation for a simple transport network model. In Section 5, we compare the analytical results with numerical simulations. Finally, we present some conclusions and perspectives.

2. Methods: Modeling Transport Network as Markov Processes

Urban traffic is the consequence of individual mobility demand to move from origins to destinations (OD) [3,4]. However, the complex structure of modern cities [31] has ubiquitously distributed activities across the urban fabric, making it extremely difficult to model human mobility using an origin–destination paradigm without collecting information at the individual level [5]. Moreover, the realization of mobility paths in the transport network [7,8] results from both physical interactions (traffic dynamics) and unpredictable individual decisions (free will) in route choice. The statistical physics approach offers a potential solution by using stochastic dynamical models related to a coarse-grained description of traffic dynamics, particularly when the macroscopic and mesoscopic states of traffic do not depend on the details of individual dynamics. The unpredictable features of individual mobility demand justify the use of random dynamical models to simulate mobility on transport networks and the application of the maximum entropy principle (MEP) [22,23] to study the properties of stationary states and the congestion transition. The existence of a mobility energy (i.e., travel time budget) is consistent with available data on individual mobility and suggests that mobility in a homogeneous transport network is characterized by an exponential path length distribution [24,25]. These observations suggest that a stochastic model for urban traffic can accurately describe the statistical properties of stationary states but may be inadequate for modeling a system out of equilibrium, where the complexity of urban mobility related to individual behavior could emerge.

The congestion formation on an urban road network depends on two main local features of traffic dynamics: the existence of a finite traffic flow rate at intersections [18] and of the maximal road capacity. The finite flow rate causes the onset of queues at intersections when the incoming flow increases and the maximal road capacity creates gridlock. A microscopic approach that simulates the dynamics of the single mobility paths requires a significant amount of information on the individual mobility demand, or it is based on the synthetic data obtained according to a priori assumptions (e.g., using optimized algorithms to compute the mobility paths from the OD mobility demand). Data on the individual mobility paths are difficult to obtain due to privacy concerns, and the use of synthetic data has the problem of not intruding bias in the system. Being interested in the reconstruction of the traffic flow dynamics along the roads, the necessary information is the measure of the probabilities, ${\pi }_{ij}$, to observe a path moving between the roads, $j\to i$. The transition probabilities, ${\pi }_{ij}$, do not require the whole knowledge of the path distribution, but they could be used to directly measure at the intersections. Then, by using a continuity argument one can compute the average traffic flows and highlight the existence of critical situations. However, if traffic distributes on the road network to avoid critical conditions (i.e., we have a balanced condition between average incoming and outgoing flows on each road), then the congestion can be triggered by traffic fluctuations when the traffic load is near to a critical value [32]. To study this effect, one introduces a stochastic process at the intersections that distributes the flows according the transition probabilities, neglecting the correlations among the not directly connected roads. The basic assumption is that each individual could be considered an independent particle and the correlation of traffic fluctuations decays rapidly with distance. This condition appears to be satisfied in modern cities, where mobility has strong random components [33]. If the only knowledge on the traffic dynamics is represented by the transition probabilities, ${\pi }_{ij}$, and one simulates the particle dynamics according to a random walks on graphs, the Markov condition is consistent with a maximum entropy principle, since this is the stochastic process that maximizes the information entropy of the particle trajectories (i.e., codes by the sequence of nodes) distribution with the constraints that the conditional probabilities, ${\pi }_{ij}$, for a single evolution step are given [34]. Therefore, random walks on graphs are the models that simulate traffic dynamics, maximizing the information entropy when the available knowledge is the transition probability matrix, ${\pi }_{ij}$, at road intersections. We also remark that the stationary distribution is completely determined by the transition probabilities so that the random walk models are able to characterize the onset of congestion and the stationary congested road distribution due to traffic fluctuations when traffic load increases. However, the traffic evolution during transient states may depends on long range spatial correlations. The application of a random walk model to study transient states is justified only if the microdynamics is strongly chaotic. Since the daily traffic load has strong periodic components, it is reasonable to assume that each individual tries to optimize their daily mobility cost. Then, one expects that mobility paths distributed on the transport network to create a Wardrop equilibrium with respect to travel time cost. Moreover, under normal conditions, mobility paths can be considered independent of traffic load (i.e., individuals do not change their mobility strategies).

Using a reductionist perspective, we introduce a mathematical description of a transport network using a graph where the nodes $i=1,\dots ,M$ represent either roads or intersections (in our case), and sometimes bus/train stations or generic transport facilities, and the links distribute the traffic load. Each node i is characterized by an internal state $n_i$ (the traffic load), and the flow ${\mathrm{\Phi }}_{ij}(n_i,n_j)$ between node j and node i depends on the node states (Markov field model) [20]. The existence of a maximum flow rate implies that the total flow, ${\mathrm{\Phi }}_j$, outgoing from the j node satisfies

$$\sum _i{\mathrm{\Phi }}_{ij}={\mathrm{\Phi }}_j\le {\mathrm{\Phi }}_j^{max}.$$

Let $\Delta t$ be the evolution time scale. If, for a given traffic load, $n_j\le {\mathrm{\Phi }}^{max}\Delta t$, we have free flow dynamics; conversely, we have queue formation at node j, and the crossing time increases proportionally to the queue length (i.e., the traffic load $n_j$). Queue formation depends on the total traffic load on the transport network, resulting in a simple fundamental diagram, as the average velocity, $\overline{v}$, decreases with the total traffic load, N. When the load, $n_j$, of a node approaches the maximum capacity, $n^{max}$, unavoidable traffic fluctuations can induce gridlock, drastically reducing incoming flows. In this case, the travel time may increase non-linearly with the total traffic load, as it depends on how the congestion spreads in the network [8]. Assuming the existence of a stationary state in a transport network, congestion formation is triggered by traffic fluctuations.

If there are no long-range spatial correlations in the mobility path distribution, the traffic dynamics at a stationary state can be approximated by a Markov random process with transition probabilities ${\pi }_{ij}$ (i.e., a non-linear random walk on a graph). The OD nature of urban mobility is not relevant for defining the stationary state when the correlation between the paths of two ’particles’ in the same node is negligible (i.e., the probability that two randomly chosen particles in the same node have similar destinations is very small). Moreover, one can justify the relaxation process to a stationary state and study the effect of fluctuations by approximating the dynamics with a Markov process [21,27]. The average dynamics of the transport network reads

$${\dot{n}}_i=\sum _j\left[{\mathrm{\Phi }}_{ij}(n_i,n_j)-{\mathrm{\Phi }}_{ji}(n_j,n_i)\right]+s_i\left(t\right),$$

(1)

where $s_i\left(t\right)$ represent the particle sources or the sinks present in the system that modulate the total traffic load. To study the stationary states we set $s_i=0$ and the traffic load ${\sum }_i{n}_i=N$ is constant. An average equilibrium solution satisfies

$$\sum _j\left[{\mathrm{\Phi }}_{ij}(n_i^{\ast },n_j^{\ast })-{\mathrm{\Phi }}_{ji}(n_j^{\ast },n_i^{\ast })\right]=0\sum _j{n}_j^{\ast }=N,$$

(2)

and the congestion transition occurs when the stationary solution becomes unstable and a new solution emerges with $n_j=n^{max}$ for a subset of nodes (the congested nodes). When the congested nodes form a giant cluster in the transport network, we say that the whole network is in a congested state.

Characterizing the Stationary State

The equilibrium states are determined by the stationary flows among the connected nodes. By definition, ${\pi }_{ij}$ is the probability that a particle performs the transition $j\to i$ and does not depend on the traffic load fluctuations. ${\pi }_{ij}$ is a stochastic matrix, and in a free traffic condition one can set ${\mathrm{\Phi }}_{ij}\propto {\pi }_{ij}{n}_j$. However, a maximum flow rate, ${\varphi }_i^{max}$, and a maximum node capacity, $n_i^{max}$, exist (i.e., ${\mathrm{\Phi }}_{ij}=0$ if $n_i\ge n_i^{max}$). A possible definition for the flows ${\mathrm{\Phi }}_{ij}(n_i,n_j)$ is

$${\mathrm{\Phi }}_{ij}(n_i,n_j)={\pi }_{ij}{\varphi }_j^{max}c(n_i/n_i^{max})\varphi \left(n_j\right),$$

(3)

where the function $\varphi \left(n_j\right)\in [0,1]$ is assumed to be monotonic, increasing with an initial linear dependence and an asymptotic limit, ${lim}_{n\to \infty }\varphi \left(n\right)=1$ (i.e., we do not consider a reduction of the outgoing flow when the node is congested, assuming that the dynamics at intersections is weakly affected by the road congestion), and the capacity function $c(n_i/n_i^{max})\in [0,1]$ is a threshold function that drops down to zero when $n_i\ge n_i^{max}$. The flow and the capacity functions simulate the effect of particle interactions that affect the traffic dynamics. In the case of urban road network, the flow function, $\varphi \left(n\right)$, is simulated by the traffic dynamics at intersections. As a consequence of definition (3), the transport network has a flow-density fundamental diagram when we increase the average traffic load for each node (cfr. Section 5). The equilibrium solution (2) follows from the condition

$$\sum _j{\pi }_{ij}{\varphi }_j^{max}c(n_i/n_i^{max})\varphi \left(n_j\right)=\sum _j{\pi }_{ji}{\varphi }_i^{max}c(n_j/n_j^{max})\varphi \left(n_i\right).$$

(4)

If there exists a solution with $n_i\le n_i^{max}$, the traffic load is sustainable, otherwise we have congested states with some nodes at maximum capacity. In the case of a low traffic load, we expect $\varphi \left(n_j\right)=\alpha n_j$ and $c(n_i/n_i^{max})=1$ and Equation (4) simplifies

$$\sum _j{\pi }_{ij}{\varphi }_j^{max}{n}_j^{\ast }={\varphi }_i^{max}{n}_i^{\ast }.$$

$n_i^{\ast }$ turns out to be the stationary eigenvector of the stochastic matrix ${\pi }_{ij}$, and the equilibrium state is stable since all the other eigenvalues, $\lambda$, of the Laplacian matrix of the network,

$$L_{ij}={\delta }_{ij}{\varphi }_i^{max}-{\pi }_{ij}{\varphi }_j^{max},$$

(5)

have a negative real part if the network is connected [35]. When the traffic load increases, the flow function $\varphi \left(n\right)\to 1$, and we have to consider a self-consistent approach assuming $c(n_i/n_i^{max})=1$ (i.e., no congestion in the network). We compute the null eigenvector ${\varphi }_i^{\ast }$,

$$\sum _j{\pi }_{ij}{\varphi }_j^{max}{\varphi }_j^{\ast }-{\varphi }_i^{max}{\varphi }_i^{\ast }=0,$$

and we look for the solutions $n_i^{\ast }$ of the system,

$$\varphi \left(n_i^{\ast }\right)\propto {\varphi }_i^{\ast }\sum _i{n}_i^{\ast }=N.$$

We study the stability of the equilibrium solution in the presence of perturbations, $\delta n_j$, such that ${\sum }_j\delta n_j=0$ (i.e., the total traffic load is constant). The stability character depends on the derivative ${\varphi }^{\prime }=d\varphi /dn$: if ${\varphi }^{\prime }\left(n_j^{\ast }\right)>0$$\forall j$, then the linearized system,

$$\delta {\dot{n}}_i=-\sum _j{L}_{ij}{\varphi }^{\prime }\left(n_j^{\ast }\right)\delta n_j,$$

is still associated with a Laplacian matrix, and the eigenvalues all have a negative real part on the invariant subspace, ${\sum }_j\delta n_j=0$. When ${\varphi }^{\prime }\left(n_j^{\ast }\right)\to 0$, the system tends to a neutral stability. However, if the equilibrium traffic load, $n_i^{\ast }$, increases up tp $n_i^{max}$, congested states appear due to the the congested function $c(n_i/n_i^{max})$ in definition (3). Assuming $n_i/n_i^{max}\simeq 1$, we can simplify Equation (4), letting $\varphi \left(n_j^{\ast }\right)\simeq {\varphi }_j^{max}$ (i.e., before the congestion the roads express their maximum flow rate), the equilibrium condition reads

$$\sum _j{\pi }_{ij}{\varphi }_j^{max}c(n_i/n_i^{max})=\sum _j{\pi }_{ji}{\varphi }_i^{max}c(n_j/n_j^{max}).$$

Therefore, the solution has the form $c(n_i^{\ast }/n_i^{max})\propto c_i$, where $c_i$ is the null eigenvector of the Laplacian matrix,

$$L_{ij}^{\prime }=\sum _k{\pi }_{ik}{\varphi }_k^{max}{\delta }_{ij}-{\pi }_{ji}{\varphi }_i^{max},$$

(6)

(cfr. Equation (5)) with the constraints

$$\sum _i{n}_i^{\ast }=N\mathrm{and}{n}_i^{\ast }<n_i^{max}.$$

By increasing the traffic load N, $n_i^{\ast }\to n_i^{max}$ for the nodes corresponding to the smallest values of $c_i$, and the congestion will start from these nodes. To study the stability of the solution, we consider the linearized system, i.e.,

$$\delta {\dot{n}}_i=\sum _j{L}_{ij}^{\prime }{c}^{\prime }\left(n_j^{\ast }\right)\delta n_j{c}^{\prime }\left(n_j^{\ast }\right)={\displaystyle \frac{dc}{d{n}_j}}\left(n_j^{\ast }\right),$$

and when $c^{\prime }\left(n_j^{\ast }\right)<0$, we have a stable congested state. We observe that if ${\mathrm{\Phi }}_j^{max}$ is the stationary eigenvector of the stochastic matrix ${\pi }_{ij}$, we have $L_{ij}=L_{ij}^{\prime }$, and the vectors ${\varphi }_i^{\ast }$ and $c_i$ are constant (i.e., at the equilibrium states, ${\mathrm{\Phi }}_i/{\varphi }_i^{max}$ and $n_i/n_i^{max}$ have the same value for all the nodes).

A possible interpretation is that individuals use the transport network in such a way that the traffic load is distributed, making the nodes equivalent from a dynamical point of view. This condition could reflect the emergence of a Wardrop equilibrium, where the paths distribution evolves to define a transition matrix, ${\pi }_{ij}$, which makes all nodes equivalent in terms of congestion formation. In other words, any change in the path distribution would make a node more susceptible to congestion.

Since the nodes are equivalent, congestion formation depends on the presence of traffic fluctuations that prevent the transport system from reaching its maximum flow rate, making congestion a dynamical stationary state. When a node i becomes congested ($n_i\ge n_i^{max}$), its incoming flow is null, rendering some links, ${\pi }_{ij}$, ineffective. Consequently, the connected nodes, j, have a higher probability of becoming congested, whereas node i can exit the congested state since its flow ${\mathrm{\Phi }}_i>0$. Thus, we do not have an equilibrium state with fixed congested nodes, but rather a stationary dynamical state where congestion moves across the network in regions of almost congested nodes.

The rise of congestion due to traffic load fluctuations is not described by the average dynamics (1), and the fluctuations define the traffic load’s stationary distribution on the transport network.

3. Random Walk on Graphs as Models of Transport Networks

To study the physics of traffic fluctuations, we propose simple models based on random walks on graphs [35]. According to the previous section, we associate a graph with a transport network where the transition probabilities, ${\pi }_{ij}$, define the link weights. This allows us to define a Markov process in which the ’particles’ move randomly according to the transition probabilities. In this way, it is possible to study the fluctuation statistics at stationary states and congestion formation.

To perform an analytical approach, we simplify model (3) by defining the flow function $\varphi \left(n\right)=\mathrm{\Theta }\left(n\right)$, where $\mathrm{\Theta }\left(n\right)$ is the Heaviside function (i.e., $\mathrm{\Theta }\left(n\right)=1$ if $n\ge 1$ and $\mathrm{\Theta }\left(0\right)=0$) and the capacity function $c\left(n\right)=\mathrm{\Theta }(n^{max}-n)$, where all the nodes have the same maximum capacity, $n_{max}$. We refer to the vector $\overrightarrow{n}$, whose components $n_i\ge 0$$i=1,...,M$ give the number of particles at node i as the dynamical state of the network, and we define $|\overrightarrow{n}|={\sum }_i{n}_i=N$ as the traffic load of the transport network. In the one-step process approximation, the evolution of the distribution function $\rho (\overrightarrow{n},t)$ is given by the master equation

$${\displaystyle \frac{1}{M}}\dot{\rho }(\overrightarrow{n},t)=\sum _{(i,j)}\left[E_i^{-}{E}_j^{+}{\mathrm{\Phi }}_{ij}(n_i,n_j)-{\mathrm{\Phi }}_{ji}(n_j,n_i)\right]\rho (\overrightarrow{n},t),$$

(7)

where we set the transition rates,

$${\mathrm{\Phi }}_{ij}(n_i,n_j)=\mathrm{\Theta }\left(n_j\right)\mathrm{\Theta }(n^{max}-n_i){\pi }_{ij},$$

(8)

and the sum runs over all the possible ordered couples $(i,j)$.

We remark that the one-step process implies an instantaneous update of information after each movement. This is a non-physical assumption for a transport network, and we also consider synchronous dynamics where all nodes evolve simultaneously, and the information of the network state is updated after the movement of all nodes. In this case, the master equation becomes more complicated and an analytical approach is not feasible. In Appendix A, we show how to write the incoming and outgoing flows for a generic node (see Equation (A5)).

We observe that if $j\ne i$, we have the identity

$$E_i^{-}{E}_j^{+}{\pi }_{ij}\mathrm{\Theta }\left(n_j\right)\mathrm{\Theta }(n^{max}-n_i)\rho \left(\overrightarrow{n}\right)=\mathrm{\Theta }\left(n_i\right)\mathrm{\Theta }(n^{max}-n_j){\pi }_{ij}\rho (\overrightarrow{n}+{\widehat{e}}_j-{\widehat{e}}_i),$$

where we set $\rho \left(\overrightarrow{n}\right)=0$ for any non-physical state and we write the master equation in the form

$${\displaystyle \frac{1}{M}}\dot{\rho }(\overrightarrow{n},t)=\sum _{(i,j)}\mathrm{\Theta }\left(n_i\right)\mathrm{\Theta }(n^{max}-n_j)\left[{\pi }_{ij}\rho (\overrightarrow{n}+{\widehat{e}}_j-{\widehat{e}}_i,t)-{\pi }_{ji}\rho (\overrightarrow{n},t)\right].$$

(9)

3.1. Equilibrium State in the Case of Detailed Balance

We are interested in the stationary distribution, ${\rho }_s\left(\overrightarrow{n}\right)$, of the master Equation (9). Let $\overrightarrow{p}$ be the null right eigenvector of the Laplacian matrix of the weighted network (i.e., $p_i$ is the probability to observe a particle in node i according to the transition rates (8)),

$$\sum _j{\pi }_{ji}{p}_i-{\pi }_{ij}{p}_j=0\sum _j{p}_j=1p_j>0.$$

(10)

We look for a stationary distribution of the form

$${\rho }_s\left(\overrightarrow{n}\right)=\left\{\begin{array}{cc}{\left[C_N^{n^{max}}\right]}^{-1}\left(\overrightarrow{p}\right)\prod _{k=1}^M{p}_k^{n_k}\hfill & |\overrightarrow{n}|=N{n}^{max}\ge n_i\ge 0\hfill \\ 0\hfill & \mathrm{otherwise}.\hfill \end{array}\right.$$

(11)

where $C_N^{n^{\ast }}\left(\overrightarrow{p}\right)$ is the normalizing constant (i.e., the partition function), and it is related to the Helmholtz Free Energy,

$$F\left(\overrightarrow{p}\right)=-ln{C}_N^{n^{max}}\left(\overrightarrow{p}\right)=-ln\sum _{|\overrightarrow{n}|=N}^{n^{max}}\prod _{i=1}^M{p}_i^{n_i},$$

(12)

where the sum runs over the physical states $n_i\le n^{max}$ (see Equation (A1) in Appendix A for some analytical estimates in the case of $n^{max}\gg 1$). The stationary condition reads

$$\sum _{(i,j)}\mathrm{\Theta }\left(n_i\right)\mathrm{\Theta }(n^{\ast }-n_j)\left[{\pi }_{ij}{\displaystyle \frac{p_j}{p_i}}-{\pi }_{ji}\right]{\rho }_s\left(\overrightarrow{n}\right)=0,$$

(13)

which has to be satisfied for all the physical states, $\overrightarrow{n}$. The main difficulty of obtaining an analytical solution of Equation (13) is the presence of empty and congested nodes at the same time. However, if the micro-dynamics satisfies a detailed balance (DB) condition,

$${\pi }_{ij}{p}_j={\pi }_{ji}{p}_i,$$

(14)

the distribution (11) is the stationary solution of the master Equation (9). The DB condition means that, in the stationary state, the probability to observe a path between the nodes $j\to i$ is the same if we consider the reverse displacement $i\to j$. We remark that the stationary distribution, $\overrightarrow{p}$, depends on the transition matrix, ${\pi }_{ij}$, but different transition matrices have the same distribution. The DB condition associates uniquely the transition matrix to the distribution $\overrightarrow{p}$. If we introduce the Gibbs Entropy,

$$\mathcal{S}\left[\rho \left(\overrightarrow{n}\right)\right]=-\sum _{|\overrightarrow{n}|=N}^{n^{max}}\rho \left(\overrightarrow{n}\right)ln\rho \left(\overrightarrow{n}\right),$$

(15)

the DB condition allows one to apply a maximum entropy principle (MPE) to compute the stationary distribution. Let ${\overline{n}}_i$ be the mean load per node,

$${\overline{n}}_i=\sum _{|\overrightarrow{n}|=N}^{n^{max}}{n}_i{\rho }_s\left(\overrightarrow{n}\right)i=1,...,M.$$

(16)

The Gibbs entropy (15) is maximum for the distribution (11) when the probabilities, $p_i$, satisfy the constraints (16). The thermodynamic approach allows one to characterize the statistical properties of the stationary state without considering the dynamics (7). An explicit computation of the entropy for a distribution of the form (11) gives

$$S(\overrightarrow{p},N)=ln{C}_N^{n^{max}}\left(\overrightarrow{p}\right)-\sum _{|\overrightarrow{n}|=N}^{n^{max}}{\rho }_s\left(\overrightarrow{n}\right)\sum _i{n}_{i}ln{p}_i=-F\left(\overrightarrow{p}\right)-\sum _i{\overline{n}}_{i}ln{p}_i,$$

and we get

$${\displaystyle \frac{\mathsf{\partial }S}{\mathsf{\partial }{p}_i}}=-{\displaystyle \frac{\mathsf{\partial }}{\mathsf{\partial }{p}_i}}F\left(\overrightarrow{p}\right)-{\displaystyle \frac{{\overline{n}}_i}{p_i}}.$$

(17)

The extremality condition reads

$${\overline{n}}_i=-p_i{\displaystyle \frac{\mathsf{\partial }}{\mathsf{\partial }{p}_i}}F\left(\overrightarrow{p}\right),$$

and it is verified by the distribution (11). When one increases the total traffic load, N, we observe that the distribution (11) is peaked on the marginal states where the nodes with the greater $p_i$ are congested, $n_i=n^{max}$, and the other nodes are empty. This means that the nodes with greater $p_i$ are hot spots for the congestion spreads and the mitigation policies have to reduce the traffic load on these nodes by redistributing the mobility paths. The situation is different when one considers a homogeneous transport network where $p_i=M^{-1}$ for all the nodes. Among all the possible distributions of the form (11), this provides the maximum value for the Gibbs entropy and the congestion formation is triggered by the traffic fluctuations since ${\overline{n}}_i=N/M$ for any node.

To study the relation between entropy and traffic load fluctuations, we compute the Hessian matrix of the entropy with respect to the variation of $p_i$,

$$\begin{array}{cc}\hfill & p_j{p}_i{\displaystyle \frac{{\mathsf{\partial }}^{2}S}{\mathsf{\partial }{p}_j\mathsf{\partial }{p}_i}}=p_j{p}_i{\displaystyle \frac{\mathsf{\partial }}{\mathsf{\partial }{p}_j}}{p}_i^{-1}\left(-p_i{\displaystyle \frac{\mathsf{\partial }}{\mathsf{\partial }{p}_i}}F\left(\overrightarrow{p}\right)-{\overline{n}}_i\right)=\hfill \\ \hfill & =-p_j{\displaystyle \frac{\mathsf{\partial }}{\mathsf{\partial }{p}_j}}{\left[C_N^{n^{\ast }}\right]}^{-1}\sum _{|\overrightarrow{n}|=N}^{n^{max}}{n}_i\prod _{k=1}^M{p}_k^{n_k}=-\left(〈n_i{n}_j〉-{\overline{n}}_j{\overline{n}}_i\right),\hfill \end{array}$$

(18)

where we have used Equation (17). When the system is in a stationary state, the covariance matrix of the traffic load, $n_i$, is the sensitivity of the entropy to the changes of the probability distribution, $p_i$. Since the distribution, $p_i$, depends on the transition probabilities, ${\pi }_{ij}$, the previous equation also gives a measure of the sensitivity to the perturbations of the transition matrix. Using the entropy as a measure of the disorder in the system, the result (18) means that when the fluctuations of the traffic load, $n_i\left(t\right)$, are large, the system may pass through ordered and disordered states during the evolution. We also have the reciprocity relations that holds in the DB condition,

$$p_j{\displaystyle \frac{\mathsf{\partial }{\overline{n}}_i}{\mathsf{\partial }{p}_j}}=p_i{\displaystyle \frac{\mathsf{\partial }{\overline{n}}_j}{\mathsf{\partial }{p}_i}}=-p_i{p}_j{\displaystyle \frac{\mathsf{\partial }F}{\mathsf{\partial }{p}_i\mathsf{\partial }{p}_j}}i\ne j,$$

to define how the average load of different nodes is affected by the change of the stationary probabilities; these relations correspond to the Onsager reciprocal relations.

In the case of a homogeneous network, the fluctuation variance at each node is a key indicator for understanding the rise of congestion. One expects that the fluctuations reduce when the traffic load is low ($n_i\ll n^{max}$) or highly congested ($n_i\simeq n^{max}$), so that one get a critical value for the total traffic load when the average fluctuation variance is maximum. Using the covariance matrix with $N{p}_i=n_i$,

$${\sigma }_{ij}\left(N\right)=-n_i{n}_j{\displaystyle \frac{{\mathsf{\partial }}^{2}S}{\mathsf{\partial }{n}_j\mathsf{\partial }{n}_i}}\left(N\right),$$

we introduce the following definition,

$$\overline{\sigma }\left(N\right)={\displaystyle \frac{1}{M}}\mathrm{Tr}\sigma ={\displaystyle \frac{1}{M}}\sum _i{n}_i^2{\displaystyle \frac{{\mathsf{\partial }}^{2}S}{\mathsf{\partial }{n}_i^2}}\left(N\right),$$

and the critical load $N_c$ satisfies

$${\displaystyle \frac{\mathsf{\partial }\overline{\sigma }}{\mathsf{\partial }N}}\left(N_c\right)=0.$$

(19)

Using numerical simulations on a simple traffic network model, we have computed the traffic load variance for the traffic load of single nodes as a function of the average traffic load (see Section 5). The simulations show that the maximum value of the standard deviation is achieved at the critical value of the flow-density fundamental diagram. The critical load, $N_c$, corresponds both to the maximum flow in the transport network and to the maximum uncertainty in traffic load distribution. When $N\ge N_c$, the congested nodes start to merge in clusters until a macroscopic large cluster emerges in the network. The congestion degree can be related to the dimension of congested clusters, which is the fingerprint of a percolation phase transition for the congestion formation [11].

3.2. Non-Equilibrium Stationary States

When the detailed balance condition (14) does not hold, the Markov process to model the vehicle dynamics realizes a non-reversible random walk, i.e., the statistics of the reverse paths on the transport network is not equivalent to the statistics of the original paths in the stationary states. As a consequence, we have the presence of probability stationary currents on the links $j\to i$, defined by

$$J_{ij}\left(\overrightarrow{n}\right)=\mathrm{\Theta }\left(n_i\right)\mathrm{\Theta }(n^{max}-n_j)\left[{\pi }_{ij}{\rho }_s(\overrightarrow{n}+{\widehat{e}}_j-{\widehat{e}}_i,t)-{\pi }_{ji}{\rho }_s(\overrightarrow{n},t)\right],$$

(20)

that correspond to net traffic flows moving on the loops of the transport network. In such a case, the MEP cannot be applied [28], but it is possible to maximize a local entropy or apply an entropy production principle [29]. If we consider the ensemble of states $\overrightarrow{n}$ with $0\le n_j<n^{max}$ so that $\mathrm{\Theta }(n^{max}-n_j)=1$, then the stationary solution satisfies

$$\sum _j\left[{\pi }_{ij}{\rho }_s(\overrightarrow{n}+{\widehat{e}}_j-{\widehat{e}}_i,t)-{\pi }_{ji}{\rho }_s(\overrightarrow{n},t)\right]=0,$$

and we recover a solution of the form (11) that maximizes the entropy if restricted to this ensemble (cfr. Equation (13)). Conversely, if one considers the ensemble of states $\overrightarrow{n}$ with $0<n_i\le n^{max}$ and $\mathrm{\Theta }\left(n_i\right)=1$ and we get the condition

$$\sum _i\left[{\pi }_{ij}{\rho }_s(\overrightarrow{n}+{\widehat{e}}_j-{\widehat{e}}_i,t)-{\pi }_{ji}{\rho }_s(\overrightarrow{n},t)\right]=0,$$

then we have a solution of the form

$${\rho }_s\left(\overrightarrow{n}\right)\propto \prod _{k=1}^M{q}_k^{-n_k},$$

(21)

where $\overrightarrow{q}$ is the stationary eigenvector of the Laplacian matrix associated with the reverse transition matrix, ${\pi }_{ij}^T$. The solution (21) approximates the stationary distribution in the case of high traffic loads, and it is a maximum entropy solution for the congested states. The distribution (21) is the distribution for the gap dynamics on a congested transport network. In the case of DB, $q_i=p_i^{-1}$ and the two distributions coincide. Changing the traffic load, the stationary distribution, ${\rho }_s\left(\overrightarrow{n}\right)$, changes its form, interpolating the two limit distributions. The congestion transition depends on the probability of intermediate states, where one has both empty and congested nodes.

4. Results: Congestion Formation in Balanced Transport Networks

To understand the role of traffic fluctuations in the congestion transition, we consider a homogeneous transport network where all the roads are equivalent, $p_i=const.$ Using the simple model (8), this condition corresponds to the balance condition

$$\sum _i{\pi }_{ij}=\sum _j{\pi }_{ij}.$$

(22)

For each node, the expected incoming flow is balanced by the expected outgoing flow. This condition is consistent with the existence of Wardrop’s equilibrium in transport systems. We also note that condition (22) does not imply detailed balance, which requires the symmetry of the transition rates, ${\pi }_{ij}$. The stationary solution can be approximated by a uniform distribution (cf. Equation (11)) for the majority of the microstates, but not for the network states $\overrightarrow{n}$, where empty and congested connected nodes are simultaneously present. Due to traffic fluctuations, these states are more probable when the average traffic load is $N/M=n^{max}/2$, where the probability of empty and congested nodes is equal, so they can play a significant role in the congestion transition.

To study how the node states are distributed in the network, we consider the single-node distribution $p\left(n\right)$, defined as the probability to observe a node in the state n in the stationary state. In an homogeneous transport network, all the nodes are equivalent so that the marginal distribution for a single node of the stationary distribution, ${\rho }_s\left(\overrightarrow{n}\right)$, approximates the traffic load distribution on the network in the thermodynamic limit. The distribution $p\left(n\right)$ can be related to the congestion degree by computing the probability of the almost congested states,

$$C(\tau ,N)\simeq \sum _{n\ge n^{max}-n\left(\tau \right)}p\left(n\right),$$

where a suitable $n\left(\tau \right)$ is a decreasing function of $\tau$.

In the case of low traffic load, $N/M\ll n^{max}$, it is possible to prove that in the thermodynamics limit the single-node distribution, $p\left(n\right)$, is approximated by an exponential distribution,

$$p\left(n\right)\simeq {\displaystyle \frac{1}{\overline{n}}}exp\left(-{\displaystyle \frac{n}{\overline{n}}}\right)\overline{n}={\displaystyle \frac{N}{M}},$$

(23)

and the congestion probability is exponentially small. Analogously, in the case of high traffic load (i.e., $N\simeq n^{max}M$), it is approximated by

$$p\left(n\right)\propto exp\left(-{\displaystyle \frac{n^{max}-n}{n^{max}-\overline{n}}}\right)n\le n^{max}$$

(24)

(see Appendix A for the proof), which explains why the number of congested nodes, $C(\tau ,N)$, increases with the traffic load, N. The true distribution interpolates the two approximations as the traffic load varies, but this is a continuous process and the distribution does not cross any singularity (see Section 5).

In the case of synchronous dynamics, each node moves one particle but can receive a variable number of particles depending on its connectivity degree. An analytic expression for the stationary distribution ${\rho }_s\left(\overrightarrow{n}\right)$ is not available, and we note that synchronous dynamics introduces correlations among the states of connected nodes. Indeed, a node connected to empty nodes has a reduced incoming flow on average, causing its population to decrease. Consequently, network states where empty nodes cluster are favored by synchronized dynamics. Similarly, if a node is connected to a congested node, its outgoing flows are reduced on average, causing its population to increase, and thus congested nodes also tend to cluster. This phenomenon can be interpreted by introducing entropic forces that attract empty and congested nodes in synchronous dynamics. In the last part of Appendix A, we discuss the relationship between synchronous dynamics and the appearance of entropic forces. In Section 5, we study the effect of entropic forces using numerical simulations.

Single-Node Dynamics

We now consider the following problem:

Is it possible to introduce an effective model for the node dynamics to quantify the effect of entropic forces?

The node state evolution can be modeled by an effective master equation that depends on the neighbor node states. In the one-step process model, it is possible to look for a self-consistent master equation for the distribution $p\left(n\right)$. One can consider the balance between the incoming and outgoing flows of a representative node i given its state $n\le n^{max}$, averaged over all the states of the network. The incoming flow is the transition probability from the state $n-1\to n$ so that $n\ge 1$ and we get

$$〈\sum _j{\pi }_{ij}\mathrm{\Theta }\left(n_j\right)-\sum _j{\pi }_{ji}\mathrm{\Theta }(n^{max}-n_j)〉=0.$$

One applies the following estimates

$$\langle \mathrm{\Theta }\left(n_j\right)\rangle =1-\pi \left(0\right|n-1)p(n-1)\langle \mathrm{\Theta }(n^{max}-n_j)\rangle =1-\pi \left(n^{max}\right|n)p\left(n\right),$$

where $\pi \left(0\right|n-1)$ is the conditional probability that a neighbor is empty when the node state is $n_i=n-1$ and $\pi \left(n^{max}\right|n)$ is the conditional probability that a neighbor is congested when the node state $n_i=n$. Using condition (22), we have the equilibrium

$$(1-\pi \left(0\right|n-1))p(n-1)-(1-\pi \left(n^{max}\right|n))p\left(n\right)=00<n\le n^{max}.$$

(25)

The conditional probabilities may depend on the degree of the node, i.e., the average is the result of an averaging on all the nodes. Equation (25) can be solved recursively using $p\left(0\right)$ to normalize the distribution. We remark that the conditional probabilities are estimated from the global dynamical properties of the transport network and measure the correlation between the states of connected nodes. In the one-step process case, the distribution is derived from an MEP, and we expect that the node states are independent so that in a thermodynamic limit $\pi \left(0\right|n)=p(0)$ and $\pi \left(n^{max}\right|n)=p\left(n^{max}\right)$. In synchronous dynamics, the stationary distribution does not maximize the entropy, and we expect that $\pi \left(0\right|n)$ is a decreasing function on n since the empty nodes tend to cluster and $\pi \left(n^{max}\right|n)$ is an increasing function. These effects are a consequence of the stochastic dynamics when the network state contains empty and congested nodes at the same time and it creates entropic forces that explains the correlation among the connected nodes.

5. Results: Numerical Simulations of the Transport Network Model

We have checked the applicability of the analytical estimate by simulating the simplified model (8) using a random network of 500 nodes with average degree $d=3$ (but the minimum degree is $d=2$). For a given transition matrix, ${\pi }_{ij}$, we have computed the stationary solution, $p_i$, and we have define the stationary flows ${\mathrm{\Phi }}_{ij}={\pi }_{ij}{p}_j$, assuming ${\varphi }^{max}=1$; we have the balance condition

$$\sum _j{\mathrm{\Phi }}_{ij}=\sum _j{\pi }_{ij}{p}_j=p_i=\sum _j{\pi }_{ji}{p}_i=\sum _i{\mathrm{\Phi }}_{ij}.$$

As previously discussed, this condition means that the flows distribute on the transport network so that the average incoming flows equal the average outgoing flows for each node, simulating an optimal use of the network. The maximum average flow for the whole network is defined by

$${\overline{\mathrm{\Phi }}}^{max}={\displaystyle \frac{1}{M}}\sum _{ij}{\mathrm{\Phi }}_{ij},$$

(26)

and the maximum node capacity is fixed at $n^{max}=10$. In the sequel, we refer to this model as the transport network model. The balance condition (22) for the flows highlights the role of fluctuations in the emergence of the congested states, but we remark that the transition matrix, ${\pi }_{ij}$, does not satisfy the DB condition (i.e., ${\mathrm{\Phi }}_{ij}$ is not symmetric). This model is used in all the simulations presented in the sequel.

In Figure 1, we plot the average flow nodes (26) with respect to the average traffic load per node $\overline{n}$ to get the fundamental diagram for the network [14] that points out that the maximum flow is reached at $\overline{n}=n^{max}/2$, which is lower than the maximum theoretical value since the presence of the congested nodes and empty nodes due to the fluctuations for any traffic load prevents the network transport capacity from reaching its maximum value. The fundamental diagram shows that the average flow is almost constant for a large fraction of traffic load $3\le \overline{n}\le 7$, and it quickly reduces to zero at the limit values $\overline{n}=0$ and $\overline{n}=10$. We observe that the traffic load $\overline{n}=7$, at which the flow rate begins to drop down, coincides with the critical value at which the clusters of congested nodes start to merge (see Figure 8). This value can be considered a precursor of a macroscopic congestion formation.

To explain the fundamental diagram, we compute the traffic load distribution, $p\left(n\right)$, on the network (see Section 4 for the definition) as a function of the average traffic load. The simulation results are shown in Figure 2, where we consider both the one-step process (top pictures) and the synchronous dynamics (bottom pictures). To average out the finite size fluctuations, the empirical stationary distributions are computed by applying an ergodic principle and averaging over ${10}^5$ evolution steps. For the one-step process, we compare the simulation results with the equilibrium distribution for the single-node dynamics (25) and find perfect agreement between the two approaches. The empirical distributions for $\overline{n}=3,7$ are also well approximated by the exponential distributions (23) and (24), as expected by an entropy principle. This is also true in the more realistic case of synchronous dynamics, where we compare the stationary empirical distribution computed by averaging over ${10}^5$ iterations with the exponential approximation for $\overline{n}=3,7$.

The distribution $p\left(n\right)$ indicates that, for any traffic load, the probability of observing both empty and congested nodes reduces the traffic flow on the network. Over the critical threshold $\overline{n}=5$, the distribution is peaked at the congested nodes. However, the emergence of congestion in the network is a continuous process, and no singular behavior is observed in the distribution.

In the synchronous evolution, we observe that, for a low traffic load, the empty state is underrepresented compared to the exponential interpolation, and states with traffic load $n>n^{max}$ are present.

Both these effects are a consequence of the fact that the node degree satisfies $d\ge 2$; indeed, the incoming flows can change the node load by up to d particles in one time step, so the frequency of the empty state is depressed compared to the one-step process. At the same time, if we have a node with a traffic load $n\in [n^{max}-1,n^{max}-d+1]$, at the next time step, its state will be $n\in [n^{max}+d-1,n^{max}]$, resulting in overloaded nodes that cannot receive any particles until their state is again $n<n^{max}$. The frequency of the overloaded node is fast decaying due to the dynamical rules. Even in the synchronous case, the exponential interpolation is a good approximation of the true stationary distribution according to a local maximum entropy principle. To illustrate how the emergence of congestion affects the single-node distribution, we have considered the dependence of the standard deviation on the mean traffic load for single nodes. The results are reported in Figure 3 for the average loads $\overline{n}=3,4,5$. In the cases $\overline{n}=3,4$, we propose a linear interpolation that corresponds to an exponential-like distribution. The point distribution in the pictures reflects the heterogeneity of node behavior in the transport network, which is maximum at the low traffic load $\overline{n}=3$. This correlation is completely lost at the critical load $\overline{n}=5$, when the standard deviation is maximum but the heterogeneity of the nodes is minimum.

However, to understand the emergence of congestion, one must also consider the spatial distribution of congested nodes to detect the presence of congested clusters. As previously discussed, synchronous dynamics introduces a correlation between the states of connected nodes when they are empty or congested. This correlation gives rise to entropic forces that tend to cluster the empty and congested nodes. To highlight this effect, we compute the stationary distribution when $\overline{n}=n^{max}/2$. The results for the considered random network are shown in Figure 4, both for the single-step process and the synchronous dynamics. We observe that, in the first case, the distribution is almost flat, whereas, in the second case, the presence of entropic forces induces a bimodal distribution with peaks at the empty and congested nodes.

To quantify the entropic forces, we have computed the conditional probabilities $\pi \left(0\right|n)$ and $\pi \left(n^{max}\right|n)$ that denote the probability that a node with load n has an empty neighbor or a congested neighbor, respectively. In the case of the single-step process, these probabilities are

$$\pi \left(0\right|n)\simeq p\left(0\right)\pi \left(n^{max}\right|n)\simeq p\left(n^{max}\right)$$

since the nodes evolve independently in the thermodynamic limit. Conversely, in the synchronous case we expect that the correlation among the node states will increase the values $\pi \left(0\right|0)$ and $\pi \left(n^{max}\right|n^{max})$ near the boundary states $n=0$ and $n=n^{max}$. The numerical results are reported in Figure 5 and Figure 6, where the effect of the entropic force is clearly visible at low and high traffic load and an interpolation of the numerical results by a power law distribution is proposed.

$$\begin{array}{cc}\hfill \pi \left(0\right|n)& \propto {\displaystyle \frac{p\left(0\right)}{{(1+n)}^{\alpha }}}\hfill \\ \hfill \pi \left(n^{max}\right|n)& \propto {\displaystyle \frac{p(n^{max}}{{(n^{max}+1-n)}^{\alpha }}}\hfill \end{array}$$

(27)

In the numerical simulations, we have initially considered an average traffic load of $\overline{n}=n^{max}/2$ (Figure 5), where the interpolation (27) of the numerical results highlights the symmetry of the entropic forces when we consider the effect on the empty and the congested nodes. Then, we have simulated a low traffic load, ${\overline{n}}_l\ll n^{max}$, and a high traffic load, ${\overline{n}}_h\simeq n^{max}$, but with the symmetric condition ${\overline{n}}_h=n^{max}-{\overline{n}}_l$ (Figure 6). The interpolation (27) provides almost the same exponent in the two cases consistently with the equivalence of the gap dynamics and the particle dynamics in the case of a balanced network. The case $\overline{n}=n^{max}/2$ is critical since the exponent is maximum (so one has a faster decaying of the correlation), suggesting that the effect of the two entropic forces tends to balance when the traffic load of the nodes is far from the limit values. This effect is illustrated by Figure 7, where we compute the normalized conditional probabilities $\pi \left(0\right|n\left)\right|p\left(0\right)$ and $\pi \left(n^{max}\right|n)/p\left(n^{max}\right)$ for the critical traffic load $\overline{n}=5$. The results also point out the symmetrical behavior of the particle and gap dynamics (except near the boundary values $n=0$ and $n\ge n^{max}$).

To see how the entropic forces depend on the traffic load, we have computed the normalized conditional probabilities $\pi \left(0\right|n)/p(0)$ and $\pi \left(n^{max}\right|n)/p\left(n^{max}\right)$ for a traffic load $\overline{n}=3$ (see Figure 6 left) and $\overline{n}=7$ (see Figure 6 right), respectively. The results suggest that the maximum effect of the entropic forces is achieved at the critical load $\overline{n}=n^{max}/2$, whereas the connected node states tend to be weakly correlated at low or high traffic loads. This can be understood since the probability $p\left(0\right)$ increases for low traffic load, decreasing the effect of the entropic force, and similarly it happens to $p(n^{max}$ at high traffic load. This explains the quite good exponential interpolation of the distribution $p\left(n\right)$ for the synchronous case (see Figure 2) that is computed by maximizing the entropy (i.e., considering the nodes independent in the thermodynamic limit).

To illustrate this effect, we compare the normalized conditional probabilities $\pi \left(0\right|n)/p(0)$ and $\pi \left(n^{max}\right|n)/p\left(n^{max}\right)$ for traffic loads $\overline{n}=3,4,5$ and $\overline{n}=7,6,5$, respectively. The results are reported in Figure 7, where we observe the entropic forces have a maximum effect at $\overline{n}=5$ and the symmetry of the particles and gap dynamics for the homogeneous transport network.

According to the previous results, the traffic load $\overline{n}=n^{max}/2$ is a critical value for system fluctuations. However, the emergence of macroscopic congestion depends on the distribution of congested nodes in the transport network as the traffic load increases. The formation of a giant cluster is the fingerprint of the percolation transition, which has been proposed to describe congestion on an urban road network [11]. We have computed the number of congested clusters, the largest cluster size, and the second largest cluster size for the transport network model, considering different traffic loads, both in the case of the one-step process approximation and in the synchronous dynamics. The simulation results are shown in Figure 8, where the percolation transition occurs when the largest cluster size steepens, while both the number of clusters and the second largest cluster size decrease due to the formation of a large congested cluster. The results indicate a percolation threshold when the average traffic load is $7<\overline{n}<8$ for both the one-step process and the synchronous dynamics. The main difference is observed in the number of congested clusters, which is lower in the case of synchronous dynamics. This is consistent with the presence of entropic forces that tend to cluster the congested nodes (thus reducing their number) and result in a larger dimension of the largest cluster in the synchronous case. The coincidence of the percolation threshold in both cases can be explained by the weak effect of entropic forces at high traffic loads. The formation of a giant congested cluster can destroy the connectivity of the transport network (percolation transition), even if the cluster is not static in the model, as the outgoing flow from the congested nodes moves the cluster within the network.

6. Conclusions

We adopted a reductionist approach to highlight the universal properties of the congestion transition in an urban road network. Our approach uses a Markov random process on a graph to simulate two fundamental features of a transport network: the existence of a maximum flow rate and the existence of a maximum capacity. Two different dynamics are considered: the one-step process, where a randomly chosen node moves each time step, and the synchronous dynamics, where all nodes move together using the same information on the network state. Our main goal is to highlight universal characteristics of the congestion transition and, in particular, to study the effects of traffic load fluctuations.

Throughout our analysis, we maintained the assumption that the transport network is in a balanced condition (cf. def. (22)), where the average incoming flow at each node equals the average outgoing flow. This condition means that the average flows are in equilibrium, and this is related to the concept of Wardrop equilibrium, where paths are optimized across the network to avoid congestion. In the case of non-balanced transport networks, there exist nodes that inevitably become congested as the traffic load increases and play the role of hot spots for congestion formation. Conversely, if the balance condition is satisfied, congestion is driven by traffic load fluctuations. Our analytical methods, particularly the application of a maximum entropy principle, allowed us to characterize the distribution of fluctuations across nodes and identify how congestion emerges in a continuous way as the overall traffic load increases. Notably, we study the entropic forces in synchronous dynamics, leading to the clustering effect of congested and empty nodes.

The stationary distribution of traffic load fluctuations can be approximated by a local entropy principle, and the congestion transition is characterized by a peak in the traffic load variance of the nodes, yet without singularity in the thermodynamic limit. Numerical simulations further demonstrated the formation of a macroscopic congested cluster, resembling a percolation-like transition within the network, requiring traffic loads nearing maximum capacity for its emergence. The simulations highlight the effect of the entropic forces in the synchronous process on the congestion transition as a greater dimension of the congested clusters compared to the one-step process, but these slightly affect the critical traffic load for the percolation transition.

The considered model is certainly a simplification of the microscopic urban traffic complexity, but it allows us to underscore some universal macroscopic features. Our results explain the emergence of congestion as a continuous process that creates clusters of congested nodes mimicking the percolation transition (i.e., the formation of a macroscopic congested cluster at high traffic load). The existence of entropic forces depends on the network dynamics and favors the cluster formation of empty and congested nodes. This research not only provides insights into the nature of traffic congestion but also contributes to the theoretical framework necessary for enhancing urban transportation systems.

The application of control strategies to mitigate the congestion effect requires two further steps. On one hand, the time scales of congestion formation related to the relaxation time scales of the system need to be studied to estimate congestion diffusion in the transport network. On the other hand, strategies must be developed to reduce the fluctuation variance of the local traffic load across the whole network using intersection dynamics. A quantitative application of our results to real transport networks certainly requires more complex models, but we do not expect new physical phenomena. The possibility of performing a validation process needs datasets containing information on traffic load fluctuations at a road scale. This could be possible using GPS datasets with a representative sample of vehicle trajectories on homogeneous road networks, but at the moment they are not publicly available.