Influence of the number of predecessors in interaction within acceleration-based flow models

TL;DR

This paper analyzes how the number of predecessors in acceleration-based flow models affects stability, revealing non-monotonic relations and providing conditions for stability in microscopic pedestrian and car-following models.

Contribution

It introduces a stability analysis considering multiple predecessors and uncovers unexpected non-monotonic effects on stability in flow models.

Findings

Stability depends on small relaxation times.

Non-monotonic relation between number of predecessors and stability.

General conditions for linear stability on ring geometries.

Abstract

In this paper, the stability of the uniform solutions is analysed for microscopic flow models in interaction with $K\ge1$ predecessors. We calculate general conditions for the linear stability on the ring geometry and explore the results with particular pedestrian and car-following models based on relaxation processes. The uniform solutions are stable if the relaxation times are sufficiently small. The analysis is focused on the relevance of the number of predecessors in the dynamics. Unexpected non-monotonic relations between $K$ and the stability are presented.

Figures

Tables

Table 1. Table 1 : Summary of the influence of the number of neighbors in interaction K 𝐾 K on the stability for the exponential and algebraic, the optimal velocity and the generalized centrifugal force acceleration-based models.

Models	Range	Positive equil. speed	Stability for $K$	Cut-off radius	Unstable wave
Exponential	Short	$A\tau {\sum }_{k=1}^K{e}^{-kd/B}\le v_0$	Bounded,	Implicit	Longest
SFM [14]			non-monotonic
Algebraic	Long	$A\tau {\sum }_{k=1}^K{k}^{-q}\le u^q{v}_0$	Bounded,	Implicit	Characteristic
[15, 32]	if $q\le 1$		non-monotonic		for $K$ large
OV [20]	Long	Given by $V(\cdot )$	Unbounded,	To be	Characteristic
$q\le 1$			non-decreasing	calibrated	for $K$ large
$q>1$	Short	Given by $V(\cdot )$	Bounded,	Implicit	Characteristic
			non-decreasing		for $K$ large
GCF [16]	Long	Solution of $K+1$	Unbounded,	To be	Longest
$\delta =1$		degree polynomials	non-decreasing	calibrated
$\delta =0$	Long	Solution of $K+1$	Bounded,	Implicit	Characteristic
		degree polynomials	non-decreasing		for $K$ large

Equations

$$x_1(0)\le x_2(0)\le \dots \le x_N(0)\le L+x_1(0).$$

$$\left\{\begin{array}[]{lcl}\dot{x}_{n}(t)&=&v_{n}(t),\\[2.84526pt] \dot{v}_{n}(t)&=&A\left(v_{n}(t),x_{n+1}(t)-x_{n}(t),v_{n+1}(t),\ldots,x_{n+K}(t)-x_{n}(t),v_{n+K}(t)\right),\end{array}\right.$$

$$A\left(v,d,v,2d,v,\dots ,Kd,v,\right)=0.$$

$$\forall t\ge 0,\forall n=1,\dots ,N,x_{n+1}^H(t)-x_n^H(t)=d,x_n^H(t)=x_n^H(0)+vt,$$

$${\tilde{x}}_n(t)=x_n(t)-(x_n^H(0)+vt)\text{mbox}and{\tilde{v}}_n(t)=v_n(t)-v.$$

$$\textstyle\dot{\tilde{v}}_{n}(t)=\sum_{k=1}^{K}\alpha_{k}\big{(}\tilde{x}_{n+k}(t)-\tilde{x}_{n}(t)\big{)}+\sum_{k=0}^{K}\beta_{k}\tilde{v}_{n+k}(t)+R\big{(}\tilde{v}_{n}(t),\tilde{x}_{n+1}(t)-\tilde{x}_{n}(t),\,\ldots\big{)},$$

$${\alpha }_k=\frac{\partial A}{\partial d_k}\left(v_0,d_1,v_1,\dots ,d_k,v_k,\dots \right),{\beta }_k=\frac{\partial A}{\partial v_k}\left(v_0,d_1,v_1,\dots ,d_k,v_k,\dots \right),$$

$$\left\{\begin{array}[]{lcl}\dot{\bar{x}}_{n}(t)&=&{\bar{v}}_{n}(t),\\[2.84526pt] \dot{\bar{v}}_{n}(t)&=&\sum_{k=1}^{K}\alpha_{k}\big{(}\bar{x}_{n+k}(t)-\bar{x}_{n}(t)\big{)}+\sum_{k=0}^{K}\beta_{k}\bar{v}_{n+k}(t).\end{array}\right.$$

$${\lambda }^2=\sum _{k=1}^K{\alpha }_k(e^{ik\theta }-1)+\lambda \sum _{k=0}^K{\beta }_k{e}^{ik\theta }.$$

$$\sum _{k=0}^K{\beta }_k<0,$$

$$\left|\begin{array}[]{lcl}\mu_{l}&=&-\sum_{k=0}^{K}\beta_{k}\,c_{lk},\\[2.84526pt] \nu_{l}&=&\sum_{k=1}^{K}\alpha_{k}\left(1-c_{lk}\right),\end{array}\right.\qquad\quad\left|\begin{array}[]{lcl}\sigma_{l}&=&-\sum_{k=1}^{K}\beta_{k}\,s_{lk},\\[2.84526pt] \rho_{l}&=&-\sum_{k=1}^{K}\alpha_{k}\,s_{lk},\end{array}\right.$$

$${\mu }_l>0\text{mbox}and\det \left[\begin{array}{ccc}{\mu }_l& 0& -{\rho }_l\\ 1& {\nu }_l& -{\sigma }_l\\ 0& {\rho }_l& {\mu }_l\end{array}\right]={\mu }_l({\nu }_l{\mu }_l+{\rho }_l{\sigma }_l)-{\rho }_l^2>0,\forall l=1,\dots ,\lceil N/2\rceil .$$

$$\left\{\begin{array}[]{lcl}\dot{x}_{n}(t)&=&v_{n}(t),\\[2.84526pt] \dot{v}_{n}(t)&=&\frac{1}{\tau}(v_{0}-v_{n}(t))-\sum_{k=1}^{K}f\left(x_{n+k}(t)-x_{n}(t)\right),\end{array}\right.$$

$$v=v_0-\tau \sum _{k=1}^{K}f(kd).$$

$$\left|\begin{array}[]{lcl}\mu_{l}&=&1/\tau,\\[2.84526pt] \nu_{l}&=&-\sum_{k=1}^{K}f^{\prime}(kd)\left(1-c_{lk}\right),\end{array}\right.\qquad\left|\begin{array}[]{lcl}\sigma_{l}&=&0,\\[2.84526pt] \rho_{l}&=&\sum_{k=1}^{K}f^{\prime}(kd)s_{lk},\end{array}\right.$$

$$0<\tau <{\tau }_K=\underset{l=1,\dots ,\lceil N/2\rceil }{\min }{\tau }_K^{(l)},\text{mbox}with{\tau }_K^{(l)}=\frac{\sqrt{-{\sum }_{k=1}^K{f}^{\prime }(kd)\left(1-c_{lk}\right)}}{\left|{\sum }_{k=1}^K{f}^{\prime }(kd)s_{lk}\right|}.$$

$$\sum _{k=1}^{\infty }|f^{\prime }(dk)|<\infty \text{mbox}foranyd>0.$$

$$f(d)=A{e}^{-d/B}.$$

$$\tau \le \frac{v_0}{A{\sum }_{k=1}^K{e}^{-kd/B}}.$$

$$u=d/B\text{mbox}and{\tilde{\tau }}_K^{(l)}=\sqrt{\frac{A}{B}}{\tau }_K^{(l)}.$$

$${\tilde{\tau }}_K^{(l)}(u)=\frac{\sqrt{{\sum }_{k=1}^K{e}^{-ku}\left(1-c_{lk}\right)}}{\left|{\sum }_{k=1}^K{e}^{-ku}{s}_{lk}\right|}.$$

$$f(d)=\frac{A}{(d/B{)}^q}.$$

$$\tau \le \frac{u^q{v}_0}{A{\sum }_{k=1}^K{k}^{-q}}.$$

$${\tilde{\tau }}_K^{(l)}(u,q)=\frac{\sqrt{\frac{u^{q+1}}{q}{\sum }_{k=1}^K{k}^{-q-1}(1-c_{lk})}}{\left|{\sum }_{k=1}^K{k}^{-q-1}{s}_{lk}\right|}.$$

$$\phi =1-\frac{{\min }_K{\tilde{\tau }}_K}{{\max }_K{\tilde{\tau }}_K}=1-\frac{{\min }_K{\tau }_K}{{\max }_K{\tau }_K},$$

$$\frac{\partial {\tilde{\tau }}_K^{(l)}}{\partial u}(u)=\frac{{\tilde{\tau }}_K^{(l)}(u)}{2}{g}_K(u)>0$$

$$\frac{\partial {\tilde{\tau }}_K^{(l)}}{\partial u}(u,q)=\frac{q+1}{u}{\tilde{\tau }}_K^{(l)}(u,q)>0.$$

$$\frac{\partial {\tilde{\tau }}_K^{(l)}}{\partial q}(u,q)=\frac{{\tilde{\tau }}_K^{(l)}(u,q)}{2}\left(\ln u-\frac{1}{q}+h_K(q)\right).$$

$$\left\{\begin{array}[]{lcl}\dot{x}_{n}(t)&=&v_{n}(t),\\[2.84526pt] \dot{v}_{n}(t)&=&\sum_{k=1}^{K}a_{k}\left\{V\left(\frac{1}{k}(x_{n+k}(t)-x_{n}(t))\right)-v_{n}(t)\right\}.\end{array}\right.$$

$$0<\tau\cdot V^{\prime}(d)<\frac{\Big{(}\sum_{k=1}^{K}k^{-q}\Big{)}^{2}\sum_{k=1}^{K}\left(1-c_{lk}\right)k^{-(q+1)}}{\left(\sum_{k=1}^{K}s_{lk}k^{-(q+1)}\right)^{2}}=:\tilde{\tau}_{K}^{(l)}.$$

Full text

Keywords: Microscopic flow models, number of neighbors in interaction, homogeneous solution, linear stability analysis

1 Introduction

1.1 Context

Self-driven many-body systems are frequently used to model pedestrian or road traffic streams [1, 2, 3]. Continuous speed and acceleration models are defined with systems of ordinary, stochastic or delayed differential equations. They admit the uniform configurations (where all the spacing and speed are constant and equal) as a stationary solution. The stability analysis allows to determine conditions on the parameters for which perturbations around the uniform solution disappear (see for instance [4, 5, 6]). Real observations of pedestrian or vehicle streams present non-homogeneous configurations, with propagation of so-called stop-and-go waves for congested density levels [7, 8]. Therefore, realistic models should have unstable uniform solutions [9].

Traffic flow models can be microscopic, mesoscopic (kinetic), or macroscopic [6]. Microscopic models describe the individual trajectories in Lagrangian coordinates. They are usually called self-driven or motorized particle systems. Microscopic models are used since the 1950’s for traffic flow modelling in 1D (see for instance [10, 11, 12]), and since the 1970’s for pedestrians in 2D (see, e.g., [13, 14, 15, 16]). Speed-based models refer to as models of the first order while acceleration-based models generally refer to as second order models, or again force-based models for pedestrian dynamics in analogy to the Newtonian mechanic. The models are investigated theoretically or as simulation tools, notably to understand the formation and propagation of stop-and-go waves [5, 17].

The number of neighbors in interaction is a central parameter in microscopic models. It is related to as an anticipation factor in traffic modelling (spatial anticipation, see [18]). Many car-following models with several predecessors in interaction exist in the literature (see, e.g., [19, 20, 21, 22, 23, 24]). Most of the approaches are extensions of the optimal velocity (OV) model [4]. All of them show that increasing the number of predecessor in interaction stabilizes the uniform solutions. However, it is unclear whether stabilization occurs for models different from extended OV models. Such a feature is fundamental for connected and autonomous vehicles (CAV) [25, 26, 27]. For the pedestrian dynamics, the models are generally defined with infinite interaction ranges (i.e. all neighbors in interaction) while finite ones (cut-off radius) are used for the simulation to limit the computational effort. Nowadays, the impact of the cut-off radius restriction on the dynamics is mainly investigated by simulation (see for instance [28]).

Technically, the stability conditions are derived by linearising the models around the homogeneous solution, determining the characteristic equation of the linear system and bounding the real part of the characteristic equation solutions. In the traffic literature, stability analysis are broadly used with particular models or even variants of existing models [29, 19, 4, 21, 30, 24]. To our knowledge, there exists no systematic analysis of generic models, especially when the number of predecessors in interaction is taken as an arbitrary parameter. Furthermore, the stability analysis are generally investigated for the longest wavelength (see, e.g., [29, 19, 21, 24]). Such an approach uses a second linearisation process in the characteristic equation allowing to simplify the calculation. However the longest wavelength is not systematically the most unstable wave in the system. Some studies notably report characteristic wavelengths in stationary states in case of instability (see for instance [9]).

In this article, we calculate general linear stability conditions of an acceleration-based model with $K\geq 1$ predecessors in interaction for any wavelength. The results are explored with particular pedestrian and car-following models including some kind of velocity-relaxation process (e.g. force-based or optimal velocity models). The analysis is focused on the relevance of the number of predecessors $K$ for the interaction on the critical relaxation time. Unexpected non-monotonic relations between $K$ and the stability are pointed out. The organization of the article is the following. The model and its uniform solution are defined in the two next subsections. The general linear stability conditions are calculated in Sec. 2. Applications to distance and speed-based models are investigated in Secs. 3, 4 and 5. Secs. 6 and 7 present summary and conclusion.

1.2 Definition of the model

We consider single-file motion with $N$ agents on a one-dimensional ring of length $L$. We denote $n$ the index of the agents and \big{(}x_{1}(t),\ldots,x_{N}(t)\big{)} their curvilinear positions at time $t$ (see Fig. 1). The initial positions are such that

[TABLE]

The dynamics of the system are supposed strictly asymmetric (i.e. forward interaction). They are described by the second order (i.e. acceleration-based) model

[TABLE]

for all $t\geq 0$ and all $n=1,\ldots,N$. The acceleration $A$ of the agent $n$ at time $t\geq 0$ depends on the speed and distance spacing, and on the speeds and distance spacings of the $K$ predecessors at the same time. We assume that $n+k$ is $n+k-N$ if $n+k>N$, and that $x_{m}-x_{n}$ is $L+x_{m}-x_{n}$ when $m<n$. Furthermore we suppose that $0<K<\!\!<N$ and that the function $A$ is at least differentiable.

1.3 Uniform solution

For a given mean spacing $d=L/N$, we suppose that there exists an equilibrium speed $v$ such that

[TABLE]

Under this assumption, the uniform (or homogeneous) configurations $H$ such that

[TABLE]

are stationary solutions for the system (2). Here the uniform configurations $(d,v)$ depend on the initial conditions \big{(}x_{1}(0),\ldots,x_{N}(0)\big{)}. We investigate in the next section the linear stability of these solutions.

2 Linear stability analysis

We determine conditions for the parameters of the acceleration-based model (2) for which the uniform solutions (4) are linearly stable. The local stability analysis, for a finite line of agents with a leader travelling at a known speed, is distinguished from the global stability, for agents on a ring or on an infinite line [5]. The global stability conditions are more restrictive since they contain as well convective perturbations that can locally vanish [6]. Here the global stability conditions are calculated for the ring geometry.

2.1 Linearisation of the system

The stability conditions are calculated by considering the evolution of the perturbation

[TABLE]

A uniform solution $H$ is said to be stable if $\tilde{x}_{n}(t)\operatorname*{\rightarrow}0$ and $\tilde{v}_{n}(t)\operatorname*{\rightarrow}0$ as $t\operatorname*{\rightarrow}\infty$ for all $n=1,\ldots,N$. By construction $\dot{\tilde{x}}_{n}(t)=\tilde{v}_{n}(t)$ and, from the Taylor’s theorem and by using (3),

[TABLE]

with

[TABLE]

the partial derivatives of $A:(v,d,v,\ldots)\mapsto A(v,d,v,\ldots)$ and where $R(v,d,v,\ldots)$ is a function tending to zero as $(v,d,v,\ldots)$ tend to zero. Neglecting the rest $R$, we obtain the linear dynamics

[TABLE]

In the following, we focus on the evolution of the variables $(\bar{x}_{n}(t),\bar{v}_{n}(t))$. A uniform configuration $H$ is said to be linearly stable if $\bar{x}_{n}(t)\operatorname*{\rightarrow}0$ and $\bar{v}_{n}(t)\operatorname*{\rightarrow}0$ as $t\operatorname*{\rightarrow}\infty$ for all $n=1,\ldots,N$.

2.2 Calculation of the characteristic equation

The linear stability analysis requires the solution of a second order polynomial characteristic equation with complex coefficients (or even a $2N$ polynomial with real coefficients). The characteristic equation of the linear system (8) can be obtained by matrix calculus (see Appendix 1), or directly by using the exponential Ansatz $\bar{x}_{n}(t)=\xi\,e^{\lambda t+in\theta}$ with $\xi\in\mathbb{R}$ and $\lambda\in\mathbb{C}$. One gets in the last case

[TABLE]

Here $\theta\in[0,2\pi]$ describes the wavelength of the system when it is unstable. We have for the ring $\theta=2\pi l/N$ and $2N$ solutions $(\lambda_{l}^{(p)},p=1,2,l=1,\ldots,N)$ for (9). For $l=N$ (or $\theta=0$), the two solutions are $\lambda^{(1)}_{N}=0$ and $\lambda^{(2)}_{N}=\sum_{k=0}^{K}\beta_{k}$. Assuming

[TABLE]

the uniform solution (4) is linearly stable if the zeros $(\lambda_{l}^{(p)},p=1,2)$ of the characteristic equation (9) have strictly negative real parts for all $l=1,\ldots,N-1$.

2.3 General stability condition

Equation (9) is a second order polynomial of the form $\lambda^{2}+w_{l}\lambda+z_{l}=0$ with complex coefficients $w_{l}=\mu_{l}+i\sigma_{l}$, $z_{l}=\nu_{l}+i\rho_{l}$ such that

[TABLE]

with $c_{lk}=\cos\left(2\pi lk/N\right)$ and $s_{lk}=\sin\left(2\pi lk/N\right)$. The sufficient and necessary conditions for that a polynomial with complex coefficients have all its zeros in the half-plane $\Re(\lambda)<0$ are given in [31, Th. 3.2]. The results are a generalization of the Hurwitz conditions for polynomials with real coefficients. These conditions are

[TABLE]

The condition is general but difficult to interpret. It will be investigated with particular distance and speed-based pedestrian and car-following models in the three next sections. Note that the inequalities are the same if we substitute $l$ by $N-l$. Therefore the study is restricted to $l=1,\ldots,\lceil N/2\rceil$ instead of $l=1,\ldots,N-1$. The condition (12) is confirmed with some particular models for which the conditions are well-known in the literature, see Appendix 2.

3 Pedestrian models with distance-based repulsive force

Many pedestrian models are described by the sum of a relaxation to the desired speed and additive repulsions with the neighbors (see for instance [14, 15, 32, 33, 34])

[TABLE]

with $v_{0}>0$ the desired speed, $\tau>0$ the relaxation time and $f(\cdot)$ a $C^{1}$ positive and non-increasing function for the repulsion. Here, the repulsive force $f(\cdot)$ solely depends on the spacing. In this model class, the equilibrium speed $v$ solution of (3) is, for a given mean spacing $d$,

[TABLE]

The equilibrium speed is positive if $\tau\leq v_{0}/\sum_{k}f(kd)$ and converges as $K\rightarrow\infty$ if $\sum_{k}f(kd)<\infty$. This is a preliminary restriction on the relaxation time parameter that a model has to satisfy. Note that long-range models cannot by construction fulfil the positivity condition of the equilibrium speed as $K\operatorname*{\rightarrow}\infty$.

With the model (13), the partial derivatives are $\beta_{0}=-1/\tau$, $\beta_{k}=0$ for all $k>0$, and $\alpha_{k}=-f^{\prime}(kd)$ for all $k$. The parameters (11) are

[TABLE]

and do not depend on the equilibrium speed $v$. The first linear stability condition in (12) is $-1/\tau<0$. It is always true and implies the preliminary assumption (10). The second condition in (12) gives the stability threshold for the relaxation time

[TABLE]

Note that $f^{\prime}(d)\leq 0$ for all $d$. The stability occurs for a relaxation time $\tau$ sufficiently small. Since $f$ is $C^{1}$ positive and non-increasing, $\lim_{x}f(x)=f(y)-\int_{y}^{\infty}|f^{\prime}(u)|\,\mbox{d}u=0$ for all $y>0$. This implies $\int_{y}^{\infty}|f^{\prime}(u)|\;\mbox{d}u<\infty$ or again, using the Cauchy criterion and changing the variable,

[TABLE]

This proves the convergence of $\tau_{K}^{(l)}$ and $\tau_{K}$ within pedestrian models with distance-based repulsive force. Therefore the impact on the stability of the number of predecessors in interaction $K$ is limited, i.e. there exists an intrinsic interaction range, as soon as the repulsion $f(x)\operatorname*{\rightarrow}0$ as $x\operatorname*{\rightarrow}\infty$ (irrespective of the convergence speed). The stability condition at the limit $K\rightarrow\infty$ may be well approximated for a finite value of $K$. Therefore the use of cut-off radius in such models can have not effect on the dynamics. Yet, the value of the range depends on the convergence speed of the series $\sum_{k}|f^{\prime}(dk)|$. This point is investigated using exponential and algebraic repulsive forces.

Exponential force

We first consider the exponential repulsive force with parameters $A,B>0$

[TABLE]

This short-range force is used within the social force model [14]. The equilibrium speed is positive with this model as soon as

[TABLE]

In the following, we use the dimensionless spacing and critical relaxation time

[TABLE]

We have with the exponential repulsive force (18)

[TABLE]

Here $u$ is a form parameter and $A$ and $B$ are scale parameters for the dimensioned $\tau_{K}^{(l)}$, while $v_{0}$ has no influence.

Algebraically decaying force

We consider the algebraic repulsive force with parameters $A,B,q>0$

[TABLE]

This model is used in [15] with $q=1$ and in [32] with $q=2$. Here, the repulsive force can be short or long range according to the value of $q$. The model is long range for $q\leq 1$ in one dimension. The equilibrium speed converge as $K\rightarrow\infty$ if and only if $q>1$ (i.e. for short range forces). More precisely, the equilibrium speed is positive if

[TABLE]

The dimensionless critical relaxation time is

[TABLE]

Here $A$ and $B$ are scale parameter for $\tau_{K}^{(l)}$, $v_{0}$ has no influence while $u$ and $q$ are form parameters.

3.1 Stability condition function of $K$

The uniform solutions (4) are linearly stable for the distance-based models with forces (18) and (22) if the relaxation time $\tau$ is strictly less than critical time $\tau_{K}=\min_{l}\tau_{K}^{(l)}$ (i.e. by using dimensionless variables, if $\sqrt{A/B}\tau<\tilde{\tau}_{K}$). Thus we have to calculate the minimum of the functions $l\mapsto\tilde{\tau}_{K}{(l)}$ to determine the stability condition. Yet, the sign of the derivative of the function $l\mapsto\tilde{\tau}_{K}{(l)}$ is hard to extract. Therefore, we investigate it numerically with $N=10^{5}$.

The critical times (21) for the exponential model are plotted as a function of $l$ in Fig. 2, top panels. Here, $K$ varies from $1$ to $25$, and we compare $u=0.4$, $1$ and $2$. We observe that for the exponential force, the $\tilde{\tau}_{K}^{(l)}$ are minimal for $l=1$ for all $K$, i.e. $\tilde{\tau}_{K}=\tilde{\tau}_{K}^{(1)}$. Further numerical investigations (not shown here) confirm it. This means that the longest wavelength is always the most unstable with the exponential repulsive force. The critical times (24) for the algebraic model are plotted as a function of $l$ in Fig. 2, bottom panels, with $q=1$, $2$ and $3$. One observes that $\tau_{K}^{(l)}$ are minimal for $l=1$ when $K$ is low. For high $K$, the minimums are reached for $l_{q,K}\propto N$. Further results show that $l_{q,K}$ converges when $K$ increases. Therefore the most unstable wavelength has a characteristic size proportional to $N$ with the algebraic force, if $K$ is sufficiently high.

The critical relaxation times $\tilde{\tau}_{K}=\min_{l}\tilde{\tau}_{K}^{(l)}$ determine the border of the linear stability for any wavelength. They are plotted as function of $K$ in Fig. 3 for the exponential force (top panels), and for the algebraic one (bottom panels). The critical time converges to a constant value through a single damped oscillation for both models. This non-monotonic relation between $K$ and the stability is unexpected. Increasing the number of pedestrians in interaction first results in a decrease of the stability (at list until $K=2$). Then increasing $K$ increases $\tau_{K}$ and so the stability. The convergence of $\tau_{K}$ is relatively smooth for the exponential force. One observes a brusque transition for the algebraic force when the stability is broken for wavelength smaller than the longest one. In general, the speed of damping of the function $K\mapsto\tilde{\tau}_{K}$ depends on the spacing parameter $u$ for the exponential model, and on the force range parameter $q$ for the algebraic model.

3.2 Variation of the stability condition

The proportion of variation of the critical relaxation time

[TABLE]

allows to quantify the influence of the number of predecessors in interaction $K$ on the stability. It is the same for initial and dimensionless critical times $\tau_{K}$ and $\tilde{\tau}_{K}$. For $\varphi\approx 0$, the model weakly depends on the interaction range and inversely for $\varphi\approx 1$. The proportion $\varphi$ does not depend on $A$ and $B$ for both exponential and algebraic forces. It depends on $u$ for the exponential model, while it depends on $q$ for the algebraic model, but not on $u$. In the Fig. 4, the variation $\varphi$ is plotted as a function of $u$ for the exponential model (left panel), and as a function of $q$ for the algebraic model (right panel). $\varphi$ tends to zero as the dimensionless mean spacing $u$ increases with the first model. This means that $\tilde{\tau}_{K}$ marginally varies for low density levels. The same phenomenon occurs as $q$ increases with the algebraic model.

3.3 Stability conditions for other parameters

The critical relaxation time $\tilde{\tau}_{K}$ depends on the parameter $u$ for the exponential force, see (21), and on $(u,q)$ for the algebraic force, see (24). We investigate here how the stability condition depends on these parameters. For the exponential model, the function $u\mapsto\tilde{\tau}_{K}(u)$ increases as $u$ increases. This means that the stability improves as the distance spacing increases, for any $K$. One has

[TABLE]

with $g_{k}(u)\geq 1$ for all $u$, $l$ and $K$. More precisely we have $g_{1}(u)=1$ for all $u$ and $l$, while $\lim_{u}g_{K}(u)=1$ for all $l$ and $K>1$. In Fig. 5, left panel, the increasing critical time $\tilde{\tau}_{K}(u)$ at the limit $K\rightarrow\infty$ is compared to the time $\tilde{\tau}_{1}(u)$ for $K=1$. One has $\lim_{K}\tilde{\tau}_{K}(u)<\tilde{\tau}_{1}(u)$ for all $u$, while, as expected since the variation tends to zero, $\lim_{u}\tilde{\tau}_{1}(u)=\lim_{u}\tilde{\tau}_{K}(u)$ for any $K$. For the algebraic force, the function $u\mapsto\tilde{\tau}_{K}(u,q)$ increases as $u$ increases since

[TABLE]

The relation $q\mapsto\tilde{\tau}_{K}^{(l)}(u,q)$ is more complicated. One has

[TABLE]

For $K=1$, $h_{1}(q)=0$ for all $q$, and the sign of $\partial\tilde{\tau}_{K}^{(l)}/\partial q$ is the sign of $\ln u-1/q$. It is negative for all $q$ if $u<1$ while it is minimum for $q=1/\ln u$ if $u>1$. Comparable properties are obtained for $K>1$ (see Fig. 5, right panel).

4 Optimal velocity model

The multi-anticipative optimal velocity model [20] with $K\geq 1$ predecessors in interaction is

[TABLE]

The equilibrium speed is given by the function $V$ for this model. We assume in the following that $V^{\prime}(d)>0$ and $1/a_{k}=\tau k^{q}$ with $\tau>0$ the relaxation time with the first predecessor and $q\geq 0$ a parameter calibrating the range of the force (in a similar way than the algebraic force). The stability condition is

[TABLE]

The mean spacing has only a role through the derivative of the optimal speed function that is a scale parameter. Only the parameter $q$ is a form one.

4.1 Stability condition function of $K$

The critical times are proportional to the square of $\sum_{k}k^{-q}$ and converges if and only if $q>1$. Oppositely to the force-based model (13), the critical times not always converge as $K$ increases within the OV model. When $q>1$, the forms of the functions $l\mapsto\tilde{\tau}_{K}^{(l)}$ are comparable to the ones obtained with the algebraic force (see Fig. 6, left panel, and Fig. 2, bottom middle panel). The critical time does not describe damped oscillation with the OV model (see Fig. 6, right panel). This means that the stability improves as the number of predecessors in interaction $K$ increases for any $V^{\prime}(d)>0$ and $q>0$. This in agreement with the results presented in the literature about extended OV models [19, 20, 21, 22, 23, 24]. The critical times converge if $q>1$. In this case, as for the previous pedestrian models, the impact of $K$ on the stability is limited (see also [35]). Note that as previously the critical relaxation times are not minimal for the longest wavelength $l=1$ if $K$ is sufficiently large (e.g. for approximately $K>17$ in Fig. 6, right panel).

4.2 Variation of the stability condition

The proportion of variation $\varphi$ of the critical time (see Eq. (25)), is not defined when $q\leq 1$ with the OV model (29) since in this case $\tilde{\tau}_{K}$ diverges. For $q>1$, the proportion tends to zero as $q$ increases (see Fig. 7, left panel). As expected, and as we observe with the algebraic model (see Fig. 4, right panel) the influence of $K$ decreases as the model becomes short range (i.e. as $q$ increases). The critical time does not depends on $q$ for $K=1$. For any $K>1$, the times decreases as $q$ increases. Therefore the stability is negatively influenced by the range $q$ (see also [36]). The constant minimal value for $K=1$ corresponds here to the limit as $q$ increases for any $K$ (see Fig. 7, right panel). This means that, oppositely to the algebraic model with $u<1$, the OV model can remain stable at the limit $q\rightarrow\infty$.

5 Generalized centrifugal force model

The repulsive force in the generalized centrifugal force (GCF) model [16] depends, as in the extended social force model given in [37], on the distance spacings and also on the speeds of the predecessors. The GCF model is

[TABLE]

with $v_{0},\tau,\eta>0$ and $h:v\mapsto a+bv$ with $a,b\geq 0$ parameters for the size of the pedestrian. Here $\delta$ is a boolean variable for the speed difference term. The model is the centrifugal force one [38] if $a=b=0$ and $\delta=1$, while it is the previous algebraic model (22) with $q=1$ if $\delta=a=b=0$.

The equilibrium speed $v$ solution of (3) is not explicit with the GCF model. For a given mean spacing $d$, it is the solution of

[TABLE]

Here $g(v)\operatorname*{\rightarrow}\infty$ as $v\rightarrow-\infty$ and $g(v)=-\infty$ as $v\rightarrow\infty$. Moreover $g(v)=-\infty$ as $v\rightarrow v_{k}^{-}$ and $g(v)=\infty$ as $v\rightarrow v_{k}^{+}$, for $v_{k}(d)=(kd/2-a)/b$ and $k=1,\ldots,K$. The function $g$ is continuous and strictly decreasing on the subsets $(-\infty,v_{1})$, $(v_{1},v_{2})$,…, $(v_{K-1},v_{K})$, $(v_{K},\infty)$. From the opposite signs of the limits, we deduce by continuity that $g(\cdot)$ admits $K+1$ solutions, one on each subset. We want that the solution $v$ satisfies $d-2h(v)>0$. It is only the case for the solution $v$ belonging to the subset $(-\infty,v_{1}(d))$. We solely consider in the following this solution. It depends on the $v_{0},K,d,\tau$, $a$ and $b$ parameters. The solution is one of the roots of a $K+1$ polynomial equation. It has no explicit definition as soon as $K>1$.

5.1 Stability condition function of $K$

We calculate the stability condition (12) within the GCF model. Denoting

[TABLE]

we have $\beta_{0}=-\frac{1}{\tau}-\sum_{k=1}^{K}A_{k}(2+bA_{k})$, $\beta_{k}=A_{k}(2-bA_{k})$, and $\alpha_{k}=A_{k}^{2}$ for all $k=1,\ldots,K$. The preliminary condition (10) as well as the first condition in (12) hold. The parameters (11) are

[TABLE]

Here the parameters depend on the equilibrium speed $v$ and therefore on $\tau$. Since the analytical form of $v$ is unknown, we can not explicitly express the stability condition on $\tau$ as we did previously. The critical relaxation times $\tau_{K}$ are investigated numerically with $N=10^{5}$. The results show that the linear stability occurs if the parameter $\tau$ is smaller than $\tau_{K}$ and inversely. As before, the stability requires a sufficiently fast relaxation to the desired speed.

The critical times as a function of $l$ are plotted Fig. 8, top panels, for $K=1,\ldots,25$. Here, the values of the parameters are: $\eta=0.3$ and $v_{0}=1.5$ m/s (see [16]). We use $a=0$ in order to obtain the algebraic distance-based model at the limit $b\rightarrow 0$ if $\delta=0$. We compare the cases $b=0.1$ s and $\delta=1$, $b=0$ and $\delta=1$, and $b=0.1$ s and $\delta=0$. The results show that the critical time is minimum for $l=1$, for all $K$ if $\delta=1$, while the critical time is minimum for $l>\!>1$ if $K$ is large when $\delta=0$. As for the exponential model, the longest wavelength is the most unstable if the speed difference term is considered. If $\delta=0$, as for the algebraic and OV models, specific wavelengths are the most unstable if $K$ is sufficiently high.

The critical time for any wavelength $K\mapsto\tau_{K}$ are drawn in Fig. 8, bottom panels. When $\delta=1$, $\tau_{K}$ linearly increases with $K$, and does not converge (further experiments done with higher $K$ confirm it). The stability always occurs if enough predecessors are taken into account. Furthermore, the speed at which $\tau_{K}$ diverges is positively correlated to $b$, i.e. increasing $b$ increases the stability. When the speed difference term is deleted, i.e. for $\delta=0$, $\tau_{K}$ converges (see Fig. 8, bottom right panel, here again further experiments done with higher $K$ confirm it). In this last case, the impact of the parameter $b$ on the stability is not straightforward. It can be positive or negative according to the value of $K$ or mean spacing $d$.

5.2 Stability condition for other parameters

The stability condition (12) is investigated with the GCF model as a function of the mean spacing $d$ and size slope $b$ parameters. As for the distance-based models (18) and (22), the stability always holds as the density tends to zero (see Fig. 9). Yet, the role of $b$ differs according to $d$, $K$ and $\delta$. When $\delta=1$, the parameter $b$ positively impacts the stability, i.e. $\tau_{K}$ increases as $b$ increases, for any $K$. For $K=1$ and $b$ sufficiently high, the function $d\mapsto\tau_{k}(d)$ admits a minimum (see Fig. 9, top left panel). This feature disappears when $K$ increases (see Fig. 9, bottom left panel), the critical time tending to be concave. Different characteristics are obtained when the speed difference term is not taken into account, i.e. for $\delta=0$. When $K$ is small, the parameter $b$ positively influences $\tau_{K}$ for small spacing $d$ (see Fig. 9, top right panel). When $d$ is high (approximately for $d>1$ m in the example), the stability decreases as $b$ increases. When the number $K$ of predecessors in interaction is high (approximately higher that $8$ here), the parameter $b$, as for the case $\delta=1$, always positively impacts the stability (see Fig. 9, bottom middle panel).

The proportion of variation of the critical time $\varphi$ (see Eq. (25) has only a sense if the speed difference term is not taken into account, i.e. for $\delta=0$ (since $\tau_{K}$ diverges for $\delta=1$). Only the parameters $d$ and $b$ are form ones for $\varphi$. It decreases as the spacing $d$ increases (see Fig. 9, bottom right panel). This means that the influence of $K$ on the stability decreases as the density decreases. When $b=0$, the model with $\delta=0$ is the algebraic model (22) for which the variation does not depend on the spacing $d$ (see dotted horizontal line in Fig. 9, bottom right panel). Increasing $b$ emphasizes the impact of $K$. Yet, the proportion tends to the constant minimal value obtained when $b=0$ as $d$ increases for any $b>0$. Oppositely to the exponential model (see Fig.4, left panel), the proportions do not tend to zero and the role of $K$ remains not negligible even if $d\rightarrow\infty$.

6 Summary

The general stability condition provided in Eq. (12) is applied to a large class of car-following and pedestrian models based on relaxation processes. The stability systematically occurs when the relaxation times are sufficiently small. However, the precise nature of the influence of the number of predecessors in interaction $K$ on stability strongly depends on the type of models.

The critical relaxation times converge as $K$ increases for the distance-based models with additive repulsion (as, e.g., the social force model [14]), for any repulsion function as soon as it tends to zero as the spacing increases. Therefore, the parameter $K$ has only limited effects on the stability. For this model class, a cut-off radius in the simulation has no influence on the dynamics if it is sufficiently large. The convergence describes a damped oscillation and the impact of the predecessor number on the stability is not systematically positive. This is not observed with optimal velocity models for which the stability always increases as the predecessor number increases (see, e.g., [20]). The role of the parameter $K$ on the stability is not negligible when the interactions are large (i.e. distance-based models with repulsion $f(d)\propto e^{-dc}$ or $d^{-q}$ or OV model such that $a_{k}=\tau k^{-q}$ with small $c,q$). On the opposite, the stability condition weakly depends on $K$ for short-range interactions (i.e. large $c,q$).

The dynamics of the generalized centrifugal force f(d,v,v_{1})\propto\big{(}\gamma+\delta(v-v_{1})\big{)}^{2}\big{(}d-b(v+v_{1})\big{)}^{-1} is different. When the speed difference term is taken into account (i.e. for $\delta=1$), the critical relaxation time $\tau_{K}$ diverges as $K$ increases. This means that the stability always holds for sufficiently high number of predecessors in interaction. Therefore the cut-off radius influences the dynamics and has to be calibrated as a parameter. The size slope parameter $b$ positively impacts the stability. However, $b$ tends to decrease the equilibrium speed. As for the distance-based models, the critical relaxation time $\tau_{K}$ converges as $K$ increases when the speed difference term is not taken into account (i.e. for $\delta=0$). Yet the speed of convergence can be slow when $b$ and $d$ are high. Here the impact of $b$ on the stability can be either positive, when $K$ is high, or negative, when $K$ is low.

We summarize in Table 1 where the models are distinguished according to: 1. Their range, that can be short or long. 2. The definition of the mean speed in homogeneous configuration, that has in one dimension to be positive. 3. The influence of the number $K$ of predecessors in interaction. The influence can be bounded, unbounded, non-decreasing or non-monotonic. 4. The impact of a cut-off radius on the dynamic. This is a consequence of the previous bounded influence of the predecessor number $K$. The cut-off is implicit when $K$ has limited effects (i.e. bounded influence) while it has to be calibrated as a parameter for unbounded influence of $K$. 5. The most unstable wavelength. The most unstable wavelength can be the longest wave. In such a case the characteristic equation of the linear stability can be solved by expanding linear approximations (as in, e.g., [29, 19, 21, 24]). In some other cases the most unstable wavelengths have specific sizes and the characteristic equation has to be solved for any wavelength.

7 Conclusion

The linear stability conditions of uniform solutions are calculated in one dimension for any wavelength for a large class of acceleration-based microscopic flow models with $K\geq 1$ predecessors in interaction. The framework is general and includes many pedestrian as well as road traffic models. The influence on stability of the number of predecessors in interaction depends strongly on the type of models. It can be both positive or negative, and also irrelevant in certain cases (see Fig. 10).

We observe that the increase of predecessors in the interaction not systematically gives increase of the stability. This result contrasts with the classical properties of OV traffic models for which increase of the interaction range systematically yields a stability improvement [19, 20, 21, 22, 23, 24], as well as to expected features of connected and autonomous vehicles [25, 26, 27].

The stability condition converges as $K$ increases for a particular class of distance-based models that notably includes the social force model [14]. The parameter $K$ is only partially relevant in the dynamics and a cut-off radius during the simulation can be reasonably used.

From a technical side, we point out that the most unstable wavelength is not systematically the longest wave. This is especially the case when an important number of predecessors is taken in account for the interaction. In such a case, linear approximations used to solve the characteristic equation of the linear stability [29, 19, 21, 24] are no more valid and the stability has to be analysed for any wavelength.

Note that the analysis is carried out in one-dimension and for strictly asymmetric interaction. Pedestrian flows and heterogeneous urban traffic stream in two dimensions, and with asymmetric interaction. Multi-dimensional streaming may have stability conditions different from the basic 1D streams and has to be investigated specifically. Furthermore, other speed-based models such that the full velocity difference model [39] remain to be analysed. The impact of the speed difference term on the stability may have different characteristics than the systematic stability improvement we observe with the GCF model.

Acknowledgements

The authors acknowledge financial support from the Deutsche Forschungsgemeinschaft (DFG) under grants SCHA 636/9-1 and SE 1789/4-1.

Appendix 1: Matrix calculation of the linear stability condition

The linear system (8) is

[TABLE]

with

[TABLE]

$T(\mathbf{x})$ being the transpose of the vector $\mathbf{x}$, and, denoting $\varphi_{K}=-\sum_{k=1}^{K}\alpha_{k}$,

[TABLE]

The size of the vectors $\mathbf{Y}$ and $\dot{\mathbf{Y}}$ is $2N$, while $M$ is a $2N\times 2N$ matrix. The solution of this linear system is

[TABLE]

We have to check if $M$ is diagonalizable. If it is, we will have $M=P\mbox{diag}(\lambda_{1},\ldots,\lambda_{2N})P^{-1}$ with $\mbox{diag}(\mathbf{x})$ the diagonal matrix with coefficients $\mathbf{x}$. This allows to obtain

[TABLE]

and to calculate a solution of the system. The results show that the system converge to a stationary state if $M$ is diagonalizable with $\Re(\lambda_{n})\leq 0$ for all $n$. We will determinate the conditions for that the system converges to the vector nil (i.e. a uniform configuration is linearly stable). The matrix $M$ has the form

[TABLE]

The matrix is invariant per circular permutation of $2\times 2$ block, i.e. we have $\theta_{2}(M)=M$ by denoting the operator shifting all the coefficients of a matrix (or a vector) of 2 ranks to the right and 2 ranks to the down

[TABLE]

If $\lambda\in\mathbb{C}$ is an eigenvalue associated to the eigenvector $\mathbf{u}\in\mathbb{C}^{2N}$, then

[TABLE]

We have

[TABLE]

By construction, $\theta_{2}(\mathbf{u})$ is also an eigenvector associated to $\lambda$. This implies

[TABLE]

with $\gamma\in\mathbb{C}$. Denoting $\mathbf{u}=T\left(u_{1},u_{2},\ldots,u_{2N}\right)$, we obtain

[TABLE]

and by iterating

[TABLE]

This directly implies

[TABLE]

as well as $u_{2m+1}=u_{1}\,\gamma_{l}^{m}$ and $u_{2m+2}=u_{2}\,\gamma_{l}^{m}$ for all $m=0,\ldots,N-1$. The eigenvectors are

[TABLE]

with $(u_{1},u_{2})\in\mathbb{C}^{2}$, $\mathbf{u}_{l}^{(1)}=T\left(1,0,\gamma_{l},0,\gamma_{l}^{2},0,\ldots,\gamma_{l}^{N-1},0\right)$ and $\mathbf{u}_{l}^{(2)}=T\left(0,1,0,\gamma_{l},0,\gamma_{l}^{2},\ldots,0,\gamma_{l}^{N-1}\right)$. The eigenvalue $\lambda_{l}$ associated to the eigenvector $\mathbf{u}_{l}$ is the solution of $M\mathbf{u}_{l}=\lambda_{l}\mathbf{u}_{l}$ or again $\lambda_{l}u_{1}=u_{2}$ and $\lambda_{l}u_{2}=u_{1}\varphi_{K}+u_{2}\sum_{k=0}^{K}\beta_{k}\,\gamma_{l}^{k}+u_{1}\sum_{k=1}^{K}\alpha_{k}\,\gamma_{l}^{k}$. This is

[TABLE]

There are two distinct complex eigenvalues $\lambda_{l}^{(1)}$ and $\lambda_{l}^{(2)}$ associated to the eigenvector $\mathbf{u}_{l}$. Since $l$ varies from $1$ to $N$, there are $2\times N$ distinct eigenvalues for $M$. This proves by construction that $M$ is diagonalizable. The equation (A14) is the characteristic equation (9) with $\theta=2\pi l/N$.

For $l=N$, $\gamma_{N}=1$, $\lambda_{N}^{(1)}=\sum_{k=0}^{K}\beta_{k}$ and $\lambda_{N}^{(2)}=0$. If we assume that

[TABLE]

then we have $\lim_{t}e^{Dt}=\mbox{diag}(0,\ldots,0,1)$ and

[TABLE]

with P=\big{[}\mathbf{u}_{1}^{(1)}\ \mathbf{u}_{1}^{(2)}\ldots\ \mathbf{u}_{N}^{(1)}\ \mathbf{u}_{N}^{(2)}\big{]}. Therefore, under the assumptions (A15), $\mathbf{Y}(t)$ converges to the vector nil if

[TABLE]

This equality allows to identify the unique uniform configuration which is linearly stable for the system.

Appendix 2: Solved examples

We illustrate here the linear stability condition (12) with two particular models for which the condition is well known in the literature.

Model with one predecessor

For $K=1$, the model is the second order model with one predecessor in interaction

[TABLE]

investigated in [40]. We suppose that for a mean spacing $d$, exists speed $v$ such that $A(v,d,v)=0$. We have for this model the three parameters $a=\alpha_{1}$, $b=\beta_{0}$ and $c=\beta_{1}$, and, denoting $c_{l}=\cos(2\pi l/N)$ and $s_{l}=\sin(2\pi l/N)$, the parameters (11) are

[TABLE]

The first condition in (12) is

[TABLE]

This condition corresponds to the first one related in [40, see Eq. (B12)]. It holds for all $l$ if $b+|c|<0$. This makes true the preliminary assumption (10). Assuming $a>0$, the second condition of (12) is

[TABLE]

This condition is the one related in [40, Eq. (B19)]. It holds for all $l$ if

[TABLE]

Optimal velocity model

We investigate the multi-anticipative optimal velocity model (29) with $K\geq 1$ predecessors in interaction [20]. With this model, the equilibrium speed $v$ corresponding to the mean spacing $d$ is $v=V(d)$, while $\beta_{0}=-\sum_{k=1}^{K}a_{k}$, $\beta_{k}=0$, for all $k\not=0$, and $\alpha_{k}=V^{\prime}(d)/k$, for all $k$. Denoting $c_{lk}=\cos(2\pi lk/N)$ and $s_{lk}=\sin(2\pi lk/N)$ the parameters (11) are

[TABLE]

The first condition (12) is

[TABLE]

It is always true and implies the preliminary assumption (10). The second condition (12) is, after rearranging

[TABLE]

This is the condition related in [20, see Eq. (14)]. Note that the case $K=1$ corresponds to the well know optimal velocity model [4]. For this model, the condition is

[TABLE]

Since $1/(1+c_{l})>1/2$, the condition holds for all $l=1,\ldots,N-1$ if

[TABLE]

(see [4]).