Properties of equilibria in transport problems with complex interactions between users

It is well known that uniqueness and stability are guaranteed properties of traﬃc equi- libria in static user-equilibrium traﬃc assignment problems, if the link travel utilities are assumed to be strictly monotonically decreasing with respect to the link traﬃc volumes. However, these preferable properties may not necessarily hold in a wide range of transport problems with complex interactions, e.g. asymmetric interactions (including dynamic traﬃc assignment), social interactions, or with economies of scale. This study aims to investigate such solution properties of transport models with complex interactions between users. Generic formulations of models are considered in this study, both for utility functions and for the evolutionary dynamics relevant to the stability analysis. Such an anal- ysis for a generic formulation is mathematically challenging due to the potential non-differentiability of the dynamical system, precluding the application of standard analyses for smooth systems. To address this issue, this study proposes a transport system with two alternatives and two user groups. While it is a simple model whose dynamics can be depicted on a plane, it also includes the core components of transport models, i.e. multiple choices and user-classes. This study classiﬁes all possible formulations into nine cases with respect to the signs (i.e. positive or negative) of interactions between users. Then, the evolutionary dynamics of each case is mathematically analysed to examine stability of equilibria. Finally, the solution properties of each case is revealed. Multiple equilibria exist in many cases. In addition, cases with no stable equilibrium are also found, yet even in such cases we are able to characterise the circumstances in which the different kinds of unstable behaviour may arise.


Introduction
It is well known that uniqueness and stability are general properties of traffic equilibria in static user-equilibrium traffic assignment problems, in which travel utility of each link is strictly monotone decreasing with respect to its traffic volume. Beckmann et al. (1956) have shown that a unique equilibrium solution exists in any case satisfying these properties. Smith (1984a) established stability (more specifically, asymptotic stability) of the equilibrium solution by a certain form of evolutionary dynamics, which is now called Smith dynamics ( Sandholm, 2010a ). These findings in the area of transport studies have led to the study of such phenomena in wider (non-transport) contexts in the field of evolutionary game theory. In this field, a game that is identical to the conventional traffic assignment model considered by Beckmann et al. (1956) is called a congestion game. Particularly, Sandholm (2001) has shown that any potential game (including any congestion game) asymptotically converges to an equilibrium solution by any evolutionary dynamics satisfying mild conditions.
The preferable properties such as uniqueness and stability of a standard traffic assignment problem, as shown above, are actually not a general property of transport problems in a wider sense. For example, the existence of multiple equilibria is well known in the case of a model with asymmetric interactions between links ( Watling, 1996 ). The asymmetry of travel cost functions has been also mentioned by Heydecker (1983) and Smith (1984b) . Recently, studies on dynamic traffic assignment (DTA) problems have shown a case in which multiple equilibria exists ( Iryo, 2011 ) or a unique equilibrium solution exists that is unstable ( Iryo, 2008;2019;Guo et al., 2018 ). These counterexamples (to the assumption of uniqueness and stability) arise because the formulation of congestion is more complex than that of conventional traffic assignment problems, especially the strong asymmetry in the utility functions. In such cases, transport analysts are no longer able to rely on the equilibrium solution with no explicit examination of its uniqueness and stability.
In addition to the complexity in formulating congestion, one may need to consider another reason that can deteriorate the nice properties of an equilibrium solution in a transport problem, that is, the existence of positive interactions. In the traditional context of transport studies, interactions between users are considered to be mostly originated by congestion phenomena, which creates negative interactions (i.e. increase of the number of users in a transport facility decreases its attractiveness for other users). In contrast, positive interaction is characterised by a situation in which the increased patronage of a transport facility increases its attractiveness for other users. There are several potential sources of positive interaction (though we are not claiming that the sources are necessarily positive in all cases, just that positive interactions may exist). Social interaction has been extensively studied as a source of positive interactions in transport studies (see e.g. Cirillo, 2016 andManess et al., 2015 andreferences therein). Carrasco et al. (2008) propose a particular theoretical frame to represent the social dimension in social activity-travel behaviour, representing the fact that 'a key cause of travel behaviour is the social dimension represented by social networks'. Brock and Durlauf (2001) explain a form of positive social interaction as being "the desire of individuals to conform to the behaviour of others in an environment of non-cooperative decision-making'. They indicate the existence of multiple equilibria in the presence of social interaction effects, a phenomenon that was examined empirically by Fukuda and Morichi (2007) . Goetzke and Rave (2011) found aggregate bicycle modal share to be a significant factor in explaining differences in individuals" probabilities to choose a bicycle mode. A different kind of positive interaction was identified by Jacobsen (2003) , who showed empirically that the collision rates between motorists and walkers/cyclists declined with increases in the numbers of people walking/cycling. A qualitatively similar effect was reflected also in studies of car-drivers ( McCarthy, 2001 ) and pedestrians ( Leden, 2002 ), and extends by implication to perceptions of safety associated with the level of other activity in an area ( Cho et al., 2009 ). Economies of scale are also a source of positive interactions. The Mohring effect ( Mohring, 1972 ) is the well-known leading example, which observes that travel cost of a public transport service decreases when the increase of its demand increases the frequency of service. In a traffic assignment problem, Ying and Yang (2005) have modelled transit lines exhibiting economies of scale and showed the existence of multiple equilibria in a numerical example. Positive interactions can also arise from demand responsive control schemes, i.e. the more people choose an alternative, the more resources are allocated to it, such as green times in a responsive traffic signal system ( Watling, 1996;Smith and Mounce, 2011 ) or bus frequencies in a demand responsive transport system ( Li and Yang, 2016;Cantarella et al., 2015 ). These authors have shown the potential existence of multiple equilibria or unstable solutions in such cases.
The aim of this study is to understand solution properties, such as uniqueness and stability, of transport models including complex forms of interactions between users. We especially try to classify a wide range of generic model specifications with respect to the solution properties. It is worthwhile to conduct a comprehensive analysis of generic models for the following reasons. First, we can identify all possible properties of equilibrium solutions in models with various types of complex interactions, such as asymmetric interactions and a mixture of positive and negative interactions. Thus far, only little is known for solution properties of such models in generic forms. Second, knowing the relationship of a model specification and properties of its equilibria comprehensively is useful to predict solution properties of a particular transport problem. With such a prediction in hand, we can more easily decide the best perspective for the analysis, e.g. whether finding equilibria is sufficient to understand model behaviour or a dynamical analysis is necessary. Lastly, checking existence of a stable equilibrium solution in generic forms of interactions will reveal whether known examples of non-existence of stable equilibrium, such as in Iryo (2008) , are very special or not. It also clarifies how far the concept of equilibrium is applicable in generic transport problems.
Investigating the dynamical adjustment process is vital to clarify which equilibrium solution is realised (or whether the system certainly reaches an equilibrium solution or not) in a transport system with the complex interactions. In the present study, we particularly use a continuous version of the dynamical adjustment processes, in which a day is described by a continuous number. A number of studies such as Smith (1984a) , Friesz et al. (1994) , Zhang and Nagurney (1996) , Mounce (2006) , Iryo (2008) , Yang and Zhang (2009) , He et al. (2010) , Han and Du (2012) , Smith and Watling (2016) , Iryo (2016) have employed this version of the dynamics. In the remaining part of the present paper, we refer to this type of dynamics as evolutionary dynamics . Several types of evolutionary dynamics have been proposed in existing studies, such as the pairwise adjustment process (or Smith dynamics; Smith (1984a) ), network tâtonnement process ( Friesz et al., 1994 ), best response dynamics, perturbed best response dynamics, replicator dynamics, and so on (see Sandholm, 2010b for a summary of them).
Although the continuous dynamical adjustment process is mathematically tractable, understanding general properties of models including general forms of interactions is a methodological challenge. The primal reason is the non-smooth property of evolutionary dynamics. Every evolutionary dynamics can be formalised as an ordinary differential equation (ODE) d x / d t = f (x ) , which is called an autonomous system (i.e. d x /d t depends on x only). The non-smoothness implies that f ( x ) is not differentiable at a certain point. A number of textbooks such as Hirsch et al. (2004) introduce how to mathematically analyse such systems, whereas they basically assume that function f ( x ) is at least first-order differentiable at all times. Particularly, calculating eigenvalues of its Jacobian matrix at a stationary point is the well-known standard technique to analyse its stability. On the other hand, most existing evolutionary dynamics have the non-smooth property around an equilibrium solution. It is caused by an important property of most evolutionary dynamics, i.e. drivers move from a costly option to a cheaper one in their day-to-day adjustment processes, while they never move in the other direction (i.e. from a cheaper one to a costly one). It is not only a natural representation of users' day-to-day adjustment behaviour but also vital to ensure that any equilibrium solution is a stationary point of the dynamics. Non-smooth autonomous dynamics have been extensively studied in mechanical engineering to analyse bouncing motions of objects. However, as far as the authors know, there exists no valid methodology for analysing stability property in general non-smooth autonomous dynamics. A textbook by Bernardo et al. (2008) (page 96) actually states that 'general progress toward understanding and classifying the dynamics of piecewise-smooth systems using such methods would appear hopeless' after showing a few examples for two-and three-dimensional problem and an approach by Lyapunov function, which has been also used by Smith (1984a) in transport studies. However, this approach is not applicable for general problems because 'finding such functions in practice is at best difficult' ( Bernardo et al., 2008 ).
To avoid this methodological barrier, we consider a simple model of a transport system consisting of two alternatives and two user groups in the present study. The proposed model can be represented by a planar dynamical system, in which classifying the dynamics of piecewise-smooth systems is still feasible. In addition, this model can include heterogeneity of users in the simplest form (i.e. two groups). Many practical transport problems naturally include heterogeneity of users, e.g. multiple origin-destination pairs, different values of time, status of car ownership (i.e. whether users own a car or not), asymmetry of between-group interactions, and so on.
The structure of the present paper is as follows. Background and motivations of the study were stated in Section 1 . The specification of the transport system and the definition of equilibrium and evolutionary dynamics are explained in Section 2 . In Section 3 , the mathematical definition of stability is shown. Then, in Section 4 , we propose a classification of the proposed transport system after a few conditions on the transport system are introduced. All possible transport systems are classified into nine cases (Cases 1 -9). Then, we analyse stability of the evolutionary dynamics for each case in Section 5 . In this section, we introduce the graphical representations of the dynamics for the analysis, which will be finally summarised in one table called the summary sheet. Section 6 proposes key examples of all nine cases and a numerical example. Section 7 concludes the study and states future directions.

Transport system
We consider a transport system called two-alternatives-two-groups transport system , which consists of two user groups A and B with a binary choice, denoted by two alternatives , e.g. different transport modes. The alternatives are denoted by 0 and 1. The number of users selecting alternative k ∈ {0, 1} for group i ∈ { A , B } is denoted by x k i . The number of users in each group is constant at any time and so without loss of generality is set to 1, i.e. x 0 A + x 1 A = 1 and x 0 B + x 1 B = 1 . We reduce the number of the variables from four (i.e. x 0 to two to simplify the analyses by introducing two variables x A and x B satisfying 0 ≤ x A ≤ 1 and 0 ≤ x B ≤ 1, respectively. The four original values are then described as follows: Variables x A and x B determine the state of the system. Note that variables x A and x 0 A are mathematically not identical because the former is correlated with x 1 . We also use the vector form of their variables, i.e., x = (x A , x B ) , referred to as the state vector . It is plotted on the two-dimensional graph called the phase diagram , in which the horizontal axis and vertical axis indicate x A and x B , respectively. An example of the phase diagram is depicted in Fig. 1 . The feasible region of x is denoted by S is divided into two parts. One is the peripheral boundary , which is denoted by the set ∂S defined by The remaining part is called the interior region , which is denoted by the set S − ∂S.

Utility function
Users' preferences of the alternatives are described by utility functions. They are denoted by u k Apart from conditions to be introduced later, we consider generic forms for these functions so that we can describe potentially complex interactions between users. Users prefer an alternative whose utility is higher than that of the other alternative and their behaviour depends only on the difference in utility between two alternatives. This difference (denoted by u i ( x )) is denoted as follows: Users in group i prefer alternative 0 when u i ( x ) > 0 and alternative 1 when u i ( x ) < 0. They are indifferent to which alternative ) contains the utility difference function for both groups. If u i (x ) = 0 is a single connected curve on the phase plane, it divides the entire feasible region into two subregions, i.e. a subregion in which users prefer alternative 0 to 1 and another in which users prefer alternative 1 to 0. On the curve u i (x ) = 0 , these alternatives are equally preferable for all users in group i . We call the curve u i (x ) = 0 the interior boundary of group i . The following two conditions are introduced as assumptions on u ( x ): [U1] u i ( x ) is Lipschitz continuous, first-order differentiable with respect to both x A and x B , and bounded.
[U2] (a) u (x These conditions are not so restrictive.
[U1] is a standard mathematical condition regarding continuity, differentiability, and boundedness.
[U2] is introduced to exclude special cases that are not so relevant for analysis, namely situations such as: (a) u (x ) = 0 and x ∈ ∂S , (b) u A (x ) = 0 or u B (x ) = 0 at a corner of the feasible region.
These are situations that are unlikely to hold in practice, which on the other hand complicate the analysis considerably. Actually, even if they arose, we could eliminate them by adding small perturbations onto the utility function. By introducing [U2], we can avoid complicated classifications in the subsequent sections.

Equilibrium
We consider Nash equilibrium as a candidate of a state that we are likely to observe in the real world. In Nash equilibrium, all users select the most preferred alternative(s) and consequently have no incentive to change their choices. That no user changes their choices implies that the state of the transport system stays at a stationary point.
A state vector x * is called an equilibrium solution if and only if it satisfies Note that u i indicates the difference in utility between two alternatives. There may exist two or more equilibrium solutions in the transport system. The set of x * satisfying Eq. (5) (i.e. the set including all equilibrium solutions) is denoted by S * . For the convenience of the subsequent analyses, we classify equilibrium solutions into two categories with respect to the value of u ( x ) in equilibrium as follows: 1. (Interior equilibrium) u (x * ) = 0 , which is identical to x * ∈ S * ∩ (S − ∂S) 2. (Peripheral equilibrium) u ( x * ) = 0 , which is identical to x * ∈ S * ∩ ∂S.
Note that the case in which both u (x * ) = 0 and x * ∈ S * ∩ ∂S are simultaneously satisfied has been already eliminated by condition [U2].

Evolutionary dynamics
To investigate stability of equilibrium, we have to analyse how the state vector moves when it is not at an equilibrium solution, and hence need to specify the model describing the motion of the state vector caused by changes of users' choice behaviour over days. It is called the evolutionary dynamics . We employ the following generic form: where ˙ x i = d x i / d t and t represents a continuous day, for the definition of the evolutionary dynamics. This formulation was proposed by Sandholm (2010a,b) . The pair of functions ρ + i (x , u i ) and ρ − i (x , u i ) , called the revision protocol , determines the dynamics of the evolution of users' choice behaviour, indicates how fast users change their alternative from 0 to 1. Note that day-to-day time is described by a continuous number in this model, although a day might be more naturally represented by a discrete number in the real world. See Sandholm (2003) and Hofbauer and Sandholm (2007) to understand that the evolutionary dynamics can be derived from the natural definition of the system comprising of a finite number of users and a discrete day variable. x is actually a function of t and hence should be described as x ( t ), while we omit t in the following descriptions for simplicity unless it needs to be explicitly shown.
We now introduce a few conditions for the revision protocol so as to keep it compatible with the definition of Nash equilibrium, in addition to a few other basic mathematical conditions. They are: , and u i , as well as bounded.
They are also partially differentiable with respect to u i > 0 and u i < 0, respectively.
[D3] The sign of u i is consistent with the values of ρ [D4] ρ + i (x , u i ) is directionally partially differentiable with respect to u i when u i → +0 and its derivative is greater than zero. Similarly, ρ − i (x , u i ) is directionally partially differentiable with respect to u i when u i → −0 and its derivative is greater than zero.
[D1] and [D2] are basic mathematical conditions. Regarding [D3], the third equation in (7) is included in order to maintain the Nash stationary property (NS), defined by which implies that any Nash equilibrium solution is a stationary solution of the evolutionary dynamics and vice versa. Other equations in (7) imply the positive correlation property (PC) ( Sandholm, 2010a;2010b ), defined by which is identical to owing to the relationships The positive correlation property in the present transport system implies that the users always switch their alternative from a less preferable one to a more preferable one. It is a natural behavioural assumption and also consistent with the concept of Nash equilibrium, which implicitly assumes that users make their decisions so as to maximise their own utility. Proof. See Appendix.
[D4] expands the differentiability of ρ + i (x , u i ) and ρ − i (x , u i ) to the neighbourhood of u i = 0 . It also guarantees that these slope in the neighbourhood of u i = 0 is not zero. We will trace the trajectory of the evolutionary dynamics for t ∈ [0, ∞ ), with x (0) given. x (0) is called the initial state vector. We will set the following mild condition to maintain simplicity in the subsequent analyses: The following proposition also holds (see Appendix for its proof), implying that we do not have to analyse anywhere not in S + whenever condition [I] holds.
Finally, we introduce the following theorem that guarantees the existence of a unique trajectory at all times: Theorem 2. (Existence of unique trajectory) Given x (0), the differential equation in (6) always has a unique solution of Proof. See Sandholm (2010a,b) . Note that the proof in this reference requires that ˙ x i must be described by a function that is locally Lipschitz continuous and bounded.
[U1] and [D1] are used to guarantee these properties.

Definition of stability and classification of instability
In this study, we consider that the term stability implies asymptotic stability (see e.g. Hirsch et al., 2004 , p. 175). An equilibrium solution x * is asymptotically stable if and only if the following two properties hold: The convergence property implies that x ( t ) converges to x * from anywhere in its neighbourhood. Lyapunov stability implies that We further classify any unstable solution into three categories: 1. Saddle equilibrium solution , which is not stable but satisfies the following property: 2. Diverging equilibrium solution , which is not a stable nor a saddle equilibrium solution.
3. Circular equilibrium solution , which is stable but does not satisfy the convergence property.
A saddle equilibrium solution is reachable from a part of its neighbourhood. A diverging equilibrium solution is completely unreachable and any state vector cannot remain in the neighbourhood of the equilibrium. A circular equilibrium solution is not reachable but the state vector can remain in the neighbourhood of the equilibrium. The classification of stable, saddle, diverging, and circular actually corresponds to one that is often introduced for a linear planar dynamical system in standard textbooks of dynamical systems (e.g. Hirsch et al., 2004 , p. 39-), namely 'sink / spiral sink', 'saddle', 'source / spiral source', and 'centre'. Two eigenvalues of the 2 × 2 matrix determine which class is applicable, e.g. two negative real parts implies stability, etc. However, we cannot use the linear approximation owing to the non-diffentiability of the evolutionary dynamics around an equilibrium solution. For example, the Smith-dynamics version of the proposed evolutionary dynamics described by ˙ implying that the right-hand side of Eq. (16) is not always differentiable with respect to u i at u i = 0 .

Global behaviour of the evolutionary dynamics
In addition to investigating the stability of equilibrium solutions, we need to analyse the global behaviour of dynamics in the entire part of the feasible region, because the existence of a stable equilibrium solution does not directly imply that the state vector converges to it at all times. We may have a case in which the state vector varies within the feasible region forever, even if there exists a stable equilibrium solution. In such a case, we need to conclude that the dynamics may or may not be stable depending on the initial state. In a planar system (i.e. a system whose feasible region is two-dimensional), we may have three types of trajectories along which the state vector moves for t → ∞ , namely a closed curve that: 1. (Limit cycle) does not pass through any equilibrium solution. 2. (Homoclinic orbit) starts from and ends at a saddle equilibrium solution. 3. (Heteroclinic orbit) starts from a saddle equilibrium solution and ends at another saddle equilibrium solution.
In Section 5 , we will examine whether the evolutionary dynamics converge to a stable equilibrium point or one of these closed orbits.

Why classification?
In this section, we will introduce a few conditions on the utility function and classify all types of transport system that meet these conditions into nine cases. This classification has a central role in our analysis of the stability of equilibrium and global behaviour of the evolutionary dynamics.
The main reason why the classification is made is that, as mentioned in Section 3.1 , the proposed evolutionary dynamics is non-differentiable at its stationary point. As mentioned in Section 1 , there exists no valid standard methodology for analysing stability in such a case. The classification is an approach to deal with this issue. Because we are considering the planar dynamics, the number of cases should not be huge. Indeed, we will show that, under the conditions to be introduced, all possible situations can be classified into nine cases.
The classification is also beneficial for analysing the global behaviour of the transport system. We will use a graphical representation to clarify properties of the global behaviour. This is especially useful to see whether the state vector will converge to a stable equilibrium solution or be moving forever.

Conditions on the utility function for classification
We introduce the following two conditions [U3] and [U4] (in addition to [U1] and [U2], which have been already introduced in Section 2.2 ) on the utility function indicates the change of group j 's utility caused by the change of group i 's choices. We call it the internal interaction when i = j and the external interaction when i = j .
The following are implications of [U3] and [U4]: [U3] implies that the sign of ∂ i u j ( x ) is fixed to being either positive or negative in the entire feasible region. We refer to ∂ i u j ( x ) > 0 as a positive interaction and ∂ i u j ( x ) < 0 as a negative interaction .

has a sign that is different from those of all others, [U4] always
holds. We refer to this as the one-different-sign case.
cates how strong the internal interactions of both groups are and | J EXT ( x )| indicates how strong the external interactions between both groups are and hence: | implies that the external interaction is stronger than the internal interaction.

Shape of interior boundary and properties of equilibrium solutions
When [U3] and [U4] hold, the following theorem guarantees the monotonicity of the interior boundaries and characterises the steepness of their slopes, where x i B ( x A ) is a function describing the interior boundary for group i , defined by

Table 1
Classification of utility functions into four types and 16 cases.

Corollary 3.2. Any equilibrium solution is isolated, i.e. no other equilibrium solution exists in the neighbourhood of any equilibrium solution.
See appendix for their proofs.

Classification of nine cases
4.4.1. Enumerating 16 cases with respect to signs of ∂ i u j ( x ) [U3] and [U4] indicate that u ( x ) can be classified by paying attention to the signs of ∂ i u j ( x ). We employ a pictogram system shown in Fig. 2 to graphically depict the classification. These symbols are defined as follows: • The filled rectangle, diamond, and circle indicate that the internal interaction of each group is positive • The arrows indicate the external interactions. A solid (dotted) line indicates the positive (negative) external interaction from the tail to the head of the arrow (e.g. ∂ A u B ( x ) > 0 if the tail is group A and the head is group B).
• Except for the one-different-sign case (Cases 7 -9), the shape of the symbols indicates whether the internal interaction is stronger than the external interaction.
-If they are rectangles, the internal interaction is stronger . Circles are used in the one-different-sign cases.
The classification of the utility functions is shown in Table 1 generates two cases except for the one-different-sign cases.
• There exists eight one-different-sign cases (either one of has a different sign from the others, which is positive or negative).
which implies that we have 4 × 4 × 2 − 8 = 24 cases. Although there are only 16 cases in Table 1 , we can generate the other 8 cases by rotating the pictograms of the asymmetric mixed signs and one-different-sign cases by 180 degrees. Hence, the pictograms shown in Table 1 cover all possible types of the utility function. Table 2 Reducing 16 cases to 9 cases by inversing signs of external interactions.

Combining 16 cases into nine cases
The 16 cases indicated in Table 1 can be further combined into nine cases. For this purpose, we replace x A and x B with x A and x B , which are defined by: We also need to redefine u i ( x ) as follows so that [D3] is still valid after the replacement: Then, we have Eq. (22) shows that using x i instead of x i reverses the signs of the external interactions between the groups. Hence, we can further combine the 16 cases shown in Table 1 into 9 cases, as shown in Table 2 . These nine cases are all the cases we need to analyse in order to fully characterise properties of the proposed system.

Stability analysis of the proposed transport system
In this section we perform a stability analysis of the proposed transport system. An overview of the analysis is shown in Fig. 3 . We first introduce a graphical representation of the dynamics in Section 5.1 , then propose the methodologies to analyse stability of the interior equilibrium solution in Section 5.2 . The methodologies for the analysis of the peripheral equilibrium solutions are introduced in Section 5.3 . Finally, all the nine cases are examined using these methodologies in Section 5.4 .

Graphical representation of dynamics on phase diagram
We first introduce a graphical scheme that depicts the dynamics on the phase plane to analyse the evolutionary dynamics. This is useful for analysing both the global behaviour in the entire region of the feasible region and the stability analysis around an equilibrium solution.
To characterise the phase diagram by the direction of motion, sets S i + and S i − are defined as for i ∈ { A , B }. Owing to PC, hold. Hence, the interior boundary of group i can be used to divide the feasible feasible into sets S i + and S i − on a phase diagram like Fig. 4 . Recall that the interior boundaries are strictly monotone increasing or decreasing owing to Theorem 3 .  The sign of the internal interaction identifies the sign of the utility function on both sides of the interior boundary be- To graphically specify whether the internal interaction is positive on the phase diagram, we employ a sign consisting of two triangles indicating the direction of the dynamics along the axis of x i as shown in Fig. 4 . This sign directly determines the direction of motion along the axis because of PC. That is, ∂ i u i ( x ) > 0 causes the movement of the state vector x towards the end of the feasible region (i.e. x i = 0 or x i = 1 ), while ∂ i u i ( x ) < 0 causes the movement towards the interior boundary of group i .
Taking the union set of S i + and S i − we define sets S , S , S , and S , referred to as subregions as Each of them corresponds to the area surrounded by interior boundaries and the peripheral boundary. Owing to Eq. (23) , the direction of motion in a subregion is characterised by We use fan-shaped symbols to denote the direction of motion on the phase diagram as shown in Fig. 5 .
Describing the motion between regions is useful to analyse global behaviour of the dynamics. The direction of motion on the interior boundary is either horizontal (i.e. ˙ x B = 0 ) or vertical (i.e. ˙ x A = 0 ), implying that the motion can be depicted by a horizontal or vertical arrow on the phase diagram like Fig. 6 . Using such arrows is a popular technique to analyse nonlinear differential equation systems. The interior boundaries actually correspond to the null clines on the phase diagram ( Hirsch et al., 2004 , p.189). A set of points on an interior boundary whose directions of motion are the same is denoted by   We now obtain the phase diagram including all symbols like Fig. 7 , in which the phase diagram is accompanied with the corresponding pictograms of the utility function. We can intuitively read at least the following properties of the dynamics from the phase diagram: • Where the interior equilibrium solution is (i.e. the intersection of the interior boundaries). • Where each peripheral equilibrium solution is and whether it is stable or not. • Possible time-series of subregions in which a state vector stays.

Stability analysis of interior equilibrium solution and examination of limit cycle
In this section, we investigate stability of the interior equilibrium solution and examine the existence of a limit cycle. Analyses in this section largely depend on whether the trajectory of the state vector forms a loop around the interior equilibrium (called loopy trajectory ) or not. To examine the existence of a loopy trajectory, we first introduce the circular phase diagram and the state diagram in Section 5.2.1 . Then, we classify the nine cases into two categories, i.e. the cases in which a loopy trajectory may exist (the loopy case ) and that in which no loopy trajectory exists at all times (the non-loopy case ). The formal definition of these two cases will be given in Section 5.2.1 . Then, the non-loopy case is examined in Section 5.2.2 , followed by Section 5.2.3 , in which the loopy case is examined.
Note that the interior equilibrium solution may or may not exist (if it exists, it is unique according to Corollary 3.1 ). In the subsequent sections, we assume that it exists except in Section 5.2.2 , in which the case without an interior equilibrium solution is examined.

State diagram and circular phase diagram
The state diagram is proposed to describe the abstract global behaviour of the system. The state diagram consists of states corresponding to subregions and transitions describing possible motions from a state to another state. Fig. 8 depicts examples of them. Tracking the states on the state diagram along with transitions, we can enumerate any time-series of subregions that can be realised by the evolutionary dynamics.   To generate the state diagram, we now introduce a graphical presentation of the phase diagram that ignores the detailed shape of the interior boundary and the situation of the peripheral boundary (i.e. how interior boundaries intersect the peripheral boundary). The state diagram is made by considering the adjacency relationship of subregions and the direction of the slope of the interior boundary separating them on the phase diagram, while it is not necessary to know how the detailed shapes of the interior boundaries are and how they cross the peripheral boundary. The reason why it is possible is explained below.
Regarding the interior boundary, the direction of motion of the state vector on an interior boundary depends on the direction of motion in the adjoining subregions only. For example, when S and S adjoin through an interior boundary, the direction of motion on it is right horizontally on the phase diagram, regardless of the shape of the interior boundary.
To determine the direction of the transition between two states, we need to specify which state is on the upstream side of the direction on the phase diagram. In the above example, if S is on the left-hand side and S is on the right-hand side, the direction of the transition is from S to S . Owing to Theorem 3 , interior boundaries are always monotone increasing or decreasing and, if both interior boundaries are monotone decreasing (i.e. Cases 1-6), the slope of the interior boundary for group A is always steeper or shallower than that for group B. In Cases 7-9, one is monotone increasing and the other is monotone decreasing. Therefore, to determine position relationships between two subregions, we do not have to identify the detailed shape of the interior boundaries. Instead, we can use straight lines in the phase diagram for this purpose. The signs of the slopes of these straight lines are consistent with those of the original interior boundaries. In addition, the slope of the straight line for the group whose interior boundary has the steeper slope is steeper than that for the other group.
Regarding the peripheral boundary, as shown on the right-hand side of Fig. 9 , there are four types of crosses between an interior boundary and the peripheral boundary. To express all of these types at once, the rectangle of the peripheral boundary is replaced by a circle on the phase diagram, as shown on the left-hand side of Fig. 9 .
We call a phase diagram using the straight interior boundaries and the circular peripheral boundary a circular phase diagram . Fig. 10 depicts examples of the circular phase diagrams and state diagrams of the two cases shown in Fig. 7 .
The state diagram, consisting of four states, may have a trapping state , which does not have an arrow from itself to another state, and hence trap the state vector in it forever. The corresponding subregion is referred to as a trapping subregion .
More precisely, x ( t 0 ) ∈ S TP implies x ( t ) ∈ S TP ∀ t > t 0 , where S TP is a trapping subregion. For example, the right-hand-side case in Fig. 10 has two trapping states located on the top right and bottom left in the state diagram. We define that: • if the state diagram has a trapping subregion, it is called a non-loopy case , • otherwise, it is called a loopy case .
Obviously, in a non-loopy case, no loopy trajectory can exist owing to the existence of the trapping subregion. On the other hand, in a loopy case, a loopy trajectory may or may not exist.

Non-loopy case
The following theorem determines a main property of non-loopy cases: Theorem 4. (Non-loopy theorem) If the state diagram has a trapping state, the system always converges to a stable or saddle equilibrium solution.
Proof. See Appendix.
Note that converging to a saddle equilibrium solution implies the possibility of the existence of a homoclinic or heteroclinic orbit. In Section 5.4 , we will prove that no such orbit exists at all times.
Theorem 4 also derives the following corollary that determines the system behaviour when there is no interior equilibrium solution in the system (i.e. there is no crossing point of the two interior boundaries):

Corollary 4.1. (No-interior-point case)
If there is no interior equilibrium solution, the system always converges to a stable or saddle unstable peripheral equilibrium solution.
Proof. See Appendix. Now we are sure that no limit cycle exists, and the remaining task for the non-loopy case is to examine stability of equilibrium solutions. Note that stability of the peripheral equilibrium solutions will be examined in Section 5.3 . In this section, we will state the following theorem for the interior equilibrium solution: Theorem 5. (Instability of non-loopy case) The interior equilibrium solution x * is unstable (either diverging or saddle equilibrium solution) if there exists a trapping subregion S TP satisfying one of the following ∀ x ∈ S TP : Proof. See Appendix.

Loopy case
In the loopy case, the state vector departing from a subregion may stay in the same region forever or may come back to the same region after passing through other subregions. If the former case is applicable, the subregion where the state vector stays forever can be recognised as a trapping subregion, and hence Theorems 4 and 5 are applicable. If the latter case is applicable, the state vector always comes back to the same region after passing through other subregions. This case is investigated in the remaining part of this section.
We call a trajectory that comes back to the same region a spiral trajectory. More precisely, a spiral trajectory is a trajectory of the state vector x ( t ) for a finite time duration t ∈ [ t 1 , t 2 ], where 0 ≤ t 1 < t 2 , x ( t 1 ) ∈ I * and x ( t 2 ) ∈ I * , and I * is either one of I → , I ← , I ↓ , or I ↑ . Note that a spiral trajectory is defined for a finite duration of time, and hence is a subset of the trajectory of the state vector x ( t ) for t ∈ [0, ∞ ). Spiral trajectories can be classified into three types according to the relationship between x ( t 1 ) and x ( t 2 ), namely converging spiral trajectory , diverging spiral trajectory , and closed orbit . In the following, the word nearer implies that the distance from the interior equilibrium solution along an interior border is smaller (see Fig. 11 ). When x ( t 2 ) is nearer than x ( t 1 ), the spiral trajectory is referred to as a converging spiral trajectory. On the other hand, when x ( t 1 ) is nearer than x ( t 2 ), it is referred to as a diverging spiral trajectory. When x ( t 2 ) = x ( t 1 ) , the spiral trajectory forms a closed orbit.
The following theorem enumerates all the possibilities of the appearance of these types in spiral trajectories: Theorem 6. (Monotonicity of spiral trajectory) In a loopy case, depending on the initial state vector, we will find either one of three cases: 1. Converging spiral trajectories repeatedly appear; the state vector converges to the interior equilibrium solution or a limit cycle. 2. Diverging spiral trajectories repeatedly appear; the state vector converges to a limit cycle or a peripheral equilibrium solution. 3. Closed orbits repeatedly appear; the state vector is always on a limit cycle, which is identical to the closed orbit itself.

Proof. See Appendix
Theorem 6 is used to prove Theorem 7 . In addition, it also states that the series of spiral trajectories always converges to either an equilibrium solution or a limit cycle. This proposition is highly correlated with the Poincaré-Bendixson theorem, while we newly state this theorem because the dynamics considered in this study can have non-smooth points on interior boundaries.
Finally, the following two theorems are stated to judge the stability of the interior equilibrium solution surrounded by the spiral trajectories: Theorem 7. (Examination of direction of spiral trajectory) Consider a spiral trajectory. Consider any t 1 and t 2 where t 2 > t 1 and x ( t 1 ) and x ( t 2 ) are on the same interior boundary (i.e. x ( t 1 ) ∈ I * and x ( t 2 ) ∈ I * ). Consider a closed curve consisting of the spiral trajectory and an additional curve that connects x ( t 2 ) and x ( t 1 ) along the interior boundary. Let S C be a region that is surrounded by the closed curve but does not include any interior boundary. Let which is an integral of the divergence of the dynamics in region S C . Then, the following three propositions hold: • D > 0 ⇔ The spiral trajectory is a diverging spiral trajectory, which emanates from the interior equilibrium solution.
• D < 0 ⇔ The spiral trajectory is a converging spiral trajectory. which approaches the interior equilibrium solution.
• D = 0 ⇔ The spiral trajectory is a closed orbit.
Proof. See Appendix ( Theorem 6 is used here).

Proof. See Appendix
The combination of Theorems 7 and 8 invokes the following corollary, which will be used to judge the stability of the interior equilibrium solution when the trajectory around it is a spiral trajectory.

Corollary 8.1. (Checking direction of spiral trajectory) For any form of evolutionary dynamics satisfying conditions [D1]-[D4], if
x is sufficiently close to the interior equilibrium solution, 1. ∂ i u i (x ) < 0 ∀ i = { A, B } implies that the spiral trajectory is a converging spiral trajectory and hence the interior equilibrium solution is asymptotically stable.
2. ∂ i u i (x ) > 0 ∀ i = { A, B } implies that the spiral trajectory is a diverging spiral trajectory and hence the interior equilibrium solution is a diverging equilibrium solution.

Stability analysis of peripheral equilibrium solutions
Before examining stability, we first show that any peripheral equilibrium solution can exist either on an interior boundary or at a corner of the feasible region. Consider a peripheral equilibrium point x * on the left side of the peripheral boundary (i.e. x * A = 0 ). If x * is not on a corner, 0 < x * B < 1 , implying u B (x * ) = 0 owing to the definition of equilibrium by Eq. (5) .
Therefore, x * must be on the interior boundary of group B. A similar procedure can be used for all sides (i.e. upper, lower, and right side). We first examine stability of a peripheral equilibrium solution at a corner. Fig. 12 depicts all possible cases (Cases C-1 to C-4) of dynamics around a corner. Apparently, only Case C-4 has an equilibrium point, which is stable. Note that [U2] implies that any internal boundary does not intersect with a corner and hence we do not have to consider the case in which an internal boundary divides the neighbourhood of a corner into two different subregions.
We then examine stability of a peripheral equilibrium solution on an interior boundary. Note that [U2] implies that the interior boundaries do not cross on the peripheral boundary, and hence we need to consider the situation in which the peripheral equilibrium solution overlaps either one of them. Fig. 13 enumerates all possible situations of such an intersection point on the phase diagram. Note that it only depicts a horizontal side. We can turn the picture 90 degrees clockwise to obtain a phase diagram of a vertical side (see Fig. 14 ). There exist eight cases, which are identified by the type of the interior boundary intersecting with the side. Among them, only two cases (Case P-7 and P-8) have an equilibrium point. The equilibrium point of Case P-8 is stable because any trajectory of the state vector asymptotically converges to the equilibrium point. Case P-7 is a saddle equilibrium point, which is unstable.
We finally introduce a diagram, called the peripheral phase diagram, that depicts all possible patterns of peripheral equilibrium solutions on the left-half part of the phase diagram as follows: • For Cases 1 -6, because sgn ( ∂ A u i ) = sgn ( ∂ B u i ) ∀ i ∈ { A, B } and Eq. (17) of Theorem 3 , the slope of any interior boundary is negative. Hence, we can draw the phase diagrams of these cases like cases a to c in Fig. 15 . Because the slope of the interior boundary is negative, it can intersect with the peripheral boundary on the upper side or the left side only. We have three cases: (i) two interior boundaries intersect on the upper side, (ii) one intersects on the upper side and the other intersects on the left side, and (iii) two intersect on the left side.   • For Cases 7 -9, the slope of one interior boundary is positive and that of the other one is negative owing to the onedifferent-sign property. In these cases, cases d -g of Fig. 15  Specifying the motion of the direction in each subregion, we can identify the stable peripheral equilibrium point by applying the classifications shown in Fig. 15 . Note that, interior boundaries in these figures are drawn as a straight line, while it does not mean that this figure only considers a special case in which any interior boundary is straight. We only have to see the intersection points of them and the peripheral boundary and hence any monotone-decreasing (and monotone-increasing for cases d -g) curve is available for the interior boundaries. In addition, we only need to draw the left half of the phase diagram on the peripheral phase diagram because we can obtain its right half by turning 180 degrees.

Summary sheet and result chart
Using the methodologies proposed in Sections 5.2 and 5.3 , we finally show the results of the complete stability analysis for the nine cases in this section. The results are shown on the summary sheet , which summarises the circular phase diagram and state diagram proposed in Section 5.2 and the peripheral diagram proposed in Section 5.3 . The format of the summary sheet is shown in Fig. 16 . It also indicates relevant section numbers.
The result of the stability analysis is shown by the result chart on the right-hand side bottom of the summary sheet. It indicates where stable peripheral equilibrium solutions are (if they exist), whether there exits a stable interior equilibrium solution or not, and whether there may exist a limit cycle or not on (in) a box representing the feasible region. A possible position of a peripheral equilibrium solution is depicted by a mark like that illustrated in Fig. 17 . If a peripheral equilibrium  point certainly exists at a corner, a circle is marked on this corner. If it may or may not exist, a dotted circle is marked.
If it exists on a corner or a side (or sides), a box is marked in the corresponding place. If there is no mark, no peripheral equilibrium point exists. Marks shown in Fig. 18 specify the existence of the stable interior equilibrium solution and the limit cycle on the result chart. The result chart of the cases combined in Table 2 is also shown for convenience. They are shown in a parentheses. Because the replacement of variables introduced in Section 2.3 (i.e x A = x A and x B = 1 − x B ) corresponds to the inversion of the vertical axis (i.e. the axis for group B), the convergence chart of a combined case can be obtained simply by flipping the convergence chart of its original case.
The summary sheets for all cases are shown in Figs. 19 -24 . To complete the summary sheet, we need to draw the three diagrams first and then make the result chart. We can draw the diagrams by simply following the methodologies proposed in Sections 5.2 and 5.3 . In addition, the existence of a stable peripheral equilibrium solution is immediately clarified by the peripheral phase diagram. On the other hand, we need to examine the existence of the stable interior equilibrium solution and a limit cycle for each case. Details of this examination for each case will be explained in the next section.

Examining existence of the stable interior equilibrium solution and limit cycle
We examine the existence of the stable interior equilibrium solution and a limit cycle for each case as follows: • Cases 1, 3, 4, and 7: These are non-loopy cases and hence no limit cycle exists. In addition, according to the circular phase diagram, either (A) or (D) in Theorem 5 is applicable in cases 1, 3, and 7 and either (B) or (C) is applicable in case 4. Therefore, the interior equilibrium solution is unstable; it is either a diverging or saddle equilibrium point. • Case 2: This is a non-loopy case and hence no limit cycle exists. Theorem 5 is not applicable in this case. Indeed, we can prove that the interior equilibrium solution is stable (see Appendix for the proof). • Case 5: Dividing subregions is necessary to show that it is a non-loopy case. The subregions S and S is divided into two sub-subregions denoted by S * , S + and S * , S + , respectively as shown on the right-hand side of Fig. 25 . Using these Case 9  sub-subregions, the state diagram will be updated to that without a loop. Consequently, no limit cycle exists. In addition, according to the circular phase diagram, either (A) or (D) of Theorem 5 is applicable in this case. Therefore, there is no stable interior equilibrium solution in this case. • Case 6: Case 6 is a special case because the local stability of the interior equilibrium point may depend on the evolutionary dynamics and the utility function. It is a loopy case because the state diagram has no trapping state. Because is not applicable to Case 6. Indeed, D in Theorem 7 can be either positive or negative depending on , and x * , implying that the property of the spiral trajectory around the equilibrium point can change depending on their settings. • Case 8: This is a loopy case. A limit cycle may exist. If the state vector is trapped, the interior equilibrium point is not stable owing to Theorem 5 . Otherwise, the spiral trajectory is a diverging spiral trajectory if it is sufficiently close to the interior equilibrium point, owing to Corollary 8.1 . Consequently, the interior equilibrium solution is not stable. • Case 9: This is a loopy case. A limit cycle may exist. If the state vector is trapped, the interior equilibrium point is stable because the trapped state vector can only go towards the interior equilibrium point. Otherwise, the spiral trajectory is a converging spiral trajectory if it is sufficiently close to the interior equilibrium point owing to Corollary 8.1 . Consequently, the interior equilibrium point is stable.

Non-existence of homo-/heteroclinic orbit
The existence of a saddle equilibrium solution is a necessary condition for the existence of such an orbit. We first show the non-existence of these orbits associated with the peripheral saddle equilibrium solution, then the interior saddle equilibrium solution.
Regarding the peripheral equilibrium solution, if there exists an equilibrium point classified into Case P-7, there may be a homo-or heteroclinic orbit associated with it. We can actually state the following theorem: Theorem 9. No peripheral equilibrium point is associated with a homo-or heteroclinic orbit that survives forever (i.e. for t → ∞ ).

Proof. See Appendix
Regarding the interior equilibrium solution, Cases 1, 3, 4, 5, and 7 may have a saddle equilibrium solution. However, the state vector converging to the saddle interior equilibrium solution can go into a trapping subregion with a positive probability by small perturbations. Because the motion within a subregion cannot form a loop, no homo-/heteroclinic orbit can exist in these cases.

Assessing signs of
This section provides a method to quickly assess the signs of ∂ j u i ( x ) when a certain two-alternatives-two-groups transport system is given. We need to know the signs of ∂ j u i ( x ) to check whether [U3] and [U4] are satisfied or not and to classify the transport system into either one of the nine cases, whereas we basically formulate u k e. utility functions of each alternative for each group) in building a model of a transport system. For simplicity, we let ∂ k i = ∂ /∂ x k i and omit the arguments of function u i ( x ) and u k The proposed method in this section is applicable for any types of two-alternatives-two-groups transport system that satisfies either one of (1a) for all i , j ∈ { A , B } and either one of (2a) where −i is the group that is not i and h , k ∈ {0, 1}. Eq. (32) implies that the utility of an alternative is (1a) more or (1b) less influenced by the number of users of this alternative than by the number of users of the other alternative. Eq. (33) implies that the interaction within the same group is (2a) stronger or (2b) weaker than that between two groups.
We first examine the signs of ∂ j u i . Considering the relationship of ∂ x 1 i /∂ x i = −1 for i ∈ { A , B } (as has been explained in Section 2.1 ), we have and hence (1a) We then examine | J INT | and | J EXT |. We now have | ∂ i u i | > | ∂ −i u i | for case (2a) and | ∂ i u i | < | ∂ −i u i | for case (2b) owing to Eqs. (33) and (34) . Therefore, (2a) Note that both cases satisfy [U4].

Example: multi-vehicle type traffic assignment problem
As a typical example, we show that any generic form of a multi-vehicle type traffic assignment problem is categorised in Case 2 or 4. First, we assume that [U2] holds at any time. As explained in Section 2.1 , this condition is set to eliminate situations that are unlikely to hold in practice only. Then, let where: i , which is link traffic volume of group i , • T i l is a travel cost function of link l for group i , which is Lipschitz continuous, differentiable with respect to y i l , bounded (i.e. [U1] is satisfied), and strictly monotone increasing.
where h ∈ {0, 1} and j ∈ { A , B }. Because δ lki δ lk j = δ lki at all times but δ lki δ lh j = 0 for some l if k = h (note that there exists a link not shared by different routes), we obtain ∂ 0 We finally consider a single-vehicle type traffic assignment problem in which the two groups represent different OD pairs and T i where T l is a travel cost function commonly applicable to both groups. In this case, (2a) is applicable because there exists a link not shared by any combination of two routes of different OD pairs. Therefore, Case 2 is applicable. This is consistent with the well-known fact, i.e. there exists a unique stable equilibrium solution in a singlevehicle type network assignment problem.

Key examples of nine cases
It is worthwhile to propose key examples of the nine cases to understand how the proposed model can describe the abstract structure of evolutionary dynamics in a transport system consisting of both positive and negative interactions. Throughout all cases: • Social interactions cause positive interactions and congestion causes negative interactions.
• In cases 1, 3, and 6, different alternatives have no social or congestion interaction with each other, i.e. ∂ h i u k j = 0 unless h = k and consequently (1a) is applicable.
• In cases 5 and 7, alternative 1 is a public transport alternative whose travel cost is constant, implying that ∂ h (1a) is applicable.
We also explain the result of the stability analysis in each case. We use the following terms to explain the property of peripheral equilibrium solutions: • Concentrated solutions : users in both groups tend to select the same alternative. • Segregated solutions : users in the same group tend to select the same alternative, while those in different groups tend to select different alternatives.
If we place the word strongly before either one of these terms, it implies that all users in a group select the same alternative.

Case 1 ( ): Social interactions in a city
Each group represents a community in a city. More patronage in a mode increases its attractiveness owing to the social interactions (i.e ∂ k i u k j > 0 for all i , j , k implies ∂ i u j > 0 for all i , j ). Users in the same community have strong social interactions. Users in different user groups also have social interactions, which is weaker than that inside the same community (i.e. (2a) is applicable). At least we have strongly concentrated solutions. We may also have strongly segregated solutions in addition to the concentrated ones depending on the details of the transport system -in such a case we need to be aware that four stable equilibria actually exist.

Case 2 ( ): Single-vehicle type traffic assignment problem
See Section 6.2 for details.

Case 3 ( ): Social interactions with weaker congestions in different cities
There exist two cities and each group represents a community in a different city. Social interactions exist uniformly between travellers regardless of where they live and hence, like Case 1, ∂ i u j > 0 for all i , j . On the other hand, the congestion effect weakens the effect between travellers in the same city only (i.e. (2b) is applicable). In this case we always have two strongly concentrated solutions; there exists no chance to have a segregated solution, owing to the presumed existence of stronger positive external interactions.

Case 4 ( ): Multi-vehicle type traffic assignment problem in a parallel network
Suppose that drivers of cars and bicycles (represented by groups A and B) share the road. Driving on a road with other vehicles of the same type increases travel cost slightly. On the other hand, driving on a road with other vehicles of a different type increases travel cost significantly, so that (2b) holds. See Section 6.2 for details.

Case 5 ( , characterised by ): Contributors and free-riders
Consider a shared-ride transport system (represented by alternative 0) and a public transport (represented by alternative 1). Consider two types of travellers, i.e. contributors (group A) and free-riders (group B). While using the shared-ride, contributors tend to use their own vehicles that provide additional seats (i.e. ∂ A u i = ∂ A u 0 i > 0 ). On the other hand, free-riders tend not to provide seats but consume them (i.e. ∂ B u i = ∂ B u 0 i < 0 ). Shortage of seats causes a congestion effect, which is much stronger between free-riders (i.e. | ∂ B u B | is so big that (2a) is applicable). We observe concentrated solutions in this case. Owing to stronger congestion effect for free-riders, their choices can vary, whereas all contributors make an identical choice, either the shared-ride or public transport option.
6.3.6. Case 6 ( , characterised by ): Leader-follower social interactions There exist two groups of travellers, called leaders (group A) and followers (group B) in a city. Followers always observe leaders behaviour because they want to imitate it and leaders also have slight social interactions between them (i.e. ∂ 0 . Followers behaviour does not have such an influence; their patronage of a mode just makes congestion worse (i.e. ∂ 0 . Leaders' effects on followers are so strong that (2b) holds.
The ultimate system behaviour in this leader-follower example is highly unpredictable without an explicit modelling of the dynamical system, because the stability of the interior solution depends on the evolutionary dynamics.

Case 7 ( , characterised by ): New mode acceptable by users in one group only
This situation can occur when a new road transport mode (e.g. car-sharing system) is introduced in a road network, while some travellers do not accept it for some reasons (e.g. they cannot walk to a depot on foot). Travellers in group A travel by the new mode (alternative 0) or a public transport mode (alternative 1), while those in group B travel by private car (alternative 0) or public transport (alternative 1). Using any mode on the road increases the travel cost on the road due to congestion. Meanwhile, travellers using the new mode can enjoy a benefit from the positive externality because more patronage of the new mode makes it more convenient owing to economies of scale. Therefore, ∂ 0 We observe segregated solutions in this case, in which all users in group A make an identical choice, whereas choices in group B can vary owing to congestion.

Cases 8 ( , characterised by ) and 9 ( , characterised by ): Simplified morning commute problem
Commuters travel via a bottleneck on a road in the morning. There exist three time slots, denoted by τ = 1, 2, and 3.
Users in group A prefer time 2 (alternative 1) to time 1 (alternative 0) and never select time 3 because it is too late. Users in group B prefer time 3 (alternative 1) to time 2 (alternative 0) and never select time 1 because it is too early. Queueing congestion at the bottleneck has a strong temporal externality towards later times. Delay at the bottleneck for each time slot, denoted by w 1 , w 2 , and w 3 , (14) where γ > 1 is the externality ratio indicating how many times the external effect towards later times is stronger than that within the same time slot. Assuming that the schedule cost is constant for each time slot, we obtain If γ > 2, Case 8 is applicable. If γ < 2, Case 9 is applicable. Both cases can have a limit cycle, implying that this problem may not have a stable solution.

Numerical example
A numerical test of the simplified morning commute problem introduced in Section 6.3.8 is shown. The utility functions are defined by u 0 group A and B, respectively. Letting γ = 2 . 5 , c A = 1 . 75 , and c B = 2 we have which implies ∂ A u A = 0 . 5 , ∂ B u A = 1 , ∂ A u B = −1 . 5 , and ∂ B u B = 0 . 5 and this setting is classified into Case 8. The phase diagram is depicted in Fig. 26 . The interior equilibrium solution is (0.5,0.5) but it is unstable owing to the property of Case 8. We can also confirm from the phase diagram that no stable peripheral solution exists in the present example. Therefore, there exists no stable equilibrium solution. To observe this property, we numerically calculated Smith dynamics (shown in Eq. (16) ) from two initial points, x (0) = (0 . 4 9999 , 0 . 4 9999) and (0.9,0.9). The results are shown in Fig. 26 , in which we can confirm that the state vector converges to a limit cycle.

Discussions and conclusions
In this paper we have analysed a generic form of dynamical process for congestion games, inspired by theoretical and empirical transport literature that points to the existence of (what we refer to as) positive interactions in travel choice decision-making. While it has previously been shown that such forms of interaction may lead to multiple equilibria or unstable behaviour, our contribution has been to classify in detail the particular forms of dynamic system behaviour that may result, and the conditions under which they occur. We have deliberately developed our analysis to be applicable to a wide range of evolutionary behaviour, given the uncertainty over which is most appropriate to describe travel adaptation. This led to the particular challenge of dealing with non-smooth systems, precluding the use of standard first-order methods for dynamical system stability analysis.
In the paper we show how the many kinds of complex dynamics that arise from analysing such a problem with two alternatives and two user groups can in fact be reduced to essentially nine cases. Cases 1-4 (symmetric signs of interactions) always have a stable equilibrium solution. The existence of multiple equilibria is actually common. Only Case 2 can have the unique interior solution, implying that having a unique (balanced) solution is not such a common property as we may have thought. The asymmetric feature of the interactions may cause the non-existence of a stable equilibrium solution (i.e. in Cases 6 and 8), while most cases (all except case 9) can have multiple equilibria. There are two types of multiple equilibria: concentration (in which users tend to select the same alternative) and segregation (in which users in different groups tend to select different alternatives). The difference between them should be practically important because the former implies that either one mode is dominated or the latter implies that the users are segregated.
Our intention is that such a detailed classification may help users of dynamical models in one of several ways. Firstly, for cases in which the context is directly applicable (e.g. the two alternatives represent two modes) then the results can be used to infer cases (e.g. of multiple equilibria) in which measures or policies could be used to influence the transportation system towards a 'more desirable' outcome. Such policies may include not only 'hard' measures such as pricing or control, but also 'soft' measures to improve the flow of information within and between social groups. Secondly, for real-life examples that have a higher dimension (in terms of the number of alternatives or number of groups), our results serve as a qualitative guide to the model-user as to what kinds of system behaviour might be expected and which aspects of the model might influence the different kinds of system behaviour to prevail. Thirdly, we believe that such models may be exploited in the future for a more nuanced consideration of implementation pathways of policies, whereby alternative implementation plans may be tested in advance (e.g. smooth gradual improvement versus a kind of 'step-change').
Future research could investigate more complicated and realistic problem settings that can be represented by Cases 1 -9. It could clarify how the proposed methodology is helpful to have a perspective of an analysis for a particular transport problem before performing a detailed analysis. Considering evolutionary dynamics that is more empirically supported is also important. For empirical analyses, we would have to prepare a data set suitable for investigating dynamical changes of transport systems, which in the case of slowly adapting systems may involve analysing long-term transport big-data for decades or more.

Proof of Eq. (12)
Consider solving the evolutionary dynamics defined by Eq. (6) by following the trajectory in reverse to a past time.
Assume that there exists t 0 > 0 satisfying x ( t 0 ) = x * , where x * is the internal equilibrium solution. Then, the right-hand side of Eq. (6) is zero at t = t 0 owing to [D3]. Hence, x (t 0 ) = x * for any 0 ≤ t < t 0 . It contradicts [I], and hence x ( t ) cannot be x * at any time. Second, assume that there exists t 0 > 0 satisfying x i ( t 0 ) = 0 . Then, Eq. (6) becomes ˙ ) . Solving this ODE by following the trajectory to a past time, we will obtain x i (t) = 0 for any 0 ≤ t < t 0 because x i ( t ) cannot be negative.
It contradicts [I], and hence x i ( t ) > 0 at any time. Similarly, we can state that x i ( t ) < 1 as ˙

Proof of Theorem 3
Taking the derivative along the interior boundary which established Eq. (18) for Cases 1-6 because sgn (J INT (x )) = sgn (J EXT (x )) in these cases.

Proof of Corollary 3.1
Theorem 3 implies that the slope of the internal boundary for group A is always steeper (when | J INT ( x )| > | J EXT ( x )|) or shallower (when | J INT ( x )| < | J EXT ( x )|) than that for group B for cases 1 -6. For cases 7 -9, the signs of d x A B / d x A and d x B B / d x A are different at any time. Therefore, for any case, the number of intersection points of the internal boundaries is not greater than one.

Proof of Corollary 3.2
Corollary 3.1 implies any interior equilibrium solution is isolated. In addition, owing to the assumed strict monotonicity of u i ( x ), all peripheral equilibrium solutions are isolated.

Proof of Theorem 4
It is apparent that the state vector will stay in a trapping subregion after a sufficient time passes. It implies that the horizontal and vertical motion of the state vector is weakly monotone increasing or decreasing. Because the feasible region is finite, the state vector always converges to a stationary point. Therefore, it is an equilibrium point, which is either a stable equilibrium point or saddle equilibrium point.

Proof of Corollary 4.1
If there is no crossing point of two interior boundaries, only the following three cases are possible: • There is no interior boundary • One interior boundary divides the feasible region into two subregions. • Two interior boundaries divide the feasible region into three subregions.
The first case apparently has no loop because there is only one state. The second case also has no loop because the direction of the transition between two states is always unidirectional. In the third case, denote three subregions by S1, S2, S3, and let S2 be the subregion adjoining both S1 and S3. Then, both S1 and S3 do not adjoin each other. Therefore, there is no transition that generates a loop connecting the states corresponding to these three subregions.

Proof of Theorem 5
Only the proof for (A) is given here; all other cases can be proven by the same procedure. Consider a state vector (14) does not hold for this x * .

Proof of Theorem 6
The following lemma is first stated: Lemma. Consider any arbitrary trajectory x ( t ) for t ∈ [0, ∞ ). If there exists a pair of t 1 and t 2 that ensures trajectory x ( t ) for t ∈ [ t 1 , t 2 ] is a converging spiral trajectory, there exists no pair of t 2 and t 3 such that trajectory x ( t ) for t ∈ [ t 2 , t 3 ] is a diverging spiral trajectory or closed orbit. Similarly, if there exists a pair of t 1 and t 2 that ensures trajectory x ( t ) for t ∈ [ t 1 , t 2 ] is a diverging spiral trajectory, there exists no pair of t 2 and t 3 such that trajectory x ( t ) for t ∈ [ t 2 , t 3 ] is a converging spiral trajectory or closed orbit .
Proof. The former case is proven here; the latter case can be proven by an obvious adaptation of the same procedure. Fig. 27 depicts the converging spiral trajectory and an additional curve that connects two ends of the trajectory along with the interior boundary. Note that, on the additional curve, the direction of motion is always the same as that of the converging spiral trajectory. This trajectory and curve form a closed curve, called a Jordan curve, on the phase plane. As shown in Fig. 27 , the state vector at t 2 + 2 is in the interior region of the closed curve if 2 is positive and sufficiently small. If the spiral trajectory between t 2 and t 3 is a diverging trajectory, the state vector at t 3 − 3 is in the exterior region of the closed curve if 3 is positive and sufficiently small. Because the diverging spiral trajectory between t 2 and t 3 includes both t 2 + 2 and t 3 − 3 , this spiral trajectory must traverse both the interior and exterior regions of the closed curve. Owing to the Jordan curve theorem (see e.g. Spanier, 1966 ), this diverging spiral trajectory must intersect with the closed curve. Because any spiral trajectory can only intersect with the additional curve connecting the points of t 1 and t 2 from the exterior side to the interior side, there must exist a point at which the diverging spiral trajectory intersects with the converging spiral trajectory. Intersecting two trajectories is not possible at any time because the evolutionary dynamics does not allow any intersection point of the solution trajectory.
This lemma implies that we have two cases only for the spiral trajectory, i.e. converging trajectories are repeatedly generated forever or diverging trajectories are repeatedly generated forever. Let t i , where i is a positive integer, be a time satisfying x ( t i ) ∈ I * and t i < t i +1 . Because the feasible region is bounded, the distance (along the interior boundary) between x ( t i ) and x ( t i +1 ) converges to zero. This implies that the spiral trajectories converge to a closed orbit or to the interior equilibrium point.

Proof of Theorem 7
This theorem can be proven by reference to Green's Theorem, while it needs to be noted that Green's Theorem can be used only when the function in the integral is differentiable. To avoid this problem, divide S C into four regions by the interior boundaries (See Fig. 28 for an example). Denote these divided regions by S 1 C , S 2 C , S 3 C , and S 4 C . Green's Theorem implies that:  where ∂S i C is the anticlockwise path along the border of S i C . Because ˙ x A and ˙ x B are continuous with respect to x , the path integrals along the interior boundaries cancel out each other in Eq. (44) . Therefore, where P 1 is the curve of the spiral trajectory and P 2 is the additional curve connecting both ends of P 1 . Note that P 1 is always anticlockwise regardless of the actual direction of motion of the spiral trajectory. Because P 1 is the curve along the evolutionary dynamics, implying that As mentioned in the proof of Theorem 6 , the direction of motion of the state vector on P 2 is always the same at any point on it. In addition, P 2 is not along the solution trajectory. Therefore, Eq. (47) implies the following relationships (where the term of right and left is defined by setting path P 2 upright and going from the bottom to the top): • D > 0 ⇔ The direction of motion is from the left to the right relative to path P 2 .
• D < 0 ⇔ The direction of motion is from the right to the left relative to path P 2 .
• D = 0 ⇔ The length of P 2 is zero.
In the following, without loss of generality, we only consider the anticlockwise spiral trajectories. For a converging spiral trajectory, we can see that the direction of motion is from the right to the left relevant to path P 2 (see Fig. 29 left). Therefore, the converging spiral trajectory implies D < 0. For the diverging spiral trajectory, we can see that the direction of motion is from the left to the right relevant to path P 2 (see Fig. 29 right). Therefore, the converging spiral trajectory implies D > 0. For the closed orbit, D = 0 is apparently established. Combining these three propositions, we can also derive their converses because they cover all cases of D (i.e. negative, positive, or zero).

Proof of Theorem 8
This theorem is proven by linear approximation of ρ + i and ρ − i , which can be calculated by Note that the partial derivative of x i is applied to the first argument of ρ + i and ρ − i only; x in u i ( x ) is considered constant when these derivatives are taken.
The first terms of the right-hand side of Eq. (48) are examined first. They can be rewritten as ∂ i ρ + i ( x * , 0) and ∂ i ρ − i ( x * , 0) by excluding the term (x i − x * i ) . Because ρ + i (x , 0) = ρ − i (x , 0) = 0 and they are Lipschitz continuous and bounded, Proof of Corollary 8.1 From Eqs. (30) and (31) , we now have for x ∈ S − ∂S − I, where I is the set of all interior boundaries. Letting x be sufficiently close to the interior equilibrium point, we can let u i ( x ) be so small that signs of the right-hand sides of Eq. (51) are determined by their second terms, i.e. sgn ∂ ˙ Note that [U2] indicates ∂ i u i ( x ) = 0. Eq. (52) can be further converted to sgn ∂ ˙

Proof of stability of Case 2 in Section 5.4.2
Assume that the state vector is sufficiently close to the interior equilibrium point so that the interior borders can be approximated by straight lines. Consider a 'box', denoted by S BOX ( t ), around the interior equilibrium point whose upper left and lower right corners overlap with the interior border of Group B (depicted by a blue dashed line). The boundary of this box always overlaps with the state vector x ( t ). See Fig. 30 for the settings described above. It implies that: • the upper side of the boundary is in the subregions S and S , • the lower side of the boundary is in the subregions S and S , • the left side of the boundary is in the subregion S , and • the right side of the boundary is in the subregion S . Therefore, all the sides of the boundary of S BOX ( t ) always converge to the interior equilibrium point (i.e. x LU ) and hence the state vector on it monotonically converges to the interior equilibrium point, implying that it is asymptotically stable.

Proof of Theorem 9
The equilibrium point of Case P-7 is not stable because the trajectories go in the opposite direction to the equilibrium point in the horizontal direction. Such trajectories always converge to a stable equilibrium point that is on a corner unless there exists another intersection point between the equilibrium point on the side and a corner, as in Fig. 31 . Considering that three or more intersection points cannot exist on the same side (otherwise, at least one interior boundary would have to be bent like in Fig. 31 , but this is not possible owing to Theorem 3 ), at least trajectories going towards one corner will converge to it. Therefore, even if there exists a homo-or heteroclinic orbit associated with the point of Case P-7, this trajectory cannot survive for t → ∞ because, with a positive probability, a perturbation to the state vector around this equilibrium point pushes it towards the stable peripheral equilibrium point at a corner.