
Citation: Marek Bodnar, Monika Joanna Piotrowska, Urszula Foryś. Gompertz model with delays and treatment: Mathematical analysis[J]. Mathematical Biosciences and Engineering, 2013, 10(3): 551-563. doi: 10.3934/mbe.2013.10.551
[1] | Matteo Piu, Gabriella Puppo . Stability analysis of microscopic models for traffic flow with lane changing. Networks and Heterogeneous Media, 2022, 17(4): 495-518. doi: 10.3934/nhm.2022006 |
[2] | Maya Briani, Rosanna Manzo, Benedetto Piccoli, Luigi Rarità . Estimation of NO$ _{x} $ and O$ _{3} $ reduction by dissipating traffic waves. Networks and Heterogeneous Media, 2024, 19(2): 822-841. doi: 10.3934/nhm.2024037 |
[3] | Paola Goatin, Chiara Daini, Maria Laura Delle Monache, Antonella Ferrara . Interacting moving bottlenecks in traffic flow. Networks and Heterogeneous Media, 2023, 18(2): 930-945. doi: 10.3934/nhm.2023040 |
[4] | Simone Göttlich, Elisa Iacomini, Thomas Jung . Properties of the LWR model with time delay. Networks and Heterogeneous Media, 2021, 16(1): 31-47. doi: 10.3934/nhm.2020032 |
[5] | Benjamin Seibold, Morris R. Flynn, Aslan R. Kasimov, Rodolfo R. Rosales . Constructing set-valued fundamental diagrams from Jamiton solutions in second order traffic models. Networks and Heterogeneous Media, 2013, 8(3): 745-772. doi: 10.3934/nhm.2013.8.745 |
[6] | Maria Teresa Chiri, Xiaoqian Gong, Benedetto Piccoli . Mean-field limit of a hybrid system for multi-lane car-truck traffic. Networks and Heterogeneous Media, 2023, 18(2): 723-752. doi: 10.3934/nhm.2023031 |
[7] | Xiaoqian Gong, Alexander Keimer . On the well-posedness of the "Bando-follow the leader" car following model and a time-delayed version. Networks and Heterogeneous Media, 2023, 18(2): 775-798. doi: 10.3934/nhm.2023033 |
[8] | Tong Li . Qualitative analysis of some PDE models of traffic flow. Networks and Heterogeneous Media, 2013, 8(3): 773-781. doi: 10.3934/nhm.2013.8.773 |
[9] | Caterina Balzotti, Simone Göttlich . A two-dimensional multi-class traffic flow model. Networks and Heterogeneous Media, 2021, 16(1): 69-90. doi: 10.3934/nhm.2020034 |
[10] | Felisia Angela Chiarello, Paola Goatin . Non-local multi-class traffic flow models. Networks and Heterogeneous Media, 2019, 14(2): 371-387. doi: 10.3934/nhm.2019015 |
Traffic jams that appear for seemingly no reason are a phenomenon that everyone has experienced. Even though there are no accidents, no lane reductions, no highway exits, so called stop-and-go waves form. This paradox has a mathematical answer: under certain conditions, traffic equilibrium states are unstable. This phenomenon has been observed and studied in many works. The consequence of such traffic instabilities are important: accelerating and decelerating increase strongly the fuel consumption and the gas emissions compared to the associated equilibrium flow.
In particular the study of [18], and the experiment associated, provided evidence that such waves form on a circular road with only about twenty drivers starting from same spacing and speed. The impact of the waves includes increased fuel consumption as measured in a similar experiment [17]. Reducing these instabilities by acting on the traffic, for instance rendering the equilibrium stable, would have a strong impact. This is why many approaches have been considered to solve the problem, from ramp-metering control [3,10,15,19,25,26,27] to the use of connected autonomous vehicles which act as "wave-dampeners'' in the traffic (see for instance [4,17,20,22,23,28] and [5] for a more detailed review). But, to be able to act on traffic efficiently, one needs to understand how these traffic waves emerge and their stability.
The main trigger for the creation of waves is the collective behavior of human agents on the road, which typically becomes unstable when the equilibrium velocity cross a certain threshold and becomes too low. They have been extensively studied in the literature, see [7,13,14], and recently explained theoretically in [4] in a ring-road framework similar to [17,18]. These works study a single-phase traffic where all drivers follow the same car-following model. However, in real life, it is likely that the instability of equilibrium flows is strongly influenced by the differences in driving characteristics among agents. These differences could occur due to the type of vehicle (trucks, SUV, small cars, etc.) or differences in behavior between drivers. For instance, trucks platooning was studied for its fuel-consumption impact in mixed traffic [6,24]. Another interesting situation is when some drivers adopt a collaborative way of driving [11]. In this case the traffic can be seen as a system with two populations: the standard drivers and those who adopt the collaborative behavior. Several interesting questions could be asked:
● Can the presence of two or several populations create some instabilities?
● Would a minority of collaborative drivers be able to render a stable traffic that would otherwise be unstable? On the contrary, could a minority of aggressive drivers make an otherwise very stable traffic unstable?
In this paper we provide some answers to these questions. Specifically, a mixed traffic on a ring-road is analyzed, where two or more populations coexists. After characterizing the equilibrium flows, we study the stability of the overall flow depending on the proportion of cars in each class of vehicles. We show that, when an instable and a stable population co-exist, there is a critical penetration rate above which the system is stable and under which the system is unstable, provided there are enough cars on the road. This penetration rate is explicitly computable. We also provide qualitative bounds on the critical penetration rate $ \tau_{0} $, in order to elucidate what makes a collaborating behavior effective from a traffic stability point of view.
For reasonable parameters, the critical penetration rate of stable vehicle is very high. This means that, not taking into account the differences in dynamics, can lead to not understanding the source of instability of the flow. Indeed, we show that a small minority of aggressive drivers, in an otherwise very stable flow, is enough to break the stability. This is something that a single-population model cannot grasp.
Surprisingly, we also show that the ability of a group of drivers to stabilize or destabilize the system only depends on its penetration rate in the traffic and not on the order of the cars. One could think for instance that three trucks in a row is less effective to stabilize the traffic than three trucks equally distributed on the ring-road. But it is not.
Finally, we also show that small experiments, as considered in [17,18], could lead to underestimate the instability of the traffic flow and the critical penetration rate. More precisely, for a small number of vehicles the system can still be stable even below the critical penetration rate, but becomes unstable when the number of cars increases.
The paper is organized as follows: in Section 2 we present the framework and the system, in Section 3 we state the main results, in Section 4 we provide an analysis of the traffic around its equilibrium flows, while in Section 5 we show our results for a two-phase traffic flow and in Section 6 for a general multi-phase traffic flow. Finally, in Section 7 we provide some numerical simulations to illustrate our results.
We study a general traffic model with $ n $ vehicles and a ring-road of length $ L $. Mathematically this means that we consider the system on the domain $ \mathbb{T} = \mathbb{R}/L \mathbb{Z} $. We denote by $ \{x_j(t)\}_{j = 1}^n $ the location of the cars. By further denoting the headway and velocity of the cars as
$ hj(t)=xj+1(t)−xj(t) and vj(t)=˙xj(t), $
|
(2.1) |
the traffic is generally described by
$ {˙hj(t)=vj+1(t)−vj(t),˙vj(t)=fj(hj(t),˙hj(t),vj(t)),∀j∈{1,2,...,n} $
|
(2.2) |
where $ f_{j} $ is the car following model for the driver $ j $ and with the convention $ x_{n+1} = x_1+ L $ on the ring road (or, equivalently, $ x_{n+1} = x_{1} $ in $ \mathbb{T} $), which means that
$ n∑j=1hj(t)=L,∀t≥0. $
|
(2.3) |
For a general car following model we usually have the following physical conditions on $ f_j(h, 0, v) $:
$ ∂∂hfj(h,0,v)>0,∂∂vfj(h,0,v)<0,∂∂˙hfj(h,0,v)>0. $
|
(2.4) |
The first condition simply means that, for a given speed, the incentive to accelerate increases with the headway. The second condition means that for a given headway the incentive to accelerate decreases with the speed. And the third condition means that for a given headway and speed, the incentive to accelerate increases if the headway is currently increasing (i.e., if the leading vehicle is moving away). These conditions can be found in nearly all car following models (e.g., Intelligent Driver Models [21], Follow-the-Leader [8], Bando–Follow-the-Leader [2,4], etc.)
While, in most traffic analysis, the car following model $ f $ is chosen to be identical for all cars on the road, in general the car-following model $ f_j $ depends on the driving habit of the $ j $-th driver, which may differ from driver to driver. It can also depends on the $ j $-th vehicle itself since different kinds of cars may result in different parameters (e.g., one can simply compare a truck and a mini cooper on the road). We summarize three typical cases below:
● Unified model
This is the case that is most commonly considered, where $ f_j $ does not depend on $ j\in \{1, 2, ..., n\} $ and there is a single type of vehicle on the road. The stability of this system has been studied in [4] and their results are recalled below. Under such a unified setting, a further special case can be the so-called Bando–Follow-the-Leader model (Bando-FTL, or equivalently OV-FTL), that combines a Bando (or Optimal Velocity) part introduced in [2] which represents the preference of a driver to reach its "own" optimal velocity (that depends on the headway), and a "Follow-the-Leader" (FTL) part introduced in [8] which represent the incentive of the driver to mimic its leader. More precisely, we have
$ fj(h,˙h,v)=a⋅(V(h)−v)+b⋅˙hh2, $
|
(2.5) |
the Bando part has a weight $ a $ and $ V $ is the optimal velocity function discussed in the next paragraph, while the FTL part has a weight $ b $. This model was studied for instance in [4,5,12].
● Mixed traffic or Collaborative driving
In a mixed traffic, the functions $ \{f_j\}_{j = 1}^n $ are chosen from a finite set. In other words, the drivers and the vehicles can be classified in finitely many categories:
$ fj∈{Fk:k=1,2,...,m},∀j∈{1,2,...,n} $
|
(2.6) |
with $ m\ll n $. A particular example is the two population traffic, where $ m = 2 $ while $ n $ might be large. This situation corresponds for instance to a collaborative driving where some drivers follow a collaborative behaviors and have a "good" function $ F_{1} $ while the rest of the traffic follows a standard function $ F_{2} $ that would lead to instabilities and stop-and-go waves. We deal with this case in Section 5. To illustrate the mixed traffic setting, we can look at the case where the $ \{f_{j}\} $ are given by the "general" Bando-FTL model, characterised by
$ fj(h,˙h,v)=aj⋅(Vj(h)−v)+bj⋅˙hh2, $
|
(2.7) |
$ (aj,bj,Vj)∈{(Ak,Bk,Vk):k=1,2,...,m}, $
|
(2.8) |
In this case, the weights $ (a_j, b_j) $ represents driving habit of the driver, and that the "Bando function" $ V_j(h) $ represents its velocity preference, which can depend for instance of the type of vehicle (cars, trucks, etc.). This optimal velocity is typically taken as a $ C^{2} $ bounded and strictly increasing function of the headway (the farther the leading vehicle, is the faster the velocity preference is likely to be, up to a limit value which is our own speed limit in the absence of traffic). We remark that in the literature for Bando model the function $ V(h) $ is fixed to a ratio of hyperbolic tangents, however, in reality different types of vehicles may have direct influence, that is the reason we call it "general" Bando and adapt it with different $ V_j(h) $ functions. We provide in Section 7 numerical simulations for this particular model.
● Mixed traffic with common velocity preference
This is a particular case of mixed traffic where all drivers have the same velocity preference. For instance for the Bando-FTL model this means that $ V_{j} = V $ and is the same across drivers, while $ a_{j} $ and $ b_{j} $ remain driver-dependent. As we will see in Section 4.1, in this case the equilibrium flows are the same as the equilibrium flows in the unified model.
We are interested in the stability of the system (2.2) around its equilibrium flow. A precise description of the equilibrium flows is given below in Section 4.1. In the unified model, where all drivers have the same car-following model $ f_{j} = f $, the stability of the system was studied in [4]. Recall the physical "common sense" condition in (2.4). For any uniform flow $ (\bar h, \bar v) $, where $ \bar{h} $ and $ \bar{v} $ are both constants that do not depend on time, we shall define
$ (α,β,γ)=(∂f∂h,∂f∂˙h−∂f∂v,∂f∂˙h)(ˉh,0,ˉv), $
|
(3.1) |
satisfying
$ α>0,β>γ>0. $
|
(3.2) |
We have used here the notation convention where $ (f_1, f_2, f_3)(x, y, z) $ denotes $ (f_1(x, y, z), ..., f_3((x, y, z))) $.
The authors of [4] showed the following results.
Theorem 3.1 ([4], unified model). A uniform flow equilibrium $ (\bar{h}, \bar{v}) $ of the system (2.2) is
● locally stable around this flow, if
$ β2−γ2−2α≥0, $
|
(3.3) |
● unstable around this flow provided $ n $ sufficiently large, if
$ β2−γ2−2α<0. $
|
(3.4) |
The definition of local stability in this context is recalled in Definition 4.1. In this paper we investigate what happens when there is not anymore a single population of vehicles but several. In particular, what happens when vehicles with an unstable behavior (i.e., which satisfies Eq (3.4)) coexists with vehicles with a stable behavior (i.e., which satisfies Eq (3.3))?
Our main results are the following: consider first a two population system where the vehicles follow either the car following model $ f^{1} $ or the car following model $ f^{2} $. In this case, any stationary state (or equilibrium flow) can be described by $ (\bar{h}_{i}, \bar{v}) $ with $ i\in\{1, 2\} $ (see Section 4.1). Define
$ (αi,βi,γi)=(∂fi∂h,∂fi∂˙h−∂fi∂v,∂fi∂˙h)(ˉhi,0,ˉv),i∈{1,2}, $
|
(3.5) |
we have the following theorem when both populations have a stable behavior
Theorem 3.2. If $ (\beta^{1})^2-(\gamma^{1})^2-2\alpha^{1}\geq 0 $ and $ (\beta^{2})^2-(\gamma^{2})^2-2\alpha^{2}\geq 0 $, then the steady-state $ (\bar{h}_{i}, \bar{v})_{i\in\{1, 2\}} $ of the ring road system (2.2) and (2.3) is locally exponentially stable.
When the population have different behaviors, we denote $ n_{1} $ and $ n_{2} $ the number of cars of each population in the road, and we have the following theorem
Theorem 3.3. Assume that $ \Delta_{1}: = (\beta^{1})^2-(\gamma^{1})^2-2\alpha^{1} > 0 $ and $ \Delta_{2}: = (\beta^{2})^2-(\gamma^{2})^2-2\alpha^{2} < 0 $. There exists a critical penetration rate $ \tau_0\in (0, 1) $ such that for any pair $ (n_1, n_2)\in \mathbb{N}^2 $ verifying
$ n1n1+n2>τ0, $
|
(3.6) |
the ring road traffic system (2.2) and (2.3) is locally exponentially stable around any equilibrium flow $ (\bar{h}_{i}, \bar{v})_{i\in\{1, 2\}} $, whatever the ordering of the cars.
On the other hand, for any fixed penetration rate
$ τ=n1n1+n2<τ0 $
|
(3.7) |
there exists some $ M > 0 $ effectively computable such that for any $ {n_1, n_2}\in \mathbb{N}^2 $ satisfying $ n_1+n_2 > M $, the ring road traffic system (2.2) and (2.3) is unstable around the equilibrium flow $ (\bar{h}_{i}, \bar{v})_{i\in\{1, 2\}} $.
Moreover, the critical penetration rate $ \tau_{0} $ is explicitly given by
$ τ0=1−(1+max{−H2(y)H1(y):y∈(0,Γ2]})−1,whereHi(y):=log((αi)2+(γi)2y(αi)2+((βi)2−2αi)y+y2),fori∈{1,2},andΓ2:=−(α2)2+√(α2)4−(α2)2(γ2)2Δ2(γ2)2∈(0,−Δ2). $
|
(3.8) |
Even though $ \tau_{0} $ can be computed easily through a minimization algorithm, we can also give some practical upper and lower bounds for qualitative studies:
$ τ0≥−Δ2(α1)2Δ1(α2)2−Δ2(α1)2, and τ0≤(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2)(β2)2(α1)2Δ1+(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2). $
|
(3.9) |
This allows a remark: for a stable class of vehicles to be efficient at stabilizing a mixed traffic flow with a small penetration rate, $ \alpha^{1} $ should be small. This means that an efficient collaborating behavior for stabilizing traffic flow will be composed of vehicles driving without taking much the headway into consideration (apart for safety reasons).
Finally, when one of the populations has a stable behavior but corresponding to the critical case $ (\beta^{1})^2-(\gamma^{1})^2-2\alpha^{1} = 0 $ and the other population has an unstable behavior, we have the following theorem
Theorem 3.4. If $ \Delta_{1}: = (\beta^{1})^2-(\gamma^{1})^2-2\alpha^{1} = 0 $ and $ \Delta_{2}: = (\beta^{2})^2-(\gamma^{2})^2-2\alpha^{2} < 0 $, then provided sufficiently many vehicles on the ring road, the traffic system (2.2)–(2.3) is unstable around the the equilibrium flow $ (\bar{h}_{i}, \bar{v})_{i\in\{1, 2\}} $.
This can be generalized for an multi-phase traffic with more than two populations (see Theorem 6.1 in Section 6).
The aim of this section is to describe the stationary states, or equilibrium flows, of the ring road traffic (2.2) and (2.3). Here, "stationary" and "equilibrium" means that the headway and the velocity do not change with respect to time, hence, for an equilibrium flow $ (\bar{h}_{j}, \bar{v}_{j})_{j\in\{1, ..., n\}} $, we have $ h_j(t) = \bar{h}_j \text{ and } v_j(t) = \bar{v}_j, \; \; \forall t\in[0, +\infty) $. Considering the fact that the headway between two cars does not change, we know that the value of $ \bar v_j $ does not depend on $ j\in \{1, 2, ..., n\} $, therefore
$ hj(t)=ˉhj and vj(t)=ˉv,∀t∈[0,+∞). $
|
(4.1) |
Before going into the details, let us introduce the "velocity preferred headway" of $ f_j $: for a given velocity $ v $, the so-called "velocity preferred headway" is given, when it exists, by $ h $ such that
$ fj(h,0,v)=0, $
|
(4.2) |
which, according to Eq (2.4), admits at most a unique value, and we denote it by $ g_j(v) $ when it exists.
Since the headway does not change, we have
$ ˙hj(t)=˙xj+1(t)−˙xj(t)=vj+1(t)−vj(t), $
|
(4.3) |
which, to be combined with Eq (4.1), imply that there exists a constant $ \bar{v} $ such that,
$ vj(t)=ˉv,∀j∈{1,2,...,n},∀t∈R+ $
|
(4.4) |
$ hj(t)=ˉhj,∀t∈R+. $
|
(4.5) |
Furthermore, thanks to Eq (2.2), we have
$ fj(ˉhj,0,ˉv)=0, $
|
(4.6) |
thus $ \bar{v} $ needs to be chosen in such a way that $ g_{j}(\bar{v}) $ exists for all $ j\in\{1, ..n\} $ and that
$ ˉhj=gj(ˉv). $
|
(4.7) |
Therefore, the equilibrium flow of the traffic is given by,
$ (hj(t),vj(t))=(gj(ˉv),ˉv) $
|
(4.8) |
on a ring road having length
$ L=n∑j=1gj(ˉv). $
|
(4.9) |
Keeping in mind the preceding characterization of equilibrium flow, we are able to address the stability of the ring road traffic around such flows.
Definition 4.1 (Local stability around the equilibrium flow). Let us consider an equilibrium flow, $ \{(\bar x_j(t), \bar v)\}_{j = 1}^n $. The ring road traffic is said to be exponentially stable around this equilibrium flow, if there exist some constant $ C > 0 $ and $ \lambda > 0 $ such that, for any initial state, $ (x_1(0), ..., x_n(0), v_1(0), ..., v_n(0))^{T} $ being sufficiently close to $ (\bar x_1(0), ..., \bar x_n(0), \bar v, ..., \bar v)^{T} $, the traffic satisfies
$ |(x1(t)−ˉx1(t)−c,...,xn(t)−ˉxn(t)−c,v1(t)−ˉv,...,vn(t)−ˉv)| $
|
(4.10) |
$ ≤Ce−λt|(x1(0)−ˉx1(0),...,xn(0)−ˉxn(0),v1(0)−ˉv,...,vn(0)−ˉv)| $
|
(4.11) |
with some constant $ c\in \mathbb{R} $ depending on the initial state. The constant $ c $ comes from the fact that an equilibrium flow defines a headway and a velocity, but the location of the cars is only defined up to a constant.
We remark here that this definition is equivalent to Definition 4.5 given below of the stability of the traffic in terms of headway-velocity.
Remark 4.1 (Ring road condition). We remark that for fixed $ \{f_j\}_{j = 1}^n $ and $ L $ there is at most one equilibrium flow, $ i.e., $ at most one value of $ \bar v $ such that Eq (4.9) holds. This means that the length of the road imposes the equilibrium flow and the steady velocity $ \bar{v} $ and reciprocally that imposing a given speed $ \bar{v} $ and a given system of $ n $ cars determines the length of the road. Since in reality we may want to study the limit where both $ n $ and $ L $ go to $ \infty $, what really matters for us is the desired steady velocity. For this reason we do not fix the value of $ L $ in this paper but we fix instead the value of $ \bar{v} $, which in turns, defines the value of $ L $ via Eq (4.9). This situation will be named as "ring road condition" in the following.
Another important property of the "ring road condition" is that the value of $ L $, Eq (4.9), does not depend on the order of $ \{f_j\}_{j = 1}^n $, namely if $ \{\tilde f_j\}_{j = 1}^n = \{f_j\}_{j = 1}^n $ up to some permutations then their lengths of ring road coincide.
Remark 4.2 (Number of parameters describing the equilibrium flow). Looking at Eq (4.6), the stationary headway $ \bar h_j = g_j(\bar v) $ only depends on the driving characterization function $ f_j $. This means that it is identical for all vehicles with the same class of parameters. In particular
● For unified model, $ g_j(\bar v) $ does not depend on $ j\in \{1, 2, ..., n\} $. Thus, only 2 parameters describe the n-vehicles equilibrium flow: $ \bar v $ and $ g(\bar{v}) $.
● For general collaborative driving as indicated in Eq (2.6) or Eq (2.7) and (2.8), $ g_j(\bar v) $ has $ m $ different choices. Hence, the $ n $-vehicles equilibrium flow is described by m+1 parameters.
● For the unified car Bando-FTL model of collaborative driving described by Eqs (2.7) and (2.8) with $ V_j = V $, where the drivers have different driving habits but the same velocity preference, $ g_j(\bar v) $ does not depend on $ j\in \{1, 2, ..., n\} $ and the equilibrium flow is still only described by two parameters, despite different driving habits. Actually, in this particular case, we observe from Eq (2.7) that $ g_j(\bar v) = V^{-1}(\bar v) $. As a direct consequence, the "ring road condition" Eq (4.9) simply becomes $ L = n V^{-1}(\bar v) $.
In this section we describe the linearized system around the equilibrium flows presented in the previous section. Let a general traffic model be given by Eqs (2.2) and (2.3). Let $ (\bar{h}_{i}, \bar v) $ be an equilibrium flow. From the previous Section, $ (\bar{h}_{i}, \bar v) $ satisfies (4.8) and (4.9). In particular, we remark that $ (\bar{h}_{i}, \bar v) $ is entirely parametrized by $ \bar{v} $.
By denoting the perturbation in headway and velocity as
$ yj(t)=hj(t)−ˉhj and uj(t)=vj(t)−ˉv, $
|
(4.12) |
the traffic flow characterized by Eqs (2.2) and (2.3) can be written in terms of $ (y_j, u_j) $:
$ {˙yj=uj+1−uj,˙uj=fj(ˉhj+yj,uj+1−uj,ˉv+uj), $
|
(4.13) |
and Eq (2.3) becomes
$ n∑j=1yj(t)=n∑j=1hj(t)−n∑j=1ˉhj=0. $
|
(4.14) |
This condition is equivalent to say that the solutions of the preceding system always stay in the $ (2n-1) $ dimensional subspace of $ \mathbb{C}^{2} $:
$ H:={(y1,y2,...,yn,u1,u2,...,un)∈C2n:n∑k=1yk=0}, $
|
(4.15) |
in particular, when the initial state takes value from $ \mathcal{H}\cap \mathbb{R}^{2n} $, the solution also stays in $ \mathcal{H}\cap \mathbb{R}^{2n} $. Here, we introduce complex valued spaces for the ease of notations when considering eigenvectors and eigenvalues.
Therefore the system (4.13) and (4.14) can be expressed as
$ {˙yj=uj+1−uj,˙uj=fj(ˉhj+yj,uj+1−uj,ˉv+uj),(yj(t),uj(t))∈H. $
|
(4.16) |
Linearized traffic around equilibrium flow
Next, standard linearization yields the linearized ring road traffic system
$ {˙yj=uj+1−uj,˙uj=αjyj−βjuj+γjuj+1,(yj(t),uj(t))∈H, $
|
(4.17) |
where $ (\alpha_j, \beta_j, \gamma_j) $ (independent of time) are given by
$ (αj,βj,γj)=(∂fj∂h,∂fj∂˙h−∂fj∂v,∂fj∂˙h)(gj(ˉv),0,ˉv) $
|
(4.18) |
which satisfies, from Eq (2.4),
$ αj>0,βj>γj>0. $
|
(4.19) |
Form now on, for the ease of notations, we will denote by $ \Lambda_{i} $ the set of parameters describing the car following model $ f_{i} $ and $ \Delta_{j} $ the quantity describing its stability around the equilibrium flow.
Definition 4.2. For any given trio
$ Λj:=(αj,βj,γj) $
|
(4.20) |
we define its discriminant as
$ Δj:=(βj)2−(γj)2−2αj. $
|
(4.21) |
The expression of $ \Delta_{j} $ comes from [4] (see Theorem 3.1), it describes the stability of $ f_{j} $ around an equilibrium flow in a single-phase traffic, i.e., when all cars have the same car-following model $ f : = f_{j} $. What we will see in the following is that $ \Delta_{j} $ is also a good indicator of the stability in a multi-phase traffic. For this reason, we introduce the following definition, inspired by Theorem 3.1,
Definition 4.3. We classify a trio, $ \Lambda : = (\alpha, \beta, \gamma)\in \mathbb{R}^3 $ satisfying the "common sense" condition (3.2), by the value of its discriminant as follows. We say that $ \Lambda $ is
● stable, if $ \Delta : = \beta^2-\gamma^2- 2\alpha > 0 $.
● critical, if $ \Delta = 0 $.
● unstable, if $ \Delta < 0 $.
Finally, for a general traffic model (2.2) with $ n $ cars satisfying the ring road condition (4.1), note that there is at most $ n $ different set of parameters $ \Lambda_j $ (one per car). However, if we restrict it into some $ m $-phase mixed traffic model (2.6), where all the vehicles can be classified in $ m $ different types, then, thanks to Eqs (4.6)–(4.8), $ \Lambda_j $ only have at most $ m $ different values.
Let us now denote
$ z(t)=(y1,...,yn,u1,...un)T(t)∈H $
|
(4.22) |
the linearized traffic system (4.17) becomes
$ ˙z(t)=Mnz(t),z(t)∈H, $
|
(4.23) |
with
$ Mn:=(OnAnCnBn) $
|
(4.24) |
where
$ On:=(000......0),An:=(−11−11−11......1−1). $
|
(4.25) |
$ Cn:=(α1α2α3......αn),Bn:=(−β1γ1−β2γ2−β3γ3......γn−βn). $
|
(4.26) |
The following definitions describe stability properties of both the linearized and nonlinear systems around equilibrium flows.
Definition 4.4 (Linearized stability around the equilibrium flow). Let us consider an equilibrium flow on the ring road: $ \{(\bar h_j, \bar v)\}_{j = 1}^n $. The linearized ring road traffic (4.17) is said to be exponentially stable around this equilibrium flow, if there exist some $ \lambda > 0 $ and $ C > 0 $ such that for any initial state, $ z(0) = (y_1(0), ..., y_n(0), u_1(0), ..., u_n(0))^{T}\in \mathcal{H} $, the solution of the Cauchy problem (4.23) (by recalling Eqs (4.17)–(4.26)) satisfies
$ |z(t)|≤Ce−λt|z(0)|,∀t∈R+. $
|
(4.27) |
The preceding definition describes the stability of the linearized system (4.17). Actually, thanks to the standard linearization argument, when the linearized system is stable in the sense of Definition 4, the original nonlinear system (4.13) is automatically locally stable in $ \mathcal{H} $ in the following sense:
Definition 4.5 (Local stability around the equilibrium flow: an alternative definition). Let us consider an equilibrium flow of the ring road: $ \{(\bar h_j, \bar v)\}_{j = 1}^n $. The ring road traffic (4.16) is said to be locally exponentially stable around this equilibrium flow, if there exist some $ \varepsilon > 0 $, $ \lambda > 0 $ and $ C > 0 $ such that for any initial state, $ z(0) = (y_1(0), ..., y_n(0), u_1(0), ..., u_n(0))^{T}\in \mathcal{H} $, satisfying
$ |z(0)|≤ε, $
|
(4.28) |
the solution of the Cauchy problem (4.13) satisfies
$ |z(t)|≤Ce−λt|z(0)|,∀t∈R+. $
|
(4.29) |
Let us remark that in the preceding two definitions, it is equivalent to express the stability in terms of position-velocity around the equilibrium flow, $ \{(\bar x_j(t), \bar v)\}_{j = 1}^n $, like Definition 4.1.
Remark 4.3. One can note that we do not look here at the well-posedness in general of the nonlinear system. The reason is that we look at the local stability around steady-states, and the well-posedness in this framework is immediate if the stability guaranteed. However, it would be interesting to know whether the nonlinear system is well-posed for any initial condition or any reasonable initial condition. When using Bando-Follow the leader model, the single population system has been shown to be well-posed as long as $ h_{j} $ are initially above a lower bound in [9] (for the IDM the situation is more complicated, one can look at [1] for instance). It would be interesting to investigate whether the same type of argument would work in the multi-population case.
Since the local exponential stability of the nonlinear system (4.16) can be directly deduced from the exponential stability of the linearized system (4.23) we are going to focus on the latter. Using Routh–Hurwitz criterion the exponential stability of the linearized system (4.23) depends on the spectrum of the operator
$ L:H→Hx↦Mnx. $
|
(4.30) |
Therefore we investigate the spectrum of the matrix $ M_{n} $ on $ \mathcal{H} $.
In this section we provide a two-step approach to find the eigenvalues of the matrix $ M_{n} $. In particular, the second approach gives a full description of the multiplicity of the eigenvalues. We first look at the spectrum of $ M_{n} $ on $ \mathbb{C}^{2n} $ and then see which eigenvalue remains when we restrict $ M_{n} $ to $ \mathcal{H} $.
The spectrum of $ M_n $ on $ \mathbb{C}^{2n} $.
A first method consists in considering the eigen-pairs $ (\omega, z) $ of the matrix $ M_n $:
$ Mnz=ωz. $
|
(4.31) |
We know from the preceding equation that, for any $ j\in \{1, 2, ..., n\} $,
$ ωyj=uj+1−uj, $
|
(4.32) |
$ ωuj=αjyj−βjuj+γjuj+1, $
|
(4.33) |
which further implies
$ (ω2+βjω+αj)uj=(γjω+αj)uj+1. $
|
(4.34) |
Therefore
$ n∏j=1Fj(ω)=1,Fj(ω):=γjω+αjω2+βjω+αj, $
|
(4.35) |
which algebraically provides $ 2n $ solutions (with multiplicity). However, this approach does not yet give any information on the multiplicity of the eigenvalues: indeed, we only know that any eigenvalue $ \omega $ should satisfy Eq (4.35). Thus, we use a second method to obtain this piece of information.
Recall that the full spectrum information of a matrix is given by the roots of the characteristic polynomial. Thus we shall this direct method to find all eigenvalues (counting multiplicity) of the matrix $ M_n $. We are going to compute
$ χ(λ):=det(λIn−An−CnλIn−Bn) $
|
(4.36) |
Notice that $ \lambda I_n C_n = C_n \lambda I_n $, we know that
$ det(λIn−An−CnλIn−Bn)=det((λIn)(λIn−Bn)−(−Cn)(−An))=det(λ2In−λBn−CnAn)=det(˜Mn), $
|
(4.37) |
where, $ \widetilde{M}_n $ is given by
$ (λ2+β1λ+α1−γ1λ−α1λ2+β2λ+α2−γ2λ−α2λ2+β3λ+α3−γ3λ−α3......−γnλ−αnλ2+βnλ+αn). $
|
(4.38) |
By considering the first column of the matrix, its determinant read as
$ χ(λ)=det(˜Mn)=n∏i=1(λ2+βiλ+αi)+(−1)n−1(−γnλ−αn)n−1∏j=1(−γjλ−αj)=n∏i=1(λ2+βiλ+αi)−n∏j=1(γjλ+αj). $
|
(4.39) |
Therefore, all the eigenvalues are given by
$ n∏j=1γjλ+αjλ2+βjλ+αj=1, $
|
(4.40) |
which is coincident with Eq (4.35).
Remark 4.4 (Independence of the stability with the order of the cars). Similarly to the length of the road for a given steady-state (see Remark 4.1), again, we remark that the spectrum information of the traffic does not depend on the order of $ \{f_j\}_{j = 1}^n $. Namely, for a given desired steady velocity $ \bar v $, if $ \{\tilde f_j\}_{j = 1}^n = \{f_j\}_{j = 1}^n $ up to some permutations, then their stability around the related equilibrium flows coincide. This means that, in the mixed-traffic setting, the stability only depends on the penetration rate of the different types of cars.
The spectrum of $ M_{n} $ on $ \mathcal{H} $
By looking at Eq (4.40), it is easy to observe that $ \lambda = 0 $ is an eigenvalue of the matrix $ M_n $. In the following, we prove that $ \lambda = 0 $ is a simple eigenvalue of the matrix $ M_n $ acting on $ \mathbb{C}^{2n} $ however, it is not a eigenvalue of $ M_n $ acting on $ \mathcal{H} $. This is an important point as, otherwise, we would not be able to deduce the stability of the system (4.23) from the eigenvalue analysis.
Let us start by showing that $ \lambda = 0 $ is a simple eigenvalue of the matrix $ M_n $. By comparing the coefficients of the characteristic polynomial $ \chi(\lambda) $ given by Eq (4.39), we immediately notice that $ \chi(0) = 0 $. Then it suffices to show that $ \chi'(0)\neq 0 $. Indeed, suppose that $ \lambda = 0 $ has (at least) multiplicity two, then the characteristic polynomial can be written as $ \chi(\lambda) = \lambda^{2}P(\lambda) $ where $ P $ is again a polynomial, and consequently $ \chi'(0) = 0 $. Proving that $ \chi'(0)\neq 0 $ is equivalent to prove that
$ n∑i=1βi(∏nk=1αkαi)≠n∑i=1γi(∏nk=1αkαi), $
|
(4.41) |
which is guaranteed by the "common sense" condition (4.19).
Next, we show that even though 0 is a simple eigenvalue of the matrix $ M_n $, it is not included in the finite spectrum of the operator $ M_n $ acting on $ \mathcal{H} $. Indeed, suppose by contradiction that there exists $ z = (y_1, ..., y_n, u_1, ..., u_n)^{T}\in \mathcal{H} $ such that $ z $ is an eigenvector of $ M_{n} $ associated to the eigenvalue $ 0 $. We have
$ Mnz=0. $
|
(4.42) |
As $ z\neq 0 $ and using Eqs (4.25) and (4.26) we deduce that there exists some $ C\neq 0 $ such that for any $ j\in \{1, 2, ..., n\} $,
$ uj=C,yj=βj−γjαjC. $
|
(4.43) |
Without loss of generality, we assume that $ C > 0 $. By recalling the "common sense" condition (4.19), this implies that
$ n∑j=1yj>0, $
|
(4.44) |
but as $ z\in\mathcal{H} $ we know from Eq (4.15) that
$ n∑j=1yj=0, $
|
(4.45) |
which leads to a contradiction. This implies that
$ SpH(Mn)⊆SpC2n(Mn)∖{0}, $
|
(4.46) |
where $ \text{Sp}_{\mathcal{H}}(M_{n}) $ is the spectrum of $ M_{n} $ on $ \mathcal{H} $ and $ \text{Sp}_{\mathbb{C}^{2n}}(M_{n}) $ is the spectrum of $ M_{n} $ on $ \mathbb{C}^{2n} $. On the other hand, for any $ \lambda\in \text{Sp}_{\mathbb{C}^{2n}}(M_{n})\setminus\{0\} $ we can find at least one related eigenvector $ z\in \mathbb{C}^{2n} $. It is clear that $ M_n z\in \mathcal{H} $, thus $ z\in \mathcal{H} $. Therefore,
$ SpH(Mn)=SpC2n(Mn)∖{0}. $
|
(4.47) |
In this section, we study the stability of the equilibrium flows in a two-phase traffic flow. This situation represents for instance two class of vehicles such as trucks and cars, or also the coexistence of vehicles with and without a collaborative driving behavior. In the following, these two classes of vehicles will be called Type 1 vehicles and Type 2 vehicles. Let $ \bar{v} $ be an equilibrium velocity, from Eq (4.8) this imposes the equilibrium headway $ \bar h_1 $ ($ resp. $ $ \bar h_2 $) of the Type 1 vehicles ($ resp. $ Type 2 vehicles). Thus, we are looking at a situation where, for every $ j\in \{1, 2, ..., n\} $, using the notation (4.20) and (4.21),
$ Λj=(αj,βj,γj)∈{Λ1,Λ2}={(α1,β1,γ1),(α2,β2,γ2)}. $
|
(5.1) |
Let us denote by $ n_{1} $ the number of Type 1 vehicles with parameters $ \Lambda^{1} $ and $ n_{2} $ the number of Type 2 vehicles with parameters $ \Lambda^{2} $ such that the total number of vehicles is $ n = n_1+n_2 $.
Suppose that the mixed traffic on road is represented by the "ordering" $ (a_1, a_2, ..., a_n) $ that belongs to
$ K:={(a1,a2,...,an):ak∈{1,2}∀1≤k≤n,n∑k=1ak=n1+2n2}, $
|
(5.2) |
where $ a_k\in \{1, 2\} $ implies that the $ k $-th vehicle on the road is of Type $ a_k $: because there is no lane changing in a single ring road, $ (a_1, a_2, ..., a_n) $ is invariant with respect to time. Consequently, there is a unique equilibrium flow corresponding to $ \bar{v} $ (or equivalently there is a unique equilibrium flow corresponding to $ L $, from Remark 4.1): the headway before the $ k $-th vehicle is given by $ \bar h_{a_k} $. Furthermore, from Section 4.2, the linearized system around this equilibrium flow is
$ {˙yj=uj+1−uj,˙uj=αjyj−βjuj+γjuj+1,(αj,βj,γj)=(αaj,βaj,γaj),(yj(t),uj(t))j∈{1,...,n}∈H, $
|
(5.3) |
which can be further represented in forms of Eqs (4.23)–(4.26).
We introduce the following Lemma
Lemma 5.1. Let given $ (n_1, n_2)\in \mathbb{N}^2 $. If the following inequality holds,
$ ((α1)2+(γ1)2x2(α1)2+((β1)2−2α1)x2+x4)n1((α2)2+(γ2)2x2(α2)2+((β2)2−2α2)x2+x4)n2<1,∀x∈R∖{0}, $
|
(5.4) |
then System (5.2) and (5.3) is exponentially stable in the sense of Definition 4.4.
On the other hand, if for some $ x\in \mathbb{R} $ we have
$ ((α1)2+(γ1)2x2(α1)2+((β1)2−2α1)x2+x4)n1((α2)2+(γ2)2x2(α2)2+((β2)2−2α2)x2+x4)n2>1, $
|
(5.5) |
then there exists some $ M_0\in \mathbb{N} $ effectively computable such that for any integer $ M\geq M_0 $, the System (5.2) and (5.3) with $ (n_1, n_2) $ being replaced by $ (M n_1, Mn_2) $ is unstable.
Proof. This proof is essentially the same as the one given by [4,Section II] in a unified models framework (namely $ n_2 = 0 $). For readers' convenience we sketch its proof as follows.
By representing System (5.2) and (5.3) in form of (4.23)–(4.26), thank to Section 4.3, the eigenvalues (counting multiplicity) are explicitely characterized by Eq (4.40):
$ (γ1λ+α1λ2+β1λ+α1)n1(γ2λ+α2λ2+β2λ+α2)n2=1. $
|
(5.6) |
Inspired by [4,Section II], we consider the following meromorphic function
$ G(z)=(γ1z+α1z2+β1z+α1)n1(γ2z+α2z2+β2z+α2)n2,∀z∈C. $
|
(5.7) |
Since all the poles are located on the left half plane, $ G(z) $ is holomorphic on the right half plane $ \mathbb{C}^{+} = \{z\in \mathbb{C}: \Re(z)\geq 0\} $. We notice that $ |G(z)| $ tends to 0 as $ |z| $ tends to $ +\infty $. Then, thanks to the maximum principle of holomorphic functions, the maximum of $ |G(z)| $ over $ \mathbb{C}^{+} $ must takes place at the imaginary axis. By considering $ z = ix $ we get
$ |G(z)|2=((α1)2+(γ1)2x2(α1)2+((β1)2−2α1)x2+x4)n1((α2)2+(γ2)2x2(α2)2+((β2)2−2α2)x2+x4)n2, $
|
(5.8) |
which explains the inequality (5.4).
If Condition (5.4) is satisfied, then we know that $ |G(z)|\leq 1 $ in $ \mathbb{C}^{+} $ with $ |G(z)| $ equals to 1 only at $ z = 0 $. Therefore, all the eigenvalues are located in $ \{z\in \mathbb{C}: \Re(z) < 0\}\cup \{0\} $, which, to be combined with Eq (4.47), yields the required exponential stability.
On the other hand, if Condition (5.5) is satisfied, then $ |G(z)| $ is not always smaller than 1. As $ n_{1}+n_{2}\neq0 $, without loss of generality we assume that $ n_{2}\neq 0 $ and we define $ \mu = n_{1}/n_{2} $. Condition (5.5) implies
$ ((α1)2+(γ1)2x2(α1)2+((β1)2−2α1)x2+x4)μ((α2)2+(γ2)2x2(α2)2+((β2)2−2α2)x2+x4)>1. $
|
(5.9) |
We denote
$ Gμ(z)=(γ1z+α1z2+β1z+α1)μ(γ2z+α2z2+β2z+α2),∀z∈C. $
|
(5.10) |
We consider the curve $ \mathcal{C}\subset \mathbb{C} $:
$ C:={z∈C:|Gμ(z)|=1}. $
|
(5.11) |
Let us further define an open subset of $ \mathcal{C} $ by
$ C+:={z∈C:ℜ(z)>0}, $
|
(5.12) |
which is not empty thanks to Condition (5.9), the fact that $ \lim_{\text{Re}(z)\rightarrow +\infty} |G_{\mu}(z)| = 0 $ and the continuity of $ G_{\mu} $. It is natural to consider the continuous function on $ \mathcal{C}^{+} $ defined as
$ F1:C+→S1z↦F1(z):=Gμ(z), $
|
(5.13) |
where $ \mathbb{S}^1 $ denotes the unit circle in the complex plane $ \mathbb{C} $. We can find a connected open set $ \mathcal{O}\subset \mathbb{S}^1 $ such that
$ O⊂F1(C+). $
|
(5.14) |
By denoting the length of $ \mathcal{O} $ as $ |\mathcal{O}| $, the value of $ M_0 $ can be chosen as
$ M0:=[2π|O|]+1. $
|
(5.15) |
Indeed, for any $ M\geq M_0 $, we know from the definition of $ M $ that
$ 2πM<|O|. $
|
(5.16) |
Therefore, the set
$ {e2kπ/M:k=1,2,...,M}∩O $
|
(5.17) |
is not empty. We assume that for some $ k $, the point $ e^{2k\pi/M} $ belongs to $ \mathcal{O} $. Thus, by the definition of $ \mathcal{O} $ there exists some $ z_0\in \mathcal{C}^+ $ such that
$ Gμ(z0)=F1(z0)=e2kπ/M. $
|
(5.18) |
Meanwhile, we recall that for the ring road traffic with $ (M \mu) $ Type 1 vehicles and $ M $ Type 2 vehicles the stability of the System (5.2) and (5.3) is determined by the solutions of $ G_M(z) = 1 $:
$ GM(z)=(γ1z+α1z2+β1z+α1)Mμ(γ2z+α2z2+β2z+α2)M. $
|
(5.19) |
The preceding equations immediately yield $ G_M(z_0) = 1 $ with $ \Re(z_0) > 0 $: the system is unstable.
This lemma allows us to show the Theorems 3.2–3.4 that we reformulate here for convenience in Theorem 5.1, 5.2:
Theorem 5.1. Let given $ \bar v > 0 $. If $ \Delta_{1}\geq 0 $ and $ \Delta_{2}\geq 0 $, then for any $ (n_1, n_2)\in \mathbb{N}^2 $ and any ordering of the vehicles on the road, the ring road traffic system (4.16) is locally exponentially stable around the equilibrium flow.
Proof of Theorem 5.1. The proof is straightforward. Indeed, by the definition of $ \Delta_{1} $ and the assumption that $ \Delta\geq 0 $ we know that for any $ x\in \mathbb{R}\setminus\{0\} $ there is
$ ((α1)2+(γ1)2x2(α1)2+((β1)2−2α1)x2+x4)n1((α2)2+(γ2)2x2(α2)2+((β2)2−2α2)x2+x4)n2<1. $
|
Similarly,
$ ((α2)2+(γ2)2x2(α2)2+((β2)2−2α2)x2+x4)n2<1, $
|
thus Inequality (5.4) holds. The proof of this theorem concludes by applying Lemma 5.1.
Remark 5.1. Note that, in the critical case, i.e. one or both of $ \Delta_{1} $ and $ \Delta_{2} $ equals to 0, the system is still exponentially stable.
Theorem 5.2. Let given $ \bar v > 0 $. We assume that $ \Delta_{1}\geq 0 $ and $ \Delta_{2} < 0 $.
(1) If $ \Delta_{1} > 0 $ then there exists some effectively computable threshold constant $ \tau_0\in (0, 1) $ depending on $ \Lambda^{1} $ and $ \Lambda^{2} $ and given by Eq (3.8) such that for any pair $ (n_1, n_2)\in \mathbb{N}^2 $ verifying
$ n1n1+n2>τ0, $
|
(5.20) |
the inequality (5.4) is satisfied. In other words, for any ordering of the vehicles on the road $ (a_1, a_2, ..., a_n)\in \mathcal{K} $, the ring road traffic system (4.16) is locally exponentially stable around the equilibrium flow associated to $ \bar{v} $.
On the other hand, for any penetration rate
$ τ:=n1n1+n2<τ0 $
|
(5.21) |
there exists $ M > 0 $ such that for any $ {n_1, n_2}\in \mathbb{N}^2 $ satisfying $ n_1+n_2 > M $, the ring road traffic system (2.2) and (2.3) is unstable around the equilibrium flow $ (\bar{h}_{i}, \bar{v})_{i\in\{1, 2\}} $.
(2) If $ \Delta_{1} = 0 $ (namely, $ \Lambda^{1} $ is critical), then for any penetration rate $ \tau = n_{1}/(n_{1}+n_{2}) $ there exists $ M > 0 $ effectively computable such that if $ n > M $ the ring road traffic system (4.16) is unstable around the equilibrium flow associated to $ \bar{v} $.
Finally, even though $ \tau_{0} $ can be easily calculated with the help of a minimization algorithm we show some simpler upper and lower bounds.
Corollary 5.3. The critical penetration rate $ \tau_{0} $ defined in Theorem 5.2 satisfies
$ τ0≥−Δ2(α1)2Δ1(α2)2−Δ2(α1)2,andτ0≤(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2)(β2)2(α1)2Δ1+(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2). $
|
(5.22) |
As we can see that Theorem 3.2–3.4 in Section 3 are direct consequences of the more detailed Theorems 5.1, 5.2, in the following we only give the proofs of the latter theorems.
Proof of Theorem 5.2. Note that it suffices to study the exponential stability of the linearized system (5.3) since the local exponential stability of the nonlinear system follows: At first we prove point (1) of this theorem. Looking at Lemma 5.1, it is thus sufficient to investigate Eq (5.4). Let us first get an intuition about what happens depending on the proportion of stable and unstable vehicle. To simplify the condition, we can set $ y = x^{2} $ and $ \mu = n_{1}/n_{2} $, namely $ \mu/(1+\mu) $ is the proportion of stable vehicle in the traffic. The condition (5.4) becomes
$ ((α1)2+(γ1)2y(α1)2+((β1)2−2α1)y+y2)μ((α2)2+(γ2)2y(α2)2+((β2)2−2α2)y+y2)<1,∀y∈R∗+. $
|
(5.23) |
We set
$ h1(y)=((α1)2+(γ1)2y(α1)2+((β1)2−2α1)y+y2),∀y∈R∗+. $
|
(5.24) |
and we $ h_{2} $ similarly. We observe that, for $ i\in \{1, 2\} $,
$ h′i(y)=−(γi)2y2−2(αi)2y−(αi)2Δi((αi)2+((βi)2−2αi)y+y2)2,∀y∈R∗+. $
|
(5.25) |
where, we recall that $ \Delta^{i} = (\beta^{i})^2-(\gamma^{i})^2-2\alpha^{i} $. This means that $ h_{i} $ has at most two points where its derivative vanishes and these potential points are given by
$ y±=−(αi)2(γi)2(1∓√1−(γi)2(αi)2Δi). $
|
(5.26) |
However, note that we are only interested in the values of $ h_{i} $ on $ [0, +\infty) $. This gives some insight about what happens when $ \Delta^{i} $ moves from a positive value (stable region) to a negative value (possibly unstable region): when $ \Delta^{i} $ is positive there is no non-negative vanishing points of $ h_{i}' $, which means that $ h_{i} $ start at the value $ h_{i}(0) = 1 $ and then decrease strictly continuously until it reaches the limit $ \lim_{y\rightarrow +\infty} h_{i}(y) = 0 $. When $ \Delta_{i} $ is negative, then $ y_{+} $ is the only positive vanishing point of $ h' $ which means that $ h_{i} $ still starts at the value $ h_{i}(0) = 1 $ but increase strictly up to $ y = y_{+} $ and becomes larger than $ 1 $. The critical case $ \Delta^{i} = 0 $ corresponds to the special case where $ y_{+} = 0 $ and therefore $ h $ still decreases strictly on $ [0, +\infty) $. Let us now prove (1) of Theorem 5.2.
Quantitative characterization of the optimal choice of $ \tau_0 $
We define
$ Hi(y)=log(hi)=log((αi)2+(γi)2y(αi)2+((βi)2−2αi)y+y2),i∈{1,2}. $
|
(5.27) |
then for any given $ y > 0 $, (5.4) is equivalent to
$ μH1(y)+H2(y)<0, ∀y>0. $
|
(5.28) |
Since $ H_1(y) < 0 \; \forall y > 0 $, the preceding condition is equivalent to having
$ μ>−H2(y)H1(y), ∀y>0. $
|
(5.29) |
Therefore, it leads us to introduce the quantity
$ K:=sup{−H2(y)H1(y):y∈(0,+∞)}, $
|
(5.30) |
and to show that this quantity is finite. If so, it suffices to choose $ \mu > K $ to conclude the exponential stability. We will show the following: define
$ N0=max{−H2(y)H1(y):y∈(0,Γ2]}<+∞, $
|
(5.31) |
$ Γ2:=−(α2)2+√(α2)4−(α2)2(γ2)2Δ2(γ2)2∈(0,−Δ2), $
|
(5.32) |
$ τ0:=N0N0+1, $
|
(5.33) |
then $ K = N_{0} < +\infty $ and if
$ n1n1+n2>τ0, $
|
(5.34) |
the condition (5.29) is satisfied and consequently the system is exponentially stable around the considered equilibrium flow.
From Eq (5.25)
$ H′1(y)=−(γ1)2y2−2(α1)2y−(α1)2Δ1((α1)2+(γ1)2y)((α1)2+((β1)2−2α1)y+y2),∀y∈R∗+, $
|
(5.35) |
$ H′2(y)=−(γ2)2y2−2(α2)2y−(α2)2Δ2((α2)2+(γ2)2y)((α2)2+((β2)2−2α2)y+y2)∀y∈R∗+. $
|
(5.36) |
And using this together with Eq (5.27), we deduce that
$ H1(0)=H2(0)=0, $
|
(5.37) |
$ H1(y)<0,H′1(y)<0,∀y∈(0,+∞), $
|
(5.38) |
$ H2(y)<0,∀y∈(−Δ2,+∞). $
|
(5.39) |
Concerning $ H_{2} $ observe that we have the following key symmetry
$ H2((α2)2(−Δ2−y)(α2)2+(γ2)2y)=H2(y),∀y∈[0,−Δ2]. $
|
(5.40) |
This implies, by recalling the definition of $ \Gamma^2 $ in Eq (5.32),
$ sup{−H2(y)H1(y):y∈(0,Γ2]}=sup{−H2(y)H1(y);y∈(0,+∞)}, $
|
(5.41) |
or equivalently $ N_{0} = K $. Considering the fact $ H_{2}/H_{1} $ is a continuous function, in order to prove that $ N_{0} $ is bounded it only remains to show that there exists a finite limit as $ y $ tends to $ 0^{+} $. From Eq (5.35),
$ limy→0+−H2(y)H1(y)=−H′2(0)H′1(0)=−Δ2(α1)2Δ1(α2)2∈(0,+∞), $
|
(5.42) |
which concludes that $ N_0 < +\infty $. As $ n_{1} = \mu n_{2} $, Eq (5.34) is equivalent to $ \mu > N_{0} $.
To show such a choice of $ \tau_{0} $ is optimal, is suffices to observe that if $ \mu < N_{0} $ (or equivalently $ n_1/(n_1+ n_2) < \tau_0 $) then by continuity there exists a subset $ [y_{1}, y_{2}]\subset(0, \Gamma^{2}] $ with $ y_{1}\neq y_{2} $ such that
$ μ<−H1(y)H2(y), for any y∈[y1,y2], $
|
(5.43) |
which implies that for any $ y\in [y_{1}, y_{2}] $, $ y > 0 $ and
$ ((α1)2+(γ1)2y(α1)2+((β1)2−2α1)y+y2)μ((α2)2+(γ2)2y(α2)2+((β2)2−2α2)y+y2)>1. $
|
(5.44) |
Setting $ x = \sqrt{y} $, Eq (5.44) together with Lemma 5.1 and Eq (5.5) allows us to conclude that there exists $ M $ large enough such that for any $ n_{1} > M $ and $ n_{2} > M $ satisfying $ n_1/(n_1+n_2) = \tau $, the system (4.16) is unstable around the equilibrium flow $ (\bar{h}_{i}, \bar v)_{i\in\{1, 2\}} $.
Proof of Corollary 5.3. We have now proved the existence of the critical penetration rate $ \tau_{0} $ and we obtained a quantitative characterization. Note that $ \tau_{0} $ can be easily computed by a minimization algorithm using Eqs (5.31)–(5.33) and provides the optimal penetration rate of stable cars to stabilize the traffic. In the following, for the qualitative study convenience, we also present some lower and upper bounds of $ \tau_0 $ (or equivalently, some lower and upper bounds of $ N_0 $).
A simple lower bound of $ \tau_0 $.
Thanks to Eq (5.42), we see that
$ N0≥−Δ2(α1)2Δ1(α2)2, $
|
(5.45) |
hence a lower bound of $ \tau_0 $ can be expressed by
$ Bl:=−Δ2(α1)2Δ1(α2)2−Δ2(α1)2. $
|
(5.46) |
Recall that $ -\Delta_{2} $ is strictly greater than 0, so $ B_{l}\in(0, 1) $.
A simple upper bound of $ \tau_0 $.
The precise value of $ N_0 $ is given by Eq (5.31), but it is rather difficult to determine by hand. Indeed, to do so we would need to compare the extreme points of $ - H_2/ H_1 $, which are given by
$ E:={z∈(0,Γ2):H′2(z)H1(z)−H2(z)H′1(z)=0}. $
|
(5.47) |
In terms of those extreme points, $ N_0 $ is further given by
$ N0=max{−H2(y)H1(y):y∈E∪{0,Γ2}}, $
|
(5.48) |
where
$ −H2(0)H1(0):=−H′2(0)H′1(0)=−Δ2(α11)2Δ1(α21)2. $
|
(5.49) |
Actually, we can exclude $ \Gamma^2 $ in $ \{0, \Gamma^2\} $ from the expression (5.48): suppose that $ \Gamma^2\in \mathcal{E} $ then $ \mathcal{E}\cup \{0, \Gamma^2\} = \mathcal{E}\cup \{0\} $; otherwise, there exists some $ y_0\in (\Gamma^2-\delta, \Gamma^2+ \delta) $ such that $ -H_2(y_0)/H_1(y_0) $ is bigger than those of $ \Gamma^2 $, then thanks to the symmetric property of $ H_2 $ and the monotonous property fo $ H_1 $, Eqs (5.37)–(5.40), there exists some $ y_1\in (\Gamma^2-\delta, \Gamma^2) $ such that $ -H_2(y_1)/H_1(y_1) $ is bigger than those of $ \Gamma^2 $. Therefore,
$ N0=max{−H2(y)H1(y):y∈E∪{0}}. $
|
(5.50) |
Next, notice that for any extreme point $ z\in \mathcal{E} $ we have
$ −H2(z)H1(z)=−H′2(z)H′1(z), $
|
(5.51) |
which, to be combined with Eq (5.49), yields
$ N0=max{−H′2(y)H′1(y):y∈E∪{0}}. $
|
(5.52) |
We remark here that, though it is relatively easy to get the maximum value of $ -H'_2/H'_1 $ in $ [0, \Gamma^2] $ (as its extreme points can be calculated explicitly via polynomials), this value is not necessarily equivalent to $ N_0 $. More precisely,
$ N0≤max{−H′2(y)H′1(y):y∈[0,Γ2]}=:˜N0. $
|
(5.53) |
The value of $ \widetilde N_0 $ can also be expressed in terms of extreme points:
$ ˜N0=max{−H′2(y)H′1(y):y∈˜E}, $
|
(5.54) |
$ ˜E:={z∈[0,Γ2]:H″2(z)H′1(z)−H′2(z)H″1(z)=0}. $
|
(5.55) |
In the following we present a simple upper bound for $ \widetilde N_0 $. Observe that $ - H'_2(y)/ H'_1(y) $ is characterized by
$ −(γ2)2y2−2(α2)2y−(α2)2Δ2((α2)2+(γ2)2y)((α2)2+((β2)2−2α2)y+y2)⋅((α1)2+(γ1)2y)((α1)2+((β1)2−2α1)y+y2)(γ1)2y2+2(α1)2y+(α1)2Δ1, $
|
(5.56) |
while for $ y\in [0, \Gamma^2] $ there are
$ (γ1)2y2+2(α1)2y+(α1)2Δ1≥(α1)2Δ1, $
|
(5.57) |
$ ((α1)2+(γ1)2y)((α1)2+((β1)2−2α1)y+y2)≤((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2), $
|
(5.58) |
$ −(γ2)2y2−2(α2)2y−(α2)2Δ2≤−(α2)2Δ2, $
|
(5.59) |
$ ((α2)2+(γ2)2y)((α2)2+((β2)2−2α2)y+y2)≥(α2)2(β2)2. $
|
(5.60) |
Consequently,
$ ˜N0≤(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2)(β2)2(α1)2Δ1, $
|
(5.61) |
which further implies the following upper bound of $ \tau_0 $:
$ Bu:=(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2)(β2)2(α1)2Δ1+(−Δ2)((α1)2+(γ1)2Γ2)((α1)2+((β1)2−2α1)Γ2+(Γ2)2). $
|
(5.62) |
Let us now prove point (2) of Theorem 5.2. We define again $ \mu = n_{1}/n_{2} $. Suppose that $ \Delta_{1} = 0 $ and $ \Delta_{2} < 0 $. Thanks to Eq (5.35), we know that
$ μH′1(0)+H′2(0)=−Δ2(α2)2>0, $
|
(5.63) |
which, to be combined with the fact that $ \mu H_1(0)+ H_2(0) = 0 $, imply the existence of $ y > 0 $ such that
$ μH1(y)+H2(y)>0. $
|
(5.64) |
As a direct consequence, the system is not stable.
In this section, analogy to the preceding Section, we study the stability of the equilibrium flows in a multi-phase mixed traffic flow. For any given equilibrium velocity $ \bar{v} $, keeping the notation $ \Lambda_j $ and Eq (4.21), we have that
$ Λj=(αj,βj,γj)∈{Λ1,...,Λm}={(α1,β1,γ1),...,(αm,βm,γm)}. $
|
(6.1) |
Again, for $ k\in \{1, 2, ..., m\} $ we denote by $ n_{k} $ the number of Type k vehicles with parameters $ \Lambda^{k} $ such that the total number of vehicles is $ n = \sum_{j = 1}^m n_j $.
For any ordering of the vehicle on road, i.e., the $ j $-th vehicle is of Type $ a_j $, the linearized system around the unique equilibrium flow is
$ {˙yj=uj+1−uj,˙uj=αjyj−βjuj+γjuj+1,(αj,βj,γj)=(αaj,βaj,γaj),(yj(t),uj(t))j∈{1,...,n}∈H. $
|
(6.2) |
Similar to Theorem 5.1 and Theorem 5.2 we have the following theorem concerning the stability of the $ m $-phase mixed ring road traffic.
Theorem 6.1. Let given $ \bar v > 0 $. We assume without loss of generality that $ {\Delta_{1}}\geq {\Delta_{2}}\geq...\geq \Delta^m $.
(1) If $ \Delta^m\geq 0 $, then for any $ (n_1, n_2, n_3, ..., n_m)\in \mathbb{N}^m $ and any ordering of the vehicles on the road, the ring road traffic system (4.16) is locally exponentially stable around the equilibrium flow.
(2) If $ {\Delta_{1}} > 0 $ and $ \Delta^m < 0 $, then there exists some effectively computable threshold constant $ \tau_1\in (0, 1) $ depending on $ \{\Lambda^{1}, ..., \Lambda^{m}\} $ such that for any m-tuple $ (n_1, n_2, ..., n_m)\in \mathbb{N}^2 $ verifying
$ n1∑mk=1nk>τ1, $
|
(6.3) |
in other words, for any ordering of the vehicles on the road $ (a_1, a_2, ..., a_n)\in \mathcal{K} $, the ring road traffic system (4.16) is locally exponentially stable around the equilibrium flow associated to $ \bar{v} $.
On the other hand, we further assume that for some $ p\in\{1, ..., m-1\} $ there is $ \Delta^p > 0\geq \Delta^{p+1} $. Then, there exists some $ \tau_2\in (0, 1) $ effectively computable such that for any fixed penetration rates
$ (n1∑mj=1nj,...,nm∑mj=1nj) satisfying ∑Pj=1nj∑mj=1nj<τ2 $
|
(6.4) |
there exists $ M > 0 $ such that for any $ {n_1, ..., n_m}\in \mathbb{N}^m $ satisfying $ \sum_{k = 1}^m n_k > M $, the ring road traffic system (4.16) is unstable around the equilibrium flow associated to $ \bar v $.
(3) If $ {\Delta_{1}} = 0 $ (namely, $ \Lambda^{1} $ is critical), then for any fixed penetration rate of the cars the ring road traffic system (4.16) is unstable around the equilibrium flow associated to $ \bar{v} $ provided sufficiently many cars on road.
Indeed, similar to the two-phase traffic, we look at the ring road traffic with $ m $ populations having respectively $ (n_1, n_2, ..., n_m) $ many vehicles. The local exponential stability of this system is still equivalent to the study of the real part of eigenvalues of the linearized system that are explicitly given by the roots of the polynomial (4.35) presented in Section 4.3:
$ n∏j=1Fj(ω)=1,Fj(ω):=γjω+αjω2+βjω+αj. $
|
(6.5) |
It becomes more delicate to calculate the exact roots of the polynomial. Instead, again, we investigate the value of the polynomial on the imaginary axis. The question becomes:
$ whether m∏k=1((αk)2+(γk)2x2(αk)2+((βk)2−2αk)x2+x4)nk<1 holds for any x∈R∖{0}? $
|
(6.6) |
Thanks to the same reasoning as in Lemma 5.1, we get the following conclusion.
● If Condition (6.6) is satisfied, then all the roots of the polynomial excluding 0 are distributed on the left complex region i.e., $ \{z\in \mathbb{C}: \Re{(z)} < 0\} $. Hence, the linearized system is exponentially stable.
● If Condition (6.6) is not verified, and if further there exists some $ x_0\in \mathbb{R} $ such that the value of the function in (6.6) is strictly larger than 1, then the ring road traffic system with these penetration rates of vehicles is not stable provided sufficiently many cars. This is obtained using the same reasoning as in the two population case (see Eqs (5.9)–(5.19)).
Recalling Definitions 4.2, 4.3 concerning the classification of $ \Delta^k $, we know that
● if $ \Delta^k > 0 $, then
$ ((αk)2+(γk)2x2(αk)2+((βk)2−2αk)x2+x4)nk<1,∀x∈R∖{0}; $
|
(6.7) |
● if $ \Delta^k = 0 $, then Inequality (6.7) is also satisfied;
● if $ \Delta^k < 0 $, then
$ ((αk)2+(γk)2x2(αk)2+((βk)2−2αk)x2+x4)nk>1, for some x∈R∖{0}. $
|
(6.8) |
This observation, together with Condition (6.6), finally leads to Theorem 6.1 in a way that is very similar to the proofs of Theorem 5.1, 5.2. The easier cases (1) and (3) are direct consequences of the preceding observations. Concerning the mixed case (2) such that both stable and unstable vehicles coexist, from a heuristic point of view the more stable cars there are on the road, the more likely it is that Condition (6.6) is satisfied, i.e., that for any ordering of the vehicles on the road $ (a_1, a_2, ..., a_n)\in \mathcal{K} $, the ring road traffic system (4.16) is locally exponentially stable around the equilibrium flow associated to $ \bar{v} $. Otherwise with fewer stable cars on road the unstable cars will dominate the traffic to prevent us from getting Condition (6.6). This comes from the fact that we get the following condition for the stability instead of Eq (5.28)
$ n∑i=1niHi(y)<0,∀y≥0, $
|
(6.9) |
where $ H_{i}(y) = \log\left(\frac{(\alpha^{i})^2+ (\gamma^{i})^2 y}{(\alpha^{i})^2+ ((\beta^{i})^2- 2\alpha^{i}) y+ y^2}\right) $ is defined similarly as in Eq (5.27); in this specific case there is
$ Hi(0)=0 and H′i(0)<0,∀i∈{1,2,...,p}, $
|
(6.10) |
$ Hi(0)=0 and H′i(0)>0,∀i∈{p+1,p+2,...,m}. $
|
(6.11) |
In this Section we present some numerical experiments to illustrate Theorems 3.2–3.4 with two type of populations described by some set of parameters $ \Lambda_{1} $ and $ \Lambda_{2} $.
We assume for this example that the vehicles are described by the nonlinear Bando-FTL model (2.5) and we consider two populations described by the parameters $ (a_{1}, b_{1}) $ and $ (a_{2}, b_{2}) $. We assume that they have the same velocity preference $ V $ given by
$ V(h)=Vmaxtanh(h−lvd0−2)+tanh(2)1+tanh(2),h∈(0,+∞), $
|
(7.1) |
which is a usual choice for Bando-FTL [2]. The quantity $ V_{\max} $ corresponds to the maximal velocity preference for a given car. It corresponds to the inner speed limit of a given driver with a given car (which can be different from the official speed limit). $ l_{v} $ is the length of the vehicle, while $ 2d_{0} $ is a safety distance below which the feeling of unsafeness would make the driver wants to stop. Here, their values are respectively taken as: $ V_{\max} = 9.25\; m.s^{-1} $, $ l_{v} = 4.5\; m $ and $ d_{0} = 2.5 $. In all the numerical simulations the cars are initially placed at an equal distance on the road with a random ordering between the two types of drivers. The initial velocity is close to half the velocity preference with a small random perturbation (between 0 and $ 0.3\; m.s^{-1} $).
We study a case where the first population has a very stable behavior on the road with $ \Delta_{1} > 0 $ while the second population is made of slightly more aggressive driver that have a slightly unstable behavior such that $ \Delta_{2} < 0 $ but $ |\Delta_{2}| < |\Delta_{1}| $. To do so, we choose the parameters $ a_{1} = 4 $, $ b_{1} = 20 $, $ a_{2} = 0.5 $, $ b_{2} = 20 $, thus the instability of the second population will simply come from a higher sensitivity to the velocity preference rather than the "Follow-the-leader" behavior. We choose $ L $ and $ N $ such that $ L/N = 10.4\; m $, which is a steady-state value similar to the setting of the real life experiment described in [17]. Note that parameters $ (a_{2}, b_{2}) $ are typical values that were obtained in [16] after calibration on data from real life experiments. From Theorem 5.2 we deduce the following:
$ Δ1=7.28,Δ2=−0.84,τ0=0.881 $
|
(7.2) |
as we can see, although $ |\Delta_{1}| $ is one order of magnitude above $ |\Delta_{2}| $, the penetration rate of stable vehicle needed for having a stable flow is very high and above $ 88\% $. This means that even a very low proportion of slightly more aggressive drivers in a large road can completely destabilize a traffic that would be very stable otherwise. In Figure 1 (left) we show the instantaneous speed variance across drivers of a traffic flow with 500 cars, a penetration rate of stable cars of $ 80\% $ ($ \tau = 0.802 $), and we see that the system becomes quickly unstable as the speed variance only increase during the entire simulation. In Figure 1 (right) we show the instantaneous speed variance of the same traffic when the penetration rate of stable car is $ 88.2\% $ instead and we see that the speed variance, already low at initial time, decreases exponentially fast. Here we define the instantaneous speed variance across drivers as the variance of the drivers' velocity at a given time $ t $, i.e., the quantity
$ Vs(t)=1NN∑i=1(vi(t)−1NN∑i=jvj(t))2. $
|
(7.3) |
However, as highlighted in Theorem 3.3, this instability might only occur with a large enough number of cars and might be missed when looking at experiments with a small number of cars such as [17,18] rather than a real freeway with sufficiently many cars. To get some insight on this phenomenon, we ran several numerical experiments where we kept the same steady state described by $ L/N = 10.4\; $m, hence the same $ \tau_{0} = 0.881 $, but we varied the number of cars $ N $ from 10 to 120. Moreover, for each number of cars $ N $ we also ran several simulations with different proportions of cars exhibiting a stable behavior. In Figure 2 we show for each fixed number of cars $ N $, the minimal penetration rate of cars with stable behavior above which we observe a terminal instantaneous speed variance of the system (i.e., the variance of speeds taken across all cars of the system at the given final time) below 0.01 $ m^2.s^{-2} $ after 2000s of simulations. Such a negligible terminal speed variance indicates the strong stability of the whole system while a non-negligible speed variance indicates that the cars are not in the uniform flow equilibrium at final time. We see that the effective penetration rate of stable cars above which the total system is stable is significantly lower than $ \tau_{0} $ for small number of cars, and this value becomes quickly close to the theoretical value $ \tau_{0} = 0.881 $ above $ N = 40 $ (to be compared with the values $ N\sim 20 $ of [17,18]).
The research is based upon work supported by the U.S. Department of Energy's Office of Energy Efficiency and Renewable Energy (EERE) under the Vehicle Technologies Office award number CID DE-EE0008872. The views expressed herein do not necessarily represent the views of the U.S. Department of Energy or the United States Government. The authors acknowledge the Office of Advanced Research Computing (OARC) at Rutgers, The State University of New Jersey for providing access to the Amarel cluster and associated research computing resources that have contributed to the results reported here. The authors acknowledge the support of the USDOT award # 69A3551747124 via the NYU-C2SMART project "Collaborative Driving, Ramp Metering and Mean-field Controls". The authors would also like to thank the IEA project SHYSTRA (CNRS) for their support.
The authors declare there is no conflict of interest.
[1] | Appl. Math. Lett., 13 (2000), 91-95. |
[2] | J. Biol. Sys., 15 (2007), 1-19. |
[3] | Springer-Verlag, New York, 1995. |
[4] | Math. Biosci., 191 (2004), 159-184. |
[5] | Math. Med. Biol., 26 (2009), 63-95. |
[6] | Math. Biosci., 222 (2009), 13-26. |
[7] | Accepted for Math. Pop. Studies. |
[8] | Philos. Trans. R. Soc. London, 115 (1825), 513-585. |
[9] | Cancer Res., 59 (1999), 4770-4775. |
[10] | Springer, New York, 1993. |
[11] | Ann. N. Y. Acad. Sci., 50 (1948), 221-246. |
[12] | SIAM J. Control Optim., 46 (2007), 1052-1079. |
[13] | J. Theor. Biol., 252 (2008), 295-312. |
[14] | Springer, Berlin-Heidelberg, 2007. |
[15] | (submitted). |
[16] | J. Math. Anal. Appl., 382 (2011), 180-203. |
[17] | Math. and Comp. Modelling, 54 (2011), 2183-2198. |
[18] | Math. Biosci. Eng., 8 (2011), 591-603. |
[19] | in "Mathematical Population Dynamics" (eds. O. Arino, D. Axelrod and M. Kimmel), Wuertz, Winnipeg, Canada, (1995), 335-348. |