Two-Agent Games on Graphs

TL;DR

This paper investigates the evolution of agent points on graphs through two-agent competitions, analyzing dynamics on general graphs and specifically on ring and torus structures, with connections to statistical physics models.

Contribution

It introduces a model of competitive dynamics on graphs and explores its properties, especially on regular lattice structures like rings and tori.

Findings

The model relates to vertex models and Potts models at zero temperature.

Dynamics exhibit specific behaviors on ring and torus graphs.

Theoretical analysis of evolution patterns on different graph types.

Abstract

We study the dynamics of evolution of points of agents placed in the vertices of a graph within the rules of two-agent units of competition where an edge is randomly chosen and the agent with higher points gets a new point with a probability $p$. The model is closely connected to generalized vertex models and anti-ferromagnetic Potts models at zero temperature. After studying the most general properties for generic graphs we confine the study to discrete d-dimensional tori. We mainly focus on the ring and torus graphs.

Equations

$$(L,S)$$

$$\sum _{i=1}^N{X}^i(k)=X_o+k,$$

$$\mathcal{P}(\vec{x},k+1)=\sum _{i=1}^N\mathcal{P}(\vec{x}-{\hat{e}}_i,k)v_i(\vec{x}-{\hat{e}}_i).$$

$$v_i(\vec{x})=\frac{1}{E}\sum _{j\in {\mathcal{N}}_i}\left[p\theta \left(\vec{x}\cdot ({\hat{e}}_i-{\hat{e}}_j)\right)+q\theta \left(\vec{x}\cdot ({\hat{e}}_j-{\hat{e}}_i)\right)\right]$$

$$\sum _{i=1}^N{v}_i(\vec{x})=1.$$

$$v_i(\vec{x})=q\frac{d_i}{E}+\frac{p-q}{E}\sum _{j\in {\mathcal{N}}_i}\theta (x_i-x_j).$$

$$\vec{v}(\alpha \vec{x})$$

$$\vec{v}(\vec{x}+\alpha \hat{\mathcal{E}})$$

$$\hat{\mathcal{E}}=\frac{1}{\sqrt{N}}\sum _{i=1}^N{\hat{e}}_i$$

$$\vec{v}(\vec{x}|x_i=x_j)=\frac{1}{2}\left[\vec{v}(\vec{x}|x_i>x_j)+\vec{v}(\vec{x}|x_i<x_j)\right].$$

$${v_{\parallel }}^{+}-{v_{\parallel }}^{-}$$

$$v_{\perp }^{+}-v_{\perp }^{-}$$

$$q(d_i-d_j)+(p-q)(l_i-l_j)+(p-q)({ϵ}_i-{ϵ}_j)/2=0,$$

$$\frac{\partial \mathcal{P}}{\partial t}=-\vec{\nabla }\cdot \left(\vec{v}\mathcal{P}\right)+\frac{1}{2}\sum _{i=1}^N{{\nabla }_i}^2\left(v_i\mathcal{P}\right)+\cdots .$$

$$\frac{\partial \mathcal{P}}{\partial t}=-\vec{v}\cdot \vec{\nabla }\mathcal{P}-\left(\vec{\nabla }\cdot \vec{v}\right)\mathcal{P}.$$

$$\frac{d\vec{x}(s)}{ds}$$

$$\frac{dt}{ds}$$

$$\frac{d\mathcal{P}(\vec{x}(s),t(s))}{ds}$$

$$\vec{\nabla }\cdot \vec{v}=\frac{p-q}{E}\sum _i\sum _{j\in {\mathcal{N}}_i}\delta (X_i-X_j),$$

$$⟨F(\vec{x})⟩(k)\equiv \sum _{\vec{x}}\mathcal{P}(\vec{x},k)F(\vec{x}).$$

$$⟨F(\vec{x})⟩(k+1)=\sum _i⟨v_i(\vec{x})F(\vec{x}+{\hat{e}}_i)⟩(k).$$

$$⟨\vec{x}⟩(k+1)-⟨\vec{x}⟩(k)$$

$$⟨\vec{v}(\vec{x})⟩(k+1)$$

$$\vec{X}(k+1)=\vec{X}(k)+\vec{v}(\vec{X}(k)).$$

$$\frac{d⟨F⟩}{dt}=⟨\vec{v}\cdot \vec{\nabla }F⟩.$$

$$\frac{d⟨\vec{x}⟩}{dt}$$

$$\frac{d⟨\vec{v}⟩}{dt}$$

$$a_i=\frac{p-q}{E}\sum _{j\in {\mathcal{N}}_i}(v_i-v_j)\delta (X_i-X_j).$$

$$\underset{t\to \infty }{\lim }\frac{d⟨\vec{v}⟩}{dt}\to 0,$$

$$\underset{t\to \infty }{\lim }⟨\vec{X}⟩\to ⟨\vec{v}⟩t+⟨\vec{X}⟩(t=0).$$

$$\underset{t\to \infty }{\lim }⟨\vec{X}⟩\to ⟨\vec{v}⟩t.$$

Full text

I Introduction

Mathematical modelling is becoming more and more important for problems not only in natural but for economical and social sciences in general rev1 ; rev2 . In this work we study the dynamics of a model which rests on two-agent competitions first introduced in red1 ; the model has attracted attention and has found various implementations rede1 -rede10 . The model was later generalized to 3-agent games 3pap ; rad1 and it also created a good laboratory for studying effects of mergers rad2 . The 4-agent version of the model was extensively studied in der1 .

The works on the mentioned model focused exclusively on the implementation of the game on the complete graph which would constitute a mean field approximation for results on generic graphs. In this work we apply the model to general graphs and study various of its properties. What one can say with generality is that the model is a generic vertex model on graphs, and hence is also related to a spin 1/2 system where spins are living on the edges of the graph. One can also say with generality that in the competitive limit the model is mapped to a constrained graph colouring problem. In that sense it is also kin to anti-ferromagnetic Potts models at zero temperature with constraints on the configurations. Furthermore because the system is a generic vertex model its dynamics is also described within the mathematical field of hyperplane arrangements; the hyperplanes being defined with neighbouring vertices.

In II we introduce the model and study its general properties. In III we study the hydrodynamical approximation and discuss its properties. In IV apply our findings to some finite graphs, including the complete graph. Section V is reserved for an in depth study of the model on the ring graph, where we show that the configurations are 3-colourings of the ring such that there are checkerboard islands of two colours separated by one colour. The size distribution of these islands and the correlation function of the colourings are also studies. We have also implemented a chemical potential approach to study the statistics of the model and show that it is quite accurate even for the details of the model. In section VI we study the model on the two dimensional discrete torus and show that the sizes and perimeters of the mentioned checkerboard islands follow a scaling law. In section VII we digress on the behaviour of the model for d-dimensional tori.

II The Model and its General Properties

Let us assume we have $N$ agents occupying vertices on a connected graph with $E$ edges and with no bubbles so that an agent is not connected to itself by an edge. We label their position with an integer $i=1,2,\cdots,N$ and endow them with integer points $X^{i}$. Let us now randomly pick an edge on an equal likelihood basis. This unambiguously defines two agents. Say, their points are $L$ and $S$ respectively and obey the ordering $L\geq S$. We give a single point to one and only one of them according to the following rules

[TABLE]

whenever $L>S$. If, on the other hand $L=S$, we shall still grant the point, but now on an equal likelihood basis;

[TABLE]

Since we have $N$ agents and they have integer points, the state of the system can be described with a vector $\vec{X}=X^{i}\hat{e}_{i}$ in an $N$ dimensional cubic lattice vector space with an orthonormal basis $\hat{e}_{i}$. Note that since we always give a single point after a single unit of competition, after $k$ such steps the state of the system obeys,

[TABLE]

where $X_{o}=\sum_{i}{X_{o}}^{i}$ is the sum of the points for the initial state.

The problem at hand is to find $P(\vec{x},k)$ which is the probability at step $k$ that the random variable $\vec{X}$ has value $\vec{x}$ given that at $k=0$ one has the distribution $P(\vec{x},0)$. If the state vector at step $k+1$ is $\vec{X}$, it must have been at a neighbouring vector at step $k$ according to the rules above. We can therefore write down the following master equation for the probability distribution

[TABLE]

Where

[TABLE]

is, determined by the microscopic rules, the probability that agent $i$ wins under the condition that one has $\vec{X}=\vec{x}$. In the above $\mathcal{N}_{i}$ is the set of agents neighbouring agent $i$ excluding agent $i$ itself. Also we define ${\rm card}(\mathcal{N}_{i})\equiv d_{i}$, the degree of the vertex associated with $i$’th agent. We use $\theta(0)=1/2$ to implement the tie breaker and also one has $\theta(u)+\theta(-u)=1$. Since $p+q=1$, using the handshaking theorem which states that $\sum_{i}d_{i}=2E$, one gets

[TABLE]

One can thus recast $v_{i}$ as follows,

[TABLE]

We note that for $p=q=1/2$ one has $v_{i}=qd_{i}/E$ and the impact of step functions which depend on the ordering of points is absent and $\vec{v}=\sum_{i}\hat{e}_{i}v_{i}(\vec{x})$ is the same constant vector everywhere in the configuration space.

There are two essential symmetries of the object $\vec{v}$ which we would like to call the velocity field. These are expressed as

[TABLE]

Here

[TABLE]

is a unit vector along the body diagonal in the $N$ dimensional configuration space.

The velocity field has other properties related to the fact that its dependence on $\vec{x}$ is only via step functions. Let us consider the hyperplane $x_{i}=x_{j}$. The microscopic rules can be used to show that

[TABLE]

One can show that this averaging property holds also for the intersections of the hyperplanes and the intersections of the intersections and so on. Of course for a given intersection one will have to average over all hyperplanes, intersections of the hyperplanes and so on that are connected to the intersection.

Let us again consider two agents $i$ and $j$ that are linked. Since a hyperplane will be of dimension $N-1$, the unique normal to it is along the vector $\hat{e}_{i}-\hat{e}_{j}$; that is the perpendicular component is $v_{i}-v_{j}$. One can thus show that the velocity field obeys the following boundary conditions

[TABLE]

where the superscripts $\pm$ refer to $x_{i}>x_{j}$ and $x_{i}<x_{j}$ respectively.

Various results are implied from the discontinuity of $v_{\perp}$. For $p>q$ we realize two possibilities: Either in both sides the velocity field is pointing away from the hyperplane or in one side it is pointing towards the hyperplane and on the other pointing away such that the one that points away has larger magnitude, which means that for both cases the hyperplane can only be a transient location. For $p<q$ things are a bit more involved. We still have the case where in one part velocity field is pointing away and the other pointing towards such that the one pointing towards has larger magnitude; the opposite situation on the first case for $p>q$. But we also have the possibility that in both regions the velocities are towards the hyperplane. So for $p<q$ we have the possibility that the hyperplane may be recurrent. However, we have to remind the reader that a hyperplane is never absolutely stable in view of the fact that we have a tie breaker rule in the microscopics of the model which ensures that one can not indefinitely stay on the hyperplane. But the velocity field on the regions between the hyperplanes might enforce an eventual return to it. So for $q>p$ we may have agents condensed around the hyperplanes, coming in and out from both regions. Pertinent to this, is the case where $v_{\perp}^{-}=-v_{\perp}^{+}$ which would mean that on the hyperplane $x_{i}=x_{j}$ the perpendicular component of the velocity field vanishes.

Another possibility is to have $v_{\perp}$ vanishing. It is clear, due to the discontinuity in $v_{\perp}$, that $v_{i}=v_{j}$ can occur only in one region separated by the $x_{i}=x_{j}$ hyperplane. The condition for $v_{i}=v_{j}$ is the requirement that

[TABLE]

where $l_{i}$ is the number of neighbours of $i$’th agent that has point strictly less than $x_{i}$. Also $\epsilon_{i}$ is the number of neighbours of agent $i$ with which agent $i$ has the same point. Surely in a region where there is an unambiguous ordering of points all $\epsilon_{i}$ are equal to zero. This condition generally depends on $q$ if the graph is rather generic. But if we have a regular graph where $d_{i}=d_{j}$ the condition does not depend on $q$ anymore  111Except for the case $q=1/2$ when it is identically satisfied.. Now if $v_{i}=v_{j}$ the velocity field in this region will be entirely parallel to the hyperplane defined by $x_{i}=x_{j}$. However due to the discontinuity in $v_{\perp}$, on the other side of the hyperplane the velocity field will definitely point away from the hyperplane for $p>q$ and toward the hyperplane for $p<q$.

One last bit of information about the velocity field is that since by construction all $v_{i}\geq 0$, nowhere the velocity field points towards any $x_{i}=0$ hyperplane. This is also evident from the fact that within the microscopic rules there is no lowering of points; at any step every agent either stays at the same point or increases its point.

II.1 A buckets-pipes analogy for the model

We can form an analogy of the model for $p=0$ or $p=1$ using a network of buckets of water connected with pipes presented in Fig.1. The buckets represent the agents, the pipes represent the links and the water level at each pipe represent the points collected. When we pick a link we shall drop a defined unit of water on the pipe and let it flow to one of the buckets. If the water level is the same for two neighbouring buckets we force the water flow randomly to only one to comply with our rules 222One could, in this case, split the water evenly to both buckets. However giving points to ties do not change the general structure of the model..

For $p=1$ when we pick a link the agent that has more water in the bucket is to surely receive the unit volume of water. So we can visualize the system as follows. Assume all buckets are connected to the ceiling with linear elastic bands and the buckets are connected to each other by pipes in a way described by the graph. Let us also assume that the pipes are connected to the top of the buckets. Now let us start with a given water level distribution for the buckets and drip water on one of the pipes chosen at random and assume that the water will flow instantaneously or that we wait until the system balances via the presence of the elastic bands until we pick another pipe. Under the action of gravity the water will flow to the bucket that has more water and this establishes the analogy. The case $p=1$ however seems to have a more natural connection to competitive systems of interacting economical entities. The dropping of water will simply stand for taking a share from the overall production.

For $p=0$ having a larger point then one’s neighbours must be punished. The way to implement this behaviour is to put the buckets on the floor and connect the neighbouring ones again with pipes. However this time do not let the pipes to have a fixed position with respect to the buckets but rest them on the water surfaces of the buckets with buoys which enforces the water flow under the action of gravity to the bucket with less water establishing the analogy.

This analogous model is very useful in assessing the qualitative aspects of our system. For instance, if there are cycles in the graph it gives an understanding that some agents might still have different water accumulation speeds even for $p=0$; since there can be a water flow trapped in a portion of the graph. If on the other hand we are talking of a tree graph the flow seems to always disperse and at $p=0$ all agents travel with the same rate of water level increase. These statements will be justified later in the manuscript. Furthermore the case for $p=0$ may find applications for water systems that are connected with each other, such as underground caves; the dropping of water in this case simply being representing rainfall. If on the other hand we still want to make an analogy to interacting entities in an economical system, the case $p=0$ resembles the case for a ”rich don’t get richer” model. Yet as we have mentioned the number of connections an agent has is also relevant, particularly in an economical system. An inequality in the number of neighbours generally tends to upset the egalitarian aspect for $p=0$ if there are cycles in the graph.

III The hydrodynamical limit

In the large time and hence large points limit one might contemplate a continuum approach akin to representing the model with a hydrodynamical fluid. In this approximation the discrete equation for $\mathcal{P}$ becomes

[TABLE]

The first term in the equation above describes a convective flow and the second is responsible for diffusion. Whenever justifiable one might ignore diffusion next to convection at very late times. Later in the manuscript, we shall elaborate on what we mean by justifiable when we discuss how the diffusion terms, and even sometimes all the terms in the expansion, may become important. For now we simply state that whenever the velocity field in a region points away from all the relevant equal points hyperplanes that bounds it the dominant term is the convective term. Focusing thus on the first term we end up with

[TABLE]

This is a linear hyperbolic equation and the general method of approach to solve it is that of characteristics. To that end one needs to solve paths in the configuration space such that along them the equation is satisfied. If those paths are described by a parameter $s$ we will have to solve the following set of equations

[TABLE]

The last equation in the list above is the material derivative of the density function in the fluid dynamics context and will be responsible in assessing the type of the flow. But since in the regions separated by relevant hyperplanes $\vec{v}$ is a constant vector field, its divergence is non-zero only along hyperplanes where neighbouring agents have equal points.

So having focused on the first non-trivial term in our master equation, what one can say is that, in the bulk regions defined between these hyperplanes the dynamics of our model is essentially that of an incompressible flow.

On the other hand, the divergence of $\vec{v}$ is given as

[TABLE]

and thus its signature is solely determined by the signature of $p-q=2p-1$. For $p>1/2$ this coefficient is positive and the density on the hyperplanes is diluted; the flow generally points away from them. If on the other hand $p<1/2$ the density on the hyperplanes is condensed; the flow generally tends toward them where the essentials depend heavily on the particulars of the graph.

So on quite general grounds one can say that for $p>1/2$ nearest neighbours have points that diverge from each other in the infinite time limit and thus the model is related to graph colouring where the colors are the asymptotic speeds. However what one can say with confidence for $p<1/2$ case is only that the nearest neighbours might have converging points and hence same speeds since this might be locally frustrated by the distribution of the number of neighbours $d_{i}$ on the graph. To make this point clear and elaborate on a point mentioned before, let us ask some questions: For $p=0$ does every agent for a given graph travel at the same asymptotic speed? If yes what was the critical value of $p$ when this is reached? If not, at $p=0$, what is the number of possible different speeds for a given graph? We shall have partial answers to these later in the manuscript.

III.1 Expectation values

Let us define the following

[TABLE]

Using the master equation it is easy to show, after a shift in the summation variables, that

[TABLE]

The expectation values of $\vec{x}$ and of $\vec{v}(\vec{x})$ are of central importance. We have

[TABLE]

Note that the equation for the expectation value of the state vector directly follows from the equation obeyed by the random variable $\vec{X}$ which is

[TABLE]

On the other hand keeping only the convective flow term in a continuum approach, and using $t$ instead of $k$ we shall get

[TABLE]

The expectation values of $\vec{x}$ and $\vec{v}$ in this approximation are,

[TABLE]

Note that the expectation value of $\vec{x}$ is already saturated by the convective approximation. Also, as in fluid dynamics, the quantity $\vec{v}\cdot\vec{\nabla}\vec{v}\equiv\vec{a}$ is the acceleration field. Its components are given as

[TABLE]

Evidently the acceleration is only non-zero when some neighbouring agents have the same points and having different speeds. Thus in the bulk regions, where no neighbouring agents have the same points, the acceleration vanishes and the expectation value of the velocity field converges to a constant value. Furthermore one can still have a vanishing acceleration on a neighbourhood hyperplane $x_{i}=x_{j}$ if $v_{i}=v_{j}$, as mentioned before.

Therefore

[TABLE]

generally in the bulk regions and also on an hyperplane if the perpendicular component vanishes.

This in turn means that

[TABLE]

As an example if one picks $\mathcal{P}(\vec{X},0)=\delta^{N}(\vec{X})$ one shall get $\langle\vec{X}\rangle(t=0)=0$. In fact if we assume that the initial distribution has finite moments we can generally say that

[TABLE]

To reiterate let us know focus on $\langle X_{a}X_{b}\rangle$. It is easy to show that

[TABLE]

Taking a further derivative and focusing on the bulk regions we can arrive at

[TABLE]

However as we have shown that the acceleration field vanishes in the bulk regions the quantity $\langle v_{a}v_{b}\rangle$ approaches a constant in the infinite time limit, and we thus expect the following if the initial distribution of points has finite moments

[TABLE]

In fact by assuming finite moments for the initial distribution the following general relation holds in the bulk regions of the configuration space for the asymptotic dynamics

[TABLE]

where $S=\sum_{k}p_{k}$.

III.2 The scaling symmetry of $\vec{v}(\vec{x})$ and the $\mathcal{P}(\vec{x}/t)$ ansatz

These facts above are closely related to the scaling symmetry of the velocity field. If we remember that $\vec{v}(\vec{x}/t)=\vec{v}(\vec{x})$ we see that the equation governing $\mathcal{P}$ will allow an asymptotic solution of the form $\mathcal{P}(\vec{x}/{t})$ provided the initial condition does not spoil it. Such types of initial conditions are for instance starting all agents from the same point, or say sampling initial conditions from inside an $N$-dimensional cube or sphere in the configuration space on an equal likelihood basis.

Defining $\vec{z}\equiv\vec{x}/t$, the equation in the bulk regions will have the form

[TABLE]

This can also be approached via the method of characteristics. The system one will have to solve is the following

[TABLE]

Since in the asymptotic regime we have found that the random variable $\vec{X}/t$ approaches a constant, which is $\vec{v}$, we arrive at the defining relation for an asymptotic solution

[TABLE]

along with

[TABLE]

and the probability that a configuration satisfying the condition above is a constant333A humourous person would write Eq.31 as $\vec{v}=\vec{v}(\vec{v})$.. Also we note that since the velocity field depends only on step functions which take only two values, the configurations defined by Eq.31 are those of a generalized vertex model on the graph.

So what we have to do has two aspects: One has to assess the list of possible configurations that satisfy Eq.31 and compute probabilities of occurrence for these configurations. The list of configurations depends only on the graph properties via the velocity field and generally on $q$ but the asymptotic occurrence probabilities of these will in general depend on the initial conditions. Thus we expect broken ergodicity on general grounds.

Let us now remember that the exact behaviour of the expectation value of $\vec{x}$ is the same as in the lowest order approximation in the hydrodynamical approximation. Therefore Eq.31 will also be valid for $p<q$ and this means we can find all the possible speeds and hence possible configurations of the system for all $p$. The generic recipie is that we look for all acyclic orientations of the graph that will yield a colouring for $p>q$; some will provide solutions and some will not. For $p<q$ however we are allowed to let neighbouring agents have the same speeds and one has to also look for those possibilities; that is say partial orientations. Again, some of these may provide solutions and some may not. We shall elaborate more on the generalalities of Eq.31 later in the manuscript.

III.3 Universality classes $p\leq 1/2$ and $p>1/2$ for regular graphs

For regular graphs there is enormous simplification of Eq.31. Since on a regular graph the number of neighbours are the same, say $d$, we have

[TABLE]

Due to the nature of the step function the following transformation

[TABLE]

will bring the equation to the form

[TABLE]

where $\epsilon=\rm{Sign}(p-q)$. This hints at two universality classes of different values of $p$ separated only by the signature of $p-q$.

To reiterate on the property just mentioned, let us go back to the original master equation. In the continuum limit keeping only the convective flow terms one can show that it can be written as

[TABLE]

where $\vec{D}=\sum\hat{e}_{i}d_{i}/E$ is the vertex degree vector and $\omega_{i}=\sum_{j\in\mathcal{N}_{i}}\theta\left(\vec{x}\cdot(\hat{e}_{i}-\hat{e}_{j})\right)/E$. Surely $\vec{\omega}$ is the same vector as $\vec{v}$ for $p=1$.

For a regular graph one has $D_{i}=d=2E/N$. Let us also remember that the velocity field enjoys a scaling symmetry and a shift symmetry if the shift is along $\hat{\mathcal{E}}$. Thus under the following transformation

[TABLE]

the velocity field is unchanged. Consequently the equation can be transformed to

[TABLE]

Which is the same equation as our original starting point had we picked $p=1$. Thus the argument for the mentioned universality classes does not depend on the $\mathcal{P}(\vec{x}/t)$ ansatz.

To summarize, for regular graphs in the $p\leq 1/2$ regime all agents have the same asymptotic speed of $1/N$ and for $p>1/2$ the problem is statistically equivalent to the case $p=1$.

Of course in an actual simulation of the system if $p\neq 1$ but still $p>1/2$ the convergence of the system to the statistics described by $p=1$ should take longer, a fact we have substantiated in simulations. The reason is that the characteristics diverge from each other in a slower way if $p\neq 1$. Thus this universality is welcome and to be exploited to make simulations run faster since for $p=1$ one does not have to pick a random number to assess which of the two neighbours that were picked will win the game.

Another property of regular graphs is that for all of them one can consider a thermodynamic limit where the number of agents is taken to infinity while keeping the properties of the graph intact.

III.4 Connection to hyperplane arrangements and graph colourings

The step functions in the velocity field divide the configuration space into regions where the velocity field is constant but have strengths and directions that depend on the region. This makes an immediate connection to a branch in combinatorial geometry called hyperplane arrangements. In particular since the arrangements of these hyperplanes depend on a given graph they fall into the category of graph arrangements. For definitions and theorems on this particular topic which we shall cite below we refer the reader to orlik and stanley .

The arrangements in our model are central ones, meaning that all hyperplanes include the origin. Because of this it is easy to show that they are all unbounded using Zaslavsky’s theorem: They bound each other but they all have an open end towards infinity in our point configuration space. This is to be expected from the general behaviour of our model since the points of agents increase indefinitely and no agent loose points.

The theory of graphic arrangements have close relations to graph colouring problem in that the defining polynomial of the arrangement is given by the chromatic polynomial of the graph. The chromatic polynomial $\chi_{G}(q)$ of a graph counts the number of ways a graph $G$ can be coloured with $q$ colours. In our context, colour means the asymptotic point gain speed of the agent in view of Eq.31, for $p>1/2$.

For a graph $G$ with $N$ vertices, the number of chambers which is the technical term used for the bulk regions defined by the velocity field where no neighbouring agent has the same points, is given by $(-1)^{N}\chi_{G}(-1)$. Also this number is equal to the number of acyclic orientations of $G$. An orientation of a graph is obtained by assigning directions on the links of the graph that joins the vertices and the total number of such orientations is simply $2^{E}$. An acyclic orientation is one that does not create oriented cycles in the diagram and thus, in our context, defines an instance of a local absolute ordering of points of the agents residing on the vertices 444If the orientation can be cyclic it can not be associated with an ordering of points for those agents residing on nodes of the cycle in question, they will have equal points in a cyclic orientation..

Using this language the meaning of Eq.31 is clearer: For $p>1/2$ it’s solutions are colourings that are defined by an acyclic orientation such that the speed, that is the colour, is equal to the number of outgoing arrows from a vertex. Note that not all acyclic orientations will therefor provide colourings that will resolve Eq.31 and conversely not all colourings obtained using the maximum number of colours defined with $d_{\rm max}$ will provide a meaningful acyclic orientation. The constraint equation Eq.31 is non-trivial in this respect.

III.5 Roles of the non-convective terms

The non-convective terms dominate the evolution for $p<1/2$ since here the points are generally drawn near the hyperplanes via the convective flow. The dynamics in this case is a complicated pattern of agents hitting a hyperplane then diffuse out only to return back in short time via the convective flow. The discussion we have below is therefore confined to the case $p>1/2$ where the convective terms tend to push the state away from the hyperplanes. However, as we have already mentioned and as we shall discuss with explicit examples the solution to the pure convective result in Eq.31 can be used to assess the configurations of asymptotic speeds even for $p<1/2$.

There are two major roles of the non-convectice terms. First they dominate the evolution of the system for early times: In a simulation where the actual configuration space is discrete one possibly spends some time early in the simulation not in the bulk regions but on the hyperplanes dividing these regions and even on their further intersections. This of course depends on the inital conditions: If one picks initial conditions from a large cube chances are one will be deep in a well defined bulk region, but if one starts all agents from the same points early times are dominated purely by the diffusion. Second, even when the state actually falls well in a bulk region the velocity field might be parallel to a hyperplane and thus for $p>1/2$ can not create a flow diverging from that hyperplane 555 Remember that we have already shown that when two connected agents $i$ and $j$ have the same speed, that is $v_{i}=v_{j}$, the velocity field is parallel to the hyperplane $x_{i}=x_{j}$..

Now, let us take two connected agents which we label as $x_{1}$ and $x_{2}$, and let us assume we are in a bulk region with a well defined ordering of points and hence one has a constant velocity field throughout. The following co-ordinate transformation

[TABLE]

will cast the hydrodynamical form of the master equation to the following

[TABLE]

where the last two summations contain no $u$ and $w$ dependence inside this chamber.

We realize that the co-ordinate $u$ which is along $\hat{e}_{1}-\hat{e}_{2}$ performs ordinary diffusion. It is well known that the eventual probability for the walker to hit $u=0$, which is $x_{1}=x_{2}$ plane, is 1. Yet the mean time to do so is infinite. Unfortunately. this is a nuissance for actual simulations which by definition must encompass a finite time.

We can thus define two types of bulk regions: Ones that we would like to call stable and the others unstable. The stable bulk regions are where, for $p>1/2$, the velocity field is pointing away from all hyperplanes bounding that region. The unstable ones are where, again for $p>1/2$, the velocity field is parallel to, possibly more than one, hyperplane.

To exemplify this behaviour let us assume the following local piece of a graph along with a given local acyclic orientation,

$\cdots$ 2110 $\cdots$

Where the labels are the would be asymptotic speeds of the agents which for simplicity we presented for $p=1$ and normalized to the number of outgoing arrows. This example is pertinent to the ring graph which we shall study at length in this manuscript. Clearly this is not a proper colouring of this portion of the graph and in the bulk region defined with the orientation above the velocity field is parallel to hyperplanes defined with the neighbouring agents which have the same speed.

Now as it stands this configuration must evolve with the diffusing terms since as we have shown the point differences of the agents with the same speed will perform an ordinary random walk. The moment when two neighbouring agents with the same speed have equal points we shall be on the hyperplane. Thus the configuration may evolve to the following

$\cdots$ 2020 $\cdots$

which constitute a stable configuration, or may fall back to the old unstable bulk region, just shy of the hyperplane. In the long run this falling back does not change the eventual probability to cross the hyperplane, albeit with an inifinite mean time to do so.

Of course if the chain with the speed $1$ agents were longer the transitions are richer and may eventually contain a stable configuration with a number of speed $1$ agents. Consider for instance the following initial configuration

$\cdots$ 21110 $\cdots$

This will eventually evolve to either the following stable configuration

$\cdots$ 20210 $\cdots$

or to the following one

$\cdots$ 21020 $\cdots$

Thus the problem here is to find the splitting probabilities of exit to either stable region. This clearly has direct dependence on the initial location inside the unstable region; the very source of broken ergodicity we have mentioned. However if the initial conditions are such that one has equal probability to start from anywhere within the unstable region of type presented above, the splitting probabilities are distributed evenly 666For instance see redbook1 and references therein. across the stable regions that surround it.

III.6 On initial conditions and interpretation of the asymptotics

The discussion on the role of non-convective terms makes it clear that the general solution, though in principle tractable, has heavy dependence on the initial conditions. If, on top of this, we remember that the convective terms dominate the long term asypmtotics in a stable chamber for $p>1/2$ we realize that a change in the questions of interest might help. So let us ask: For $p>1/2$, given that the initial condition is evenly distributed in a particular chamber how does the asymptotic state distributes itself among the stable chambers? The answer is simple. In view of the equality of splitting probabilities if one samples evenly in a bounded region one simply has to count the number of stable chambers which bound the region where the initial state resides. This is a problem in combinatorial geometry that has powerful tools. The solution however will again depend on the particulars of a given graph, and though constituting a simpler problem than the most generic initial condition may still require formidable effort. In view of this one can even generalize slightly: Say we randomly sample initial conditions from a cube cornered at the origin or a sphere centered at the origin. When the size of these grow such that the hyperplane contributions to the volume can be neglected in a discrete lattice, every chamber will sampled in accordance with the volumes they occupy. The calculation of the relative volumes of chambers in a given graph hyperplane arrangement is fairly straightforward since they are central, but a closed form is not necessarily available for regular graphs of infinite size. The final answer to the question just asked will thus have to include these ratios as well.

The second type of initial condition that stands out is one where all agents start at the same point, which can always be chosen to be zero. In this case the early stages of a simulation is complicated since on a discrete lattice the bulk regions can not be distinguished too close to the origin; there is simply not enough points on the lattice. For this reason many of the agents will be stuck at zero points for $p>1/2$ and thus the ratio of zero speed agents will be larger than the previous type of initial conditions; we are mostly on the hyperplanes defined by the co-ordinate zeros. The asymptotic properties will thus be different.

Once these transition probabilities between chambers is known the problem is mapped to a Markov process and the final distribution can in priciple be obtained brown . Yet obtaining these transition probabilities is a formidable job in general and in this work we give explicity study of them only for a few graphs.

IV Some examples on finite graphs

In this subsection we shall study the system on some finite graphs to expose the properties we have so far established.

IV.1 The smallest graph

The smallest non-trivial graph is one that has two vertices and a single edge that connects them; note that this is a regular graph. Let us represent their points on a two dimensional plane. The visualization of the velocity field is presented in Fig.2.

The simple structure of the graph divides the configuration space into two equivalent parts. We thus infer that for $p>1/2$ no matter what the initial conditions are one of the agents will have asymptotic speed $p$ and the other $q=1-p$ this property is quite visual from the figure of the vector field. The statistics of the system is equivalent to that of $p=1$ and the probabilities of the asymptotic configurations are equally likely if one starts from the origin since the geometric structure of the bulk regions enjoys a flip symmetry about the body diagonal in the configuration space. For $p\leq 1/2$ the asymptotic speeds of the agents are the same and a shock is formed as the flow converges to the body diagonal as apparent from the figure.

IV.2 Complete graph with $N$ nodes

In a complete graph every vertex has an edge to every other: That is every agent can play a game with every other and hence this constitutes the basis for the mean field approach to an arbitrary graph. It is also a regular graph. Thus for $p\leq 1/2$ all agents will have the same asymptotic speed and a shock is formed trough the entire graph. The statistical behaviour for $p>1/2$ is identical to that of $p=1$ which has the following velocity field,

[TABLE]

where we have normalized so that the sum of its componente is $N$. Since all the agents are linked none will have the same speed and the application of the consistency relation Eq.31 will mean

[TABLE]

where we have to deny zero arguments in the step function using our argumentation on the stable and unstable chambers. The only consistent solution to the equation above, within our requirements, is that there is a unique agent with $z=0$, a unique agent with $z=1$, another unique agent with $z=2$, and finally at the end a unique agent with $z=N-1$.

Thus if we label the number of agents with points less than $s$ at time $t$ by $F(s,t)$, we know that it goes to a function of $F(s/t)$ and defining a proper normalization for $z=s/t$ such that the mean velocity of the list of agents is $1/2$ we find that $F(z)=z$ between $0\leq z\leq 1$. This result is well known since the model and some of its extensions on the complete graph are extensively studied. Here the framework we established helps us understand it as the solution to the colouring problem on the complete graph with $N$ nodes. To take a note in the literature: The fact that for $p>1/2$ all agents travel at different speeds for this system was mentioned in the past, albeit not in this form; see Fig.3 in 3pap .

We see that the constraints imposed by Eq.31 are not too spectacular for this graph, since every colouring of a complete graph is a strict ordering. Also, the hyperplane arrangement associated with this graph enjoys the full permutation symmetry and every chamber has equal relative volume when intersected with a $N$-dimensional cube or sphere with the center at the origin. This means asymptotically that every ordering is equally likely to occur if the initial condition does not upset this symmetry.

IV.3 The N=3 Tree Graph

We now consider the smallest non-regular graph given as below

acb

where the labels are those that we shall use to denote the agents residing at those vertices. The velocity field is,

[TABLE]

We start by looking for all possible asymptotic speeds. To that end we performed simulations for every value of $p$, the results are presentend in Fig.3.

We shall now show that the structure of possible speeds can be fully understood with the use of Eq.31 even for $p\leq 1/2$. First of all one can show that for this graph no partial orientation solution exists where some neighbouring agents have the same speed; more on this in just a moment. We thus turn to the complete orientations of the graph which are given below

$S_{T}\;:$acb

$S_{B}\;:$acb

$U_{L}\;:$acb

$U_{R}\;:$acb

We can easily show that the acyclic orientations called $U_{L}$ and $U_{R}$ do not provide solutions to Eq.31, except for the case $p=q=1/2$ where all regions are treated equally anyways, and this justifies calling them unstable chambers.

The acyclic orientation $S_{B}$ simply demands that one must has $z_{a}=z_{b}=p/2$ and $Z_{c}=q$ and this must be a solution to Eq.31, meaning

[TABLE]

We clearly see that the solution is valid only for $p>2/3$. If $p$ drops below this value this acyclic orientation does not provide a solution and hence disappears from the list of possibly asymptotic speeds as can be clearly seen from the behaviour of the simulations presented in Fig.3.

Coming to the case of $S_{T}$ we realize that the speeds must be $z_{a}=z_{b}=q/2$ and $z_{c}=p$. Since this must solve Eq.31, we need

[TABLE]

As can be seen the solution is valid for $p>1/3$, as can be seen from Fig.3, where there are two distinct speeds for each $p$.

What happens for $p<1/3$? Well, clearly Eq.31 can not be satisfied with all equal speeds. However we also have the constraint that $\sum_{a}z_{a}=1$. Thus below the value where there are no more solutions to Eq.31 via partial orientatons, the situation is resolved by the fact that all $z_{a}$ are the same and are equal to $1/3$. This situation is in exact analogy to the resolution of a shock in hyperbolic differential equations via the Rankine-Hugoniot conditions.

We thus realize that Eq.31 along with orientations of a given graph is giving enough information about the possible asymptotic speeds of agents. The details of course will depend heavily on the topology of the graph.

These behaviours can also be understood by plotting the velocity field of the graph which we present in Fig.4. To that end we first remind the reader that at the $k$’th step of the simulation the state vector is on the plane $a+b+c=P_{o}+k$. Thus for this graph instead of providing a three dimensional plot of the velocity field we can present their projections on the plane defined above. Since the sum of the components of the velocity field is normalized to unity everywhere the component along the normal to the plane is the same and thus one must add $\hat{\mathcal{E}}/\sqrt{3}$ to the velocities presented. One word of caution is also in order. The plane where the state vector resides changes at every step of the simulation, so in reading the graphs one must also add this into account to prevent misconceptions such as inferring that the state vector will hit the co-ordinate hyperplanes which is strictly forbidden by the microscopic rules from the outset.

IV.3.1 Initial data and asymptotics

Bulk Initial Distribution: The first initial data type we shall study is the one where we sample the initial state equally likely from a hypercube or hypersphere; thus the chambers are sampled with probabilities proportional to their volume ratio. This means that the initial state will be in chamber $S_{T}$, $S_{B}$, $U_{R}$ and $U_{L}$ with probabilities $2/6$, $2/6$, $1/6$ and $1/6$ respectively. Now let us consider the case with $p=1$ for concreteness. As we have seen in this case the velocity field in $U_{R}$ and $U_{L}$ are parallel to the hyperplanes that separate them with $S_{T}$ and hence any state in these chambers that are unstable will eventually make a transition to the chamber $S_{T}$ and remain there. This means that in the asymptotic limit the probability that the state will be in chambers $S_{T}$, $S_{B}$, $U_{R}$ and $U_{L}$ are given as $2/3$, $1/3$, [math] and [math] respectively; a result that we have substantiated with numerical simulations. As another example for this initial conditions let us consider the case where $1/3<p<2/3$. It is clear from the plots of the velocity field that in this case only $S_{T}$ is a stable chamber and all the initial states will converge here; a fact we have also substansiated with simulations.

Zero Initial Distribution: The second type of initial conditions is the one where all agents start with the same point which can be chosen to be zero without loss of generality. To keep the analysis short we shall study this only for the case of $p=1$. The most important difference of this type of initial conditions is that there will always be at least one agent with zero point and hence the state of the system travels not in the chambers that constitute the bulk but is always on its edges defined by the co-ordinate hyperplanes. We thus expect rather different behaviour as opposed to the previous case. This is, as we have stressed before, due to the broken ergodicity inherent to the model.

One can, by direct enumeration, show that if we start every agent with zero points, the asymptotic state will be in the chamber $S_{T}$ with probability $3/4$ and in the chamber $S_{B}$ with probability $1/4$, and it is never in the chambers $U_{L}$ or $U_{R}$. On the other hand one can also study the infinite time behaviour by approaching the system as a Markov process with a transition matrix between relevant locations of the state vector; which is of course an equivalent method. We present this as a graph rather than as a matrix in the same presentation we used for plotting the velocity field. From the graph one can read out the transition matrix and take its infinite power which will duplicate the results just mentioned. All of these results are substantiated by direct numerical simulations.

IV.4 General properties of the solutions to Eq.31

In this section we discuss general properties of the speed constraint equation on the configurations. The listing of all possible solutions for all $p$ depends on the particulars of the graph. However the following algorithm will provide all the possible configurations.

First we remember that for $p=1/2$ all speeds are strictly given as $z_{i}=d_{i}/2E$. For $p>1/2$, remembering that the constraint depends on step functions and hence an ordering of points we must check all the acyclic orientations of the graph. For $p=1$ some of these will provide solutions and some will not and as we run $p$ to lower values some will fail to be solutions an will disappear from the list. Since the constraint equation is linear in $p$ all the curves are linear as in Fig.3 and they meet with the solutions of $p=1/2$. For $p<1/2$ we look for all orientations including the ones that have equality of neighbouring speeds and even include cyclic orientations which will mandate the equality of speeds of vertices in the cycle; which we have dubbed before as partial orientations. Again some of these will provide solutions and some will not. For example one can show that for the path graph with 4 vertices the orientation $z_{1}<z_{2}=z_{3}>z_{4}$ will provide a solution for $1/4\leq p\leq 1/2$. As one lowers $p$ there may be a critical point when no partial orientations provide any solutions and thus there can be only one solution: All agents travel at the same speed given as $1/N$, which comes from the normalization of the velocity vector.

What is the criterion for the case where all agents travel at the same speed? We have already shown that for regular graphs all agents travel at the same speed for $p\leq 1/2$. But for other cases let us consider the agent that has the minimum number of neighbours. Clearly these agents have the minimum speed at $p=1/2$ and as one decreases $p$ their speed must increase. If all the agents travel at the same speed below a critical value $p_{c}$ this speed must be given both as $1/N$ and as the speed of this lowest speed agent which is clearly $q_{c}d_{\rm min}/E$. So we get the condition

[TABLE]

What this condition means is that if it is not satisfied from above some agents will still have different speeds even for $p=0$. This behaviour was already hinted by the bucket-pipe analogy we have discussed.

Note that the above condition is validated for regular graphs where $Nd=2E$ and hence $q_{c}=1/2$. We also realize that for tree graphs where $d_{\min}=1$ and $E=N-1$ all agents travel at the same speed below $q_{c}=(N-1)/N$. This can also be very clearly understood from the bucket-pipe analogy we have mentioned before: The the amount of water disperses equally on a tree graph for $p=0$.

For graphs that are not tree yet are still planar we can use the known result $N-E+C=1$ where $C$ is the number of independent cycles in the graph; we shall get

[TABLE]

For instance if $d_{\rm min}=1$ for a planar graph; for $C=0$ all agents travel at the same speed below $q_{C}=(N-1)/N$ as the graph will be a tree; for $C=1$ all agents travel at the same speed at $q_{c}=1$ and for $C>1$ we realize that there will be more than one possible speed at $p=0$.

If $d_{\rm min}=2$ the graph can not be a tree and if planar must have $C\geq 1$, this saturates the lower bound and we shall obtain $C<N+1$ as the condition for all agents travelling at the same speed below $q_{c}$.

V The ring graph

The ring graph is a regular planar graph with a single cycle and hence $E=N$. For $p\leq 1/2$ all agents will travel with the same asymptotic speed. As we have stated for regular graphs the statistical properties $p>1/2$ cases are identical and from now on we shall take $p=1$. Let us remember that in this regime neighbouring agents will have different asymptotic speeds.

The condition on the asymptotic speeds will now read

[TABLE]

Using the scaling invariance of the step functions we shall redefine these speeds so that one will have

[TABLE]

so that $\sum_{i}z_{i}=N$ and the graph average of the speeds become $1$. This means that all $z_{i}$ are in the set $\left\{2,1,0\right\}$, which we would like to relabel as $\{L,M,S\}$, where $L$, $M$ and $S$ refer to largest, middle and smallest of possible speeds respectively.

Furthermore since no neighbouring agents will have the same $z$ value we can ignore zero argument in the $\theta$ functions and as a result of this we have the following conditions

[TABLE]

for any $a>0$.

V.1 Constraints on the configurations

So the model of two agent games on the ring looks like a 3-colouring problem. This is a well known system and is directly tied to the 3-state anti-ferromagnetic Potts model on the ring graph at zero temperature. However this analogy is altered by constraints that naturally arise form the two-agent model via Eq.46 and by the broken ergodicity inherent in the model.

For instance consider the following configuration

[TABLE]

This does not make any sense since the agent in the middle will surely win all games against its neighbours. Assuming the rest of the configuration is meaningful in this respect, this configuration should evolve to

[TABLE]

Via a similar reasoning one has

[TABLE]

Thus the configuration space of our model actually is that of a 3-state anti-ferromagnetic Potts model at zero temperature where the colors are labeled as $\{L,M,S\}$ and where the tilings $SMS$ and $LML$ are forbidden. These conditions are symmetric with respect to $L\leftrightarrow S$ and thus the full symmetry group of the configurations for a given $N$ are $G=D_{N}\times C_{2}$ as opposed to $D_{N}\times S_{3}$ which binds the unconstrained three colourings of the cycle graph. Here $C$, $D$ and $S$ mean the cyclic, dihedral and permutation groups respectively.

As the constraint equation on the speeds involve the orientations of the graph we can also see the imposibility of a tile like $LML$ or $SMS$ as follows: M-type agents mean one has $\cdots>M>\cdots$ type of ordering. Similarly L-type agents have $\cdots<L>\cdots$ and S-type agents have $\cdots>S<\cdots$ and thus no orientation can provide the mentioned tiles.

To express these observations in mathematical terms, let us define the following quantities

[TABLE]

which are simply the eigenfunctions of $z_{i}$ in that they obey

[TABLE]

These variables further satisfy $L_{i}^{2}=L_{i}$, $M_{i}^{2}=M_{i}$, $S_{i}^{2}=S_{i}$ and $S_{i}M_{i}=S_{i}L_{i}=M_{i}L_{i}=0$. Since they are exclusive one also has $L_{i}+M_{i}+S_{i}=1$. Furthermore

[TABLE]

which simply re-express that this is a colouring problem.

Now it is easy to show that

[TABLE]

as the mathematical expressions of the constraint on the configurations.

From the constraints it is clear that the M-type agents are residing over a checkerboard formed by $L$ and $S$ type agents. Note that there will be no alteration on the checkerboard in view of the constraint mentioned above. Since for an alteration to occur a tile of the forbidden type must be used. As an example this is a valid tiling where there is no alteration of the checkerboard background $\cdots LSLSMLSMLMSMLSLSLSLSMLSLM\cdots$.

Furthermore once a given configuration is found the only affect of $L\leftrightarrow S$ is to yield just a single other configuration. From this perspective one can, when convenient, see the background only as a single type of object say $o$ and the above configuration will become $\cdots ooooMooMoMoMooooooooMoooM\cdots$. Also, as a consequence of these constraints if the ring has even(odd) number of vertices there can only be an even(odd) number of M-type vertices in any configuration.

On the other hand it is also fruitful to remember the connection to hyperplane arrangements and hence to the generic vertex-model character of the system we study. Any valid configuration described above can be mapped one-to-one to an acyclic orientation on the ring graph and it must be a stable orientation from previous arguments we have presented. Not all acyclic orientations however are goind to yield a colouring commensurable with the constraints. Yet we have already shown that the exit from an unstable chamber rests on first passage times to the neighbouring stable chambers in the configuration space. For the ring graph one can show that the unstable chambers are those that contain at least three similar signs in a row. For instance the following acyclic orientation $\cdots<><<><<<<>>>>><><\cdots$ will be mapped to an unacceptable colouring in the form $\cdots SLMSLMMMSMMMMLSL\cdots$. As we have discussed before, the resolution of this will be via the first passage properties of the series of like arrow regions in the orientation. But if we sample equally likely from the acyclic orientations the splitting probabilities will be evenly distributed and thus we simply flip orientation of the links in the chain randomly. Note that the change of the orientation of a single edge in a series of like arrows will remove two M-type agents and will replace it with a checkerboard element -LS or SL- congruent with the rules. Thus if we sample equally likely from acyclic orientations the dynamics is such that given a $k$ chain of like signs one performs sign flips in a region of length $k-2$ and this is like dimer adsorbtion in a chain of length $k-1$. However, choosing initial conditions equally likely from acycic orientations is a bit unnatural for the microscopic model we have. The most natural ones are either starting all from zero points or sampling equally likely from inside a large cube. The latter will sample from acyclic orientations commensurable with their relative volumes. We shall come back to these issues in a moment.

Going back to the study of constraints the following relations are easily obtained

[TABLE]

along with

[TABLE]

and one can iterate this last property to find

[TABLE]

which simply means that the probability distribution of a single speed variable $z_{i}$ depend on a single parameter.

Last of all we would like to talk about an exact chain of relations between expectation values of the different tilings. As well known these are obtained by first starting with a single site’s configurations: $L_{i}$, $S_{i}$ or $M_{i}$. At the zeroth level we have from the probability conservation the following

[TABLE]

where we have used the symmetry group $G$ to reduce the number of unkowns.

Now by multiplying by either $1_{i\pm 1}=L_{i\pm 1}+S_{i\pm 1}+M_{i\pm 1}$ and using all the $G$ symmetry we shall arrive at the following set

[TABLE]

Thus if $\langle M\rangle$ is known so are the two tilings. However if we proceed further, the symmetry group $G$ is not enough to reduce the unknowns to the level of the distinct equations one can obtain from this procedure and one will have to introduce at each level new parameters. One can show that these can be chosen to be $\langle M\rangle$, $\langle MoM\rangle$, $\langle MooM\rangle$, and so on. Clearly this amounts to giving the chemical potential and interaction energy if one views the $M$ type agents as particles of a gas over an $LS$ checkerboard background.

V.2 The speed distribution of agents in the mean

An important property of the model is the number of agents having points less than say $s$. This information can be obtained by studying the following random variable

[TABLE]

One can use the representation

[TABLE]

perform a series expansion, swithch to the scaling ansatz $z_{i}=x_{i}/t$ and use Eq.54 to arrive at

[TABLE]

where we have defined $z=s/t$. Thus we see that the fraction of agents with speed [math] is equal to the fraction of agents with speed $2$; somewhat expected from the $L\leftrightarrow S$ symmetry of the constraints we have found. There remains to find $n_{L}$ or $n_{M}$ for that matter since $n_{M}+2n_{L}=1$.

V.2.1 Various Initial Conditions and Simulation Results

Starting from an almost cyclic orientation: Here the initial condition is such that the points of the agents are large and follow the ordering given by the orientation $\cdots>><>>>\cdots$, that is we have $N-1$ same side orientations and one opposite. If the initial points are distributed evenly inside the chamber defined by this orientation the splitting probabilities to the neighbouring chambers are evenly distributed as well. The problem therefor is that of a dimer adsorption in a path graph of length $N-1$. In the limit of $N\to\infty$ we shall therefore get flory ; redbook2

[TABLE]

which is in full accord with our simulations.

Starting from acyclic orientations sampled equally: Here we simply sample from acyclic orientations equally likely. A given acyclic orientation has a known distribution $p_{k}$ of k consecutive like signs, which in the limit $N\to\infty$ have a simple expression. Since there will be various chains of like signs, the eventual ratio of $M$-type agents will be given as

[TABLE]

where $g(k)$ is the average number of unoccupied sites, which are M-type agents in our context, in the dimer adsorbtion problem. They are given as redbook2

[TABLE]

with

[TABLE]

Also in the large $N$ limit $p_{k}$ is given as

[TABLE]

By performing the sum above with Mathematica for $N=1000$ we find

[TABLE]

which is in full accord with our simulations.

Bulk type initial conditions: This is the first of a realistic initial condition we consider: We sample points form a large $N$ dimensional cube, which in turn yields an acyclic orientation sampled by its relative volume. The relative volumes of chambers are easy to define but we do not know of any closed form formula for the ring graph. However the ratio of $M$ type agents can be expressed using the result on equal sampling of acyclic orientations modulated by the mentioned volume ratios and for small rings an exact calculation is possible. The simulation results are

[TABLE]

Zero initial condition: In this case all agents start with the same point which, in view of the mentioned shift symmetry of the velocity field, can be taken to be zero without loss of generality. The simulation results are

[TABLE]

V.3 Two meaningful approximations

In this subsection we shall provide two approximations schemes. One for the estimation of the ratio of number of agents having types of speeds pertinent to the graph and one for the statistics of the configurations.

V.3.1 Congruent minimal tiles approximation for the density of types of vertices

One can show that the following exact relations hold

[TABLE]

Note that on the first equation there are an equal number of M-type vertices on both sides of the equation. Similarly on the second(third) equation there are an equal number of L(S)-type vertices on both sides. Furthermore the 3-tiles described above are the smallest tiles that are congruent with all the constraints of the system. As this equation is exact it will also hold for the configuration averages of the tiles that takes part in it.

The approximation is to assume that the configuration averages of the 1-tiles on the left are proportional to the number of 3-tiles to the right. This does not mean that we expect all the 3-tiles on the right hand sides occur equally likely, this does not happen as we have checked in simulations. We simply assume that the sums contrive to yield a number proportional to the number of 3-tiles on the right hand sides. The normalization of the probability does the rest and we have

[TABLE]

Note that this approximation is quite good for the bulk type initial conditions and still close for zero initial conditions. We shall see that this scheme provides a rather adequate estimate for higher dimensional tori as well.

V.3.2 Chemical potential approximation for the statistics of the configurations

Due to the broken ergodicity inherent in the model not all configurations occur with the same probability. As we have previously mentioned one can in principle calculate them from the transition probabilities from unstable chambers to the stable ones. However this is a formidable task as we have laid out for the tree graph with three vertices and we do not yet have an approach that is sure to give a closed form result.

However it is also clear that the most important quantity is the ratio of M-type agents. The chemical potential approach we advocate is to assume that all configurations with a given number of M-type agents occur equally likely weighed with a chemical potential $\mu$. This potential describes the cost of removing a pair of LS agents from the full checkerboard configuration of unspecified absolute probability and replace them with two M-type agents randomly distributed over the background such that they are not next to each other.

For instance in this approach the ratio of M-type agents will be given as

[TABLE]

where for simplicity we assumed that $N$ is even and the factors of $2$ in the numerator and the denominator are simply there because of $L\leftrightarrow S$ symmetry; as we see they cancel.

The expression above has a well defined large $N$ limit given as

[TABLE]

What one can do now is to associate a particular value of $\mu$ for a given type of initial condition, from the data or from the CMT approximation. The rest will be to calculate various statistical quantities from the partition function. In the appendix we give full details of this approach.

Of course this is only an approximation and can be further ameliorated at the next level by the inclusion of a new separate chemical potential for tiles of the type $MoM$ where $o$ stands for L or S-type vertices. However we do not pursue it further here; the chemical potential approach alone provides a rather good approximation.

V.4 Separation distribution of two M type agents

One important quantity we can study is the distance distribution between two M-type agents. We know that the space between them is filled with an island of checkerboard configuration of L and S type agents. Since we are studying the ring a configuration, $N_{M}$ M-type agents will give rise to $N_{M}$ checkerboard islands. Since the total number of agents is conserved we arrive at

[TABLE]

Since when we let $N\to\infty$ so will $N_{M}$, dividing by $N_{M}$ we arrive at an exact result

[TABLE]

To estimate the distribution using the chemical potential approach we first recall that in this approach configurations with a given number of M-type agents are sampled equally likely. This means the following: Let us assume we walk on a given configuration given that we start on an M-type agent and count the number of sites to reach the next M-type agent. The mean number of these agents is given as $n_{M}$ and since in the chemical potential approximation the configurations do not care where the M-type particles reside, we expect this to be a Poisson process. As this is a conditional probability we thus have for distribution of distances

[TABLE]

Calculating $\langle\delta\rangle$ with this distribution gives us the value of $\alpha$

[TABLE]

We derive this expected result in the Appendix. The comparison with actual data and this estimate is presented in Fig.6.

V.5 The pair correlation function of asymptotic speeds

Another quantity of interest is the pair correlation function

[TABLE]

where we have used the $G$ symmetry to reduce the form of the right hand side. Note that since $\langle z_{i}\rangle=\langle z_{i+k}\rangle=1$, to obtain the pair cumulant one has to simply subtract 1 from the above equation.

For $k=0$ we use the previous result that $z_{i}^{2}=z_{i}+2L_{i}$ and hence one has

[TABLE]

Furthermore for $k=1$ we have to calculate $\langle z_{i}z_{i+1}\rangle=2\langle L_{i}M_{i+1}\rangle+2\langle M_{i}L_{i+1}\rangle$. But from our constraints $L_{i}M_{i+1}=L_{i}M_{i+1}S_{i+2}$ and $M_{i}L_{i+1}=S_{i-1}M_{i}L_{i+1}$. Remembering Eq.67a we see that

[TABLE]

In the Appendix where we solve the problem with the chemical potential approach we derive an exact expression for the pair correlation.The results of the simulation and these theoretical calculations are shown in Fig.7.

VI The torus graph

The rectangular lattice with periodic boundary conditions is topologically equivalent to a torus. Since all agents have the same number of neighbours this graph is regular. We consider the case where links are placed along the orthogonal directions so every agent has 4 neighbours. Since the graph is regular, for $p\leq 1/2$ all agents have the same asymptotic speed and all cases $p>1/2$ are equivalent to the case $p=1$. For $p=1$ the constraint of the speeds are given, with an obvious normalization similar to the one we have done for the ring graph, as

[TABLE]

So the speeds are in the set $\{0,1,2,3,4\}$. Let us denote these speeds as $\{S,M^{-},M^{0},M^{+},L\}$ and the mean fraction of agents having those as $n_{S}$, $n_{M^{-}}$, $n_{M^{0}}$, $n_{M^{+}}$ and $n_{L}$ respectively.

In analogy to the ring graph case we can introduce the variables $S_{i,j}$,$M^{-}_{i,j}$, $M^{0}_{i,j}$, $M^{+}_{i,j}$ and $L_{i,j}$ which will be the eigenvectors of $z_{i.j}$. Let us collectively call them $\sigma^{\alpha}_{i,j}$ with $\alpha\in\{0,1,2,3,4\}$. These, as before, list all possible speeds of a given agent in terms of all the 5-node tilings centered about a given agent that are commensurable with the velocity equation. Similarly one can also easily show that $(\sigma^{\alpha}_{i,j})^{2}=\sigma^{\alpha}_{i,j}$ and $\sigma^{\alpha}_{i,j}\sigma^{\beta}_{i^{\prime},j^{\prime}}=\delta^{\alpha\beta}\sigma^{\alpha}_{i,j}$ if $(i^{\prime},j^{\prime})$ denote a neighbour of the node $(i,j)$. And thus we see that the phase space of the torus model is that of an anti-ferromagnetic Potts model with 5 states at zero temperature, and hence 5-colourings, constrained with the velocity equation.

One can show that constraints similar to the ring graph case can also be constructed for the torus graph via these variables. Here instead of presenting the full list we simply write down how they will emerge.

[TABLE]

To the right of the above equation we will have all the 5-tiles colourings commensurable with the velocity equation. Since this equation is exact, it will also hold in the mean and will constitute the basis for the CMT approximation as in the ring diagram.

VI.1 The speed distribution of agents in the mean

In precise analogy to the ring graph, one can, after some tedious but straightforward calculation, show that the mean number of agents having speeds less than a given value $z$ obeys the following

[TABLE]

along with

[TABLE]

VI.1.1 Simulation Results

Bulk type initial conditions: Here we sample initial conditions equally likely from an N-dimensional cube. The size of the cube is large enough so as to discriminate different chambers in the configuration space. The simulation results are

[TABLE]

Zero initial condition: Here we start all agents from zero points. The simulation results are

[TABLE]

VI.1.2 CMT approximation

As in the ring case the CMT approximation consists of counting all the 5-tiles commensurable with constraints and assuming that the vertex type in the middle has a configuration average proportional to the number of 5-tiles that has that type of agent in the middle. The counting yields the following

[TABLE]

We see that, as was the case for the ring graph, CMT approximation is very close to the results for bulk type initial conditions.

VI.2 The geometry of the configurations

A sample configuration on a torus is shown in Fig.8 We note that as in the one dimensional case, there are LS checkerboard islands separated by M-type vertices. There are various geometrical properties one can study. We focus on the distribution of the areas and of perimeters of LS checkerboard islands. To that end let us assume that we have $N_{c}$ islands in a given configuration. One has

[TABLE]

and

[TABLE]

where $R_{i}$ is the size of an island and $\partial R_{i}$ the size of its boundary. The factor of 2 in the second equation above comes from the fact that when one sums the boundaries consisting of the M-type agents over all islands, boundaries eventually will be counted twice.

We thus immediately see the following

[TABLE]

Note that for zero initial conditions one has $n_{L}\approx 1/3$ and thus $\langle R\rangle\approx\langle\partial R\rangle$. So for this case, in the mean, the checkerboard islands have the same discrete area and boundary size. Of course such regions have a wide spectrum of shapes.

We have also studied the size distribution of the island areas and their perimeters. The results are shown in Fig.9. We see clear scaling behaviour for both the island areas and the perimeter. The scaling exponents are extracted as approximately $-1.32$ for the areas and $-1.55$ for the perimeters.

One might have imagined that one could have formed an approach similar to the chemical potential approximation we have presented for the ring graph. The basis of that approximation was of course is exact tractability. Here in the torus graph our system can be mapped to a 16-vertex model but does not admit a solution for the Yang-Baxter equation and thus it is not integrable. So we do not pursue such an avenue here.

VII Digression on d-dimensional tori

It is clear that our model on d-dimensional tori will be that of $2d+1$ colourings of the graph subject to the constraint of the velocity equation. Here we present results for the mean asymptotic speed distribution of agents in the CMT approximation. The counting consists of $2d+1$ tilings

[TABLE]

where $n_{k}$ is the fraction of agents having speed $k$ where $k\in\{0,1,\cdots,2d\}$. The functions $s(k,d)$ and $g(d)$ are given as follows

[TABLE]

where $0^{0}=1$ to be used with impunity. An analysis using the Stirling’s approximation yields for instance

[TABLE]

The term $+1$ next to $d$ in the above formula is a phenomenological amelioration of the large $d$ expansion which succeeds to yield rather close results for all $d$, including $d=1$.

Since we know $n_{k}$ from it we can calculate the distribution of agents with speeds less than $z$ where the maximum speed is normalized to $1$. This function is given as

[TABLE]

Using Stirling’s approximation and integral approximations to the sums we obtain the following

[TABLE]

It is very interesting to recover the famous lead probability distribution of the 1D random walker here. However since the meaning of $F(z)$ is basically of the same nature as the lead probabilities this result is, though maybe not outright expected, understandable. In Fig.10 we present various plots pertinent to the CMT approximation.

It is a valid question up to what dimensionality this pattern of checkerboard islands will be appreciably present. It is clear that as $d$ increases so does the fraction of L and S type agents. One can come up with a critical dimensionality by assuming that $n^{c}_{L}=1/4$ which would mean that the M type agents overall has reached half of the population. Using the CMT approximation we can calculate that for $d=3$ one has $n_{L}\approx 0.265$ and for $d=4$ one has $n_{L}=0.235$. So it seems the critical dimensionality is $4$. Motivated by this estimate we have performed an analysis on the 3D torus with zero initial conditions. The tentative results are in aggrement with the CMT approximation.

VIII Conclusions and future directions

As we have so far seen two-agent games on graphs can yield a kaleidoscope of interesting phenomena. This richness we believe is due to its simplicity. The reason for this simplicity lies in the fact that the resolution of the two-agent unit of competition depends only on step functions. However even if one generalizes the rules of engagement the basic properties we have presented, will, we believe, remain intact. To back this assesment up let us consider two agents with points $x$ and $y$. The probability that the agent with point $x$ will win is

[TABLE]

which is a function that starts out at value $q$ and makes a jump to value $p$. Now let us assume we enlarge the region of transition; that is let us take a function that starts out at $q$ and gradually increases to value $p$. The impact of this generalization will only be to somewhat enlarge the thicknesses of the hyperplanes that separates regions in the configuration space. However in the large time and hence large points limit these thicknesses will surely become irrelevant. Sure, when the hyperplane thicknesses are finite the system should take longer to reach a stable asymptotic state, yet for $p=1$ that it will reach it beyond the hyperplanes is a fact. So the model with step function captures, and in fact defines, the essentials of models of this type.

Furthermore any model of this type will have a hydrodynamical limit where one deals with a conservation law with probability density $\mathcal{P}$ and current $\vec{v}\mathcal{P}$ where $\vec{v}$ will be defined by the microscopic rules of the game. Thus Eq.31 is generic to all such type models. Hence on discrete lattices the model is a generic vertex model which has a map to a spin 1/2 model where spins live on the edges.

So for this type of models there always is a connection between, linear hyperbolic equations, hyperplane (graph) arrangements and graph colourings. One is tempted to speculate that this connection to hyperbolic equations will shed further light on the other disciplines mentioned above.

The most obvious generalization that comes to mind is to study the dynamics of 3-agent games to graphs. As mentioned this model has been extensively studies only on the complete graph. A study along this line is currently under way.

Also, since it deals with constrained colourings, our model also has connections to constrained anti-ferromagnetic Potts models. The constraints can be seen as magnetic fields that affect only a given colour. There are already some work in the litterature along these lines.

On the two dimensional torus we have said that our model is actually related to a generic 16 vertex model. However as we have seen, the symmetry $n_{S}=n_{L}$, $n_{M^{+}}=n_{M^{-}}$ along with $2n_{L}+2n_{M^{+}}+n_{M^{0}}=1$ is imposed on the system and one actually has a two parameter vertex model. Since all sixteen possible orientations of the vertices are present we could not find a solution to the Yang-Baxter equation. But, this type of vertex model might be of interest in of itself. We reserve this for further study.

As for the applications of the 2-agent model on the graphs to various phenomena one can come up with many ideas. We have already advocated a possible application to the water levels of connected water systems. Another application may be to apply it to the economical relations of a group of individuals with a given connectivity pattern; not all individuals have relations with all other individuals espescially if there are hostilities and this may define edges that will never occur in the graph the individuals define. Once such a pattern is defined the rest will follow from the model.

The author thanks M. Mungan and F. Öztürk for various ideas and feedback they provided regarding this work.

Appendix The chemical potential approach for the ring graph

In this appendix we shall develop an approximation for the statistics of our model on the ring graph. We shall assume that all configurations with the same number of M-type agents are equally probable. But we shall weigh them with respect to the full LS checkerboard configuration via a chemical potential $\mu$. For instance if a configuration contains $k$ many M-type agents we shall weigh it with $\mu^{k}$.

As well known, to obtain the statistics of the ring diagram one will need to solve the path diagram first. But since we shall ultimately be considering the limit of infinite number of vertices we do not really need to solve for the ring diagram’s partition function: Using a well known trick, whenever we need to calculate the expectation value of a given finite size subconfiguration we can place it in the middle of an infinitely long path diagram.

Now let us remember that the configurations we are interested in are all 3-colourings where the colours are labelled as $L$, $S$ and $M$ respectively and where tiles of type $LML$ and $SMS$ are forbidden.

We therefore have the following partition function for the path graph

[TABLE]

where the measure $\mathcal{Q}$ is given as

[TABLE]

Summing over the first agent we get the following recursion relation

[TABLE]

where

[TABLE]

and thus $M_{N}$ is the probability that the first agent on the path has colour $M$. Summing over the first agent we get the following reursion relation

[TABLE]

Let us also define the following, since we shall eventually need them

[TABLE]

and thus $L_{N}$ and $S_{N}$ are the probabilities that the first agent has colour $L$ and $S$ respetively. Since these are exclusive we of course have $M_{N}+L_{N}+S_{N}=1$.

One can, summing the first agent out, show that $L_{N}$ and $S_{N}$ obey the following

[TABLE]

One can by direct enumaration show that the initial conditions for $L_{N}$ and $S_{N}$ are identical and thus we must have $L_{N}=S_{N}$ which means that they are equal to $(1-M_{N})/2$. This is expected since the constraints in the measure are symmetric under $L\leftrightarrow S$.

Now one can, by rearranging the equations involving $Z_{N}$ and $M_{N}$, show that they are equivalent to

[TABLE]

[TABLE]

To resolve Eq.101 we need to provide two initial conditions. But these initial conditions must be such that they cover all the relevant constraints and thus $Z_{3}$ must be in the list. With these we get the following solution

[TABLE]

with

[TABLE]

Note that $|\lambda_{+}|>|\lambda_{-}|$ and hence we have

[TABLE]

.1 Mean number of M-type agents

To get the average number of $M$-type agents we consider the probability that the agent in the middle of a path graph is of type $M$. This means we need to calculate,

[TABLE]

By summing over the agent $\sigma_{0}$ we shall end up with

[TABLE]

In the limit $N\to\infty$ the quantity $\bar{M}_{2N+1}$ must converge to $n_{M}$ of the paper and we recover

[TABLE]

as advertized before via other means.

.2 The distance distribution between two M-type agents

To study the probability distribution $P(k)$ of distances between two M-type agents we again place the region of interest in the middle of a path graph of length $2N+k+2$ with $k\geq 1$. A straightforward calculation yields

[TABLE]

which in the large $N$ limit will give us

[TABLE]

As expected and advertized in the paper via other means this probability distribution adds up to $n_{M}$

[TABLE]

Furthermore again as advertized in the paper via other means we see that $\lambda_{+}=e^{\alpha}$ with $\alpha=\ln(2n_{L}/(4n_{L}-1))$.

.3 The pair correlation function of asymptotic speeds

The quantity of interest is

[TABLE]

But due to the flip symmetry the middle two terms are equal and thus we have

[TABLE]

Furthermore since $L_{i+k+1}+M_{i+k+1}+S_{i+k+1}=1$ we can recast this as

[TABLE]

Some tedious algebra will yield the following equations

[TABLE]

The quantity $B_{LL}(k)$ is the probability that a path graph starts with an L-type agent and ends with an L-type agent. Similarly the quantity $B_{LS}$ is the probability that a path graph starts with an L-type agent and ends with an S-type agent. They both satisfy the following equation

[TABLE]

where $f_{k}$ stands either for $B_{LS}(k)Z(k)$ or for $B_{LL}(k)Z(k)$.

However as one can easily show that their initial conditions, which one can explicitly find by direct enumeration of graphs starting with $k=2$, are not the same. Therefor they are not equal, and in fact there is no a priori expectation for them to be equal anyways.

One can readily find the particular solution to be $f^{P}_{k}=Z_{k}/2$. Therefor the general solution can be cast as $f_{k}=f^{P}_{k}+g_{k}$ where $g_{k}$ obeys

[TABLE]

which can be readily solved to give

[TABLE]

where one has

[TABLE]

Note that $|\lambda_{+}|>|\tilde{\lambda}_{\pm}|$. All one now needs are the initial conditions and we have $B_{LL}(2)=0$ and $B_{LL}(3)=1/Z_{3}$ along with $B_{LS}(2)=1/Z_{2}$ and $B_{LS}(3)=\mu/Z_{3}$.

With these and taking the $N\to\infty$ limit we finally have

[TABLE]

with

[TABLE]

We have thus ahieved an exact expression for the pair correlation functions in the limit $N\to\infty$. One can furthermore show that as $k\to\infty$ one has $C_{k}\to 1$. This expected since $\langle z_{i}z_{i+k}\rangle$ is expected to factorize as $\langle z_{i}\rangle\langle z_{i+k}\rangle$ and $\langle z_{i}\rangle=1$ via the constraints on the system.