The dynamics of language evolution is one of many interdisciplinary fields to which methods and insights from statistical physics have been successfully applied (see BIBREF0 for an overview, and BIBREF1 for a specific comprehensive review).In this work we revisit the question of language coexistence. It is known that a sizeable fraction of the more than 6000 languages that are currently spoken, is in danger of becoming extinct BIBREF2, BIBREF3, BIBREF4. In pioneering work by Abrams and Strogatz BIBREF5, theoretical predictions were made to the effect that less attractive or otherwise unfavoured languages are generally doomed to extinction, when contacts between speakers of different languages become sufficiently frequent. Various subsequent investigations have corroborated this finding, emphasising that the simultaneous coexistence of competing languages is only possible in specific circumstances BIBREF6, BIBREF7, all of which share the common feature that they involve some symmetry breaking mechanism BIBREF1. A first scenario can be referred to as spatial symmetry breaking. Different competing languages may coexist in different geographical areas, because they are more or less favoured locally, despite the homogenising effects of migration and language shift BIBREF8, BIBREF9, BIBREF10. A second scenario corresponds to a more abstract internal symmetry breaking. Two or more competing languages may coexist at a given place if the populations of speakers of these languages have imbalanced dynamics BIBREF11, BIBREF12, BIBREF13. Moreover, it has been shown that a stable population of bilinguals or multilinguals also favours the coexistence of several languages BIBREF14, BIBREF15, BIBREF16.The aim of the present study is to provide a quantitative understanding of the conditions which ensure the coexistence of two or more competing languages within each of the symmetry breaking scenarios outlined above. Throughout this paper, in line with many earlier studies on the dynamics of languages BIBREF5, BIBREF7, BIBREF8, BIBREF10, BIBREF11, BIBREF12, BIBREF13, BIBREF14, BIBREF15, BIBREF16, and with an investigation of grammar acquisition BIBREF17, we describe the dynamics of the numbers of speakers of various languages by means of coupled rate equations. This approach is sometimes referred to as ecological modelling, because of its similarity with models used in theoretical ecology (see e.g. BIBREF18). From a broader perspective, systems of coupled differential equations, and especially Lotka-Volterra equations and replicator equations, are ubiquitous in game theory and in a broad range of areas in mathematical biology (see e.g. BIBREF19, BIBREF20, BIBREF21).The plan of this paper is as follows. For greater clarity, we first consider in Section SECREF2 the situation of several competing languages in a single geographic area where the population is well mixed. We address the situation where internal symmetry is broken by imbalanced population dynamics. The relevant concepts are reviewed in detail in the case of two competing languages in Section SECREF1, and the full phase diagram of the model is derived. The case of an arbitrary number $N$ of competing languages is then considered in Section SECREF11 in full generality. The special situation where the attractivenesses of the languages are equally spaced is studied in Section SECREF22, whereas Section SECREF34 is devoted to the case where attractivenesses are modelled as random variables. Section SECREF3 is devoted to the situation where coexistence is due to spatial symmetry breaking. We focus our attention onto the simple case of two languages in competition on a linear array of $M$ distinct geographic areas. Language attractivenesses vary arbitrarily along the array, whereas migrations take place only between neighbouring areas at a uniform rate $\gamma $. A uniform consensus is reached at high migration rate, where the same language survives everywhere. This general result is demonstrated in detail for two geographic areas (Section SECREF57), and generalised to an arbitrary number $M$ of areas (Section SECREF67). The cases of ordered and random attractiveness profiles are investigated in Sections SECREF71 and SECREF84. In Section SECREF4 we present a non-technical discussion of our findings and their implications. Two appendices contain technical details about the regime of a large number of competing languages in a single geographic area (Appendix SECREF5) and about stability matrices and their spectra (Appendix SECREF6).This section is devoted to the dynamics of languages in a single geographic area. As mentioned above, it has been shown that two or more competing languages may coexist only if the populations of speakers of these languages have imbalanced dynamics BIBREF11, BIBREF12, BIBREF13. Our goal is to make these conditions more explicit and to provide a quantitative understanding of them.We begin with the case of two competing languages. We assume that language 1 is more favoured than language 2. Throughout this work we neglect the effect of bilingualism, so that at any given time $t$ each individual speaks a single well-defined language. Let $X_1(t)$ and $X_2(t)$ denote the numbers of speakers of each language at time $t$, so that $X(t)=X_1(t)+X_2(t)$ is the total population of the area under consideration.The dynamics of the model is defined by the coupled rate equationsThe above equations are an example of Lotka-Volterra equations (see e.g. BIBREF18, BIBREF19). The terms underlined by braces describe the intrinsic dynamics of the numbers of speakers of each language. For the sake of simplicity we have chosen the well-known linear-minus-bilinear or `logistic' form which dates back to Lotka BIBREF22 and is still commonly used in population dynamics. The linear term describes population growth, whereas the quadratic terms represent a saturation mechanism.The main novelty of our approach is the introduction of the parameter $q$ in the saturation terms. This imbalance parameter is responsible for the internal symmetry breaking leading to language coexistence. It allows for the interpolation between two situations: when the saturation mechanism only involves the total population, i.e., $q=1$, and when the saturation mechanism acts separately on the populations of speakers of each language, $q=0$, which is the situation considered by Pinasco and Romanelli BIBREF11. Generic values of $q$ correspond to tunably imbalanced dynamics.The last term in each of equations (DISPLAY_FORM2), () describes the language shift consisting of the conversions of single individuals from the less favoured language 2 to the more favoured language 1. In line with earlier studies BIBREF7, BIBREF11, BIBREF12, BIBREF13, conversions are triggered by binary interactions between individuals, so that the frequency of conversions is proportional to the product $X_1(t)X_2(t)$. The reduced conversion rate $C$ measures the difference of attractivenesses between the two languages.For generic values of the parameters $q$ and $C$, the rate equations (DISPLAY_FORM2), () admit a unique stable fixed point. The dynamics converges exponentially fast to the corresponding stationary state, irrespective of initial conditions. There are two possible kinds of stationary states:I. Consensus.The solutiondescribes a consensus state where the unfavoured language 2 is extinct. The inverse relaxation times describing convergence toward the latter state are the opposites of the eigenvalues of the stability matrix associated with equations (DISPLAY_FORM2), (). The reader is referred to Appendix SECREF131 for details. These inverse relaxation times readThe above stationary solution is thus stable whenever $q+C>1$.II. Coexistence.The solutiondescribes a coexistence state where both languages survive forever. This stationary solution exists whenever $q+C<1$. It is always stable, as the inverse relaxation times readFigure FIGREF9 shows the phase diagram of the model in the $q$–$C$ plane. There is a possibility of language coexistence only for $q<1$. The vertical axis ($q=0$) corresponds to the model considered by Pinasco and Romanelli BIBREF11, where the coexistence phase is maximal and extends up to $C=1$. As the parameter $q$ is increased, the coexistence phase shrinks until it disappears at the point $q=1$, corresponding to the balanced dynamics where the saturation mechanism involves the total population.The model exhibits a continuous transition along the phase boundary between both phases ($q+C=1$). The number $X_2$ of speakers of the unfavoured language vanishes linearly as the phase boundary is approached from the coexistence phase (see (DISPLAY_FORM7)), whereas the relaxation time $1/\omega _2$ diverges linearly as the phase boundary is approached from both sides (see (DISPLAY_FORM5) and (DISPLAY_FORM8)).For parameters along the phase boundary ($q+C=1$), the less attractive language still becomes extinct, albeit very slowly. Equations (DISPLAY_FORM2), () here yield the power-law relaxation lawsirrespective of initial conditions.The above setting can be extended to the case of an arbitrary number $N$ of competing languages in a given area. Languages, numbered $i=1,\dots ,N$, are more or less favoured, depending on their attractivenesses $A_i$. The latter quantities are assumed to be quenched, i.e., fixed once for all. This non-trivial static profile of attractivenesses is responsible for conversions of single individuals from less attractive to more attractive languages.Let $X(t)$ be the total population of the area under consideration at time $t$, and $X_i(t)$ be the number of speakers of language number $i=1,\dots ,N$. The dynamics of the model are defined by the rate equationsThe terms underlined by braces describe the intrinsic dynamics of the numbers of speakers of each language. The novel feature here is again the presence of the parameter $q$, which is responsible for imbalanced dynamics, allowing thus the possibility of language coexistence. The last term in (DISPLAY_FORM12) describes the conversions of single individuals. If language $i$ is more attractive than language $j$, there is a net positive conversion rate $C_{ji}=-C_{ij}$ from language $j$ to language $i$. For the sake of simplicity, we assume that these conversion rates depend linearly on the differences of attractivenesses between departure and target languages, i.e.,in some consistent units.Throughout this work we shall not pay any attention to the evolution of the whole population $X(t)$. We therefore reformulate the model in terms of the fractionsof speakers of the various languages, which sum up to unity:The reduction to be derived below is quite natural in the present setting. It provides an example of the reduction of Lotka-Volterra equations to replicator equations, proposed in BIBREF23 (see also BIBREF19, BIBREF20, BIBREF21). In the present situation, for $q<1$, which is precisely the range of $q$ where there is a possibility of language coexistence, the dynamics of the fractions $x_i(t)$ obeys the following reduced rate equations, which can be derived from (DISPLAY_FORM12):withand where attractivenesses and conversion rates have been rescaled according toIn the following, we focus our attention onto the stationary states of the model, rather than on its dynamics. It is therefore legitimate to redefine time according toso that equations (DISPLAY_FORM16) simplify toThe rate equations (DISPLAY_FORM20) for the fractions of speakers of the $N$ competing languages will be the starting point of further developments. The quantity $Z(t)$ can be alternatively viewed as a dynamical Lagrange multiplier ensuring that the dynamics conserves the sum rule (DISPLAY_FORM15). The above equations belong to the class of replicator equations (see e.g. BIBREF19, BIBREF20, BIBREF21). Extensive studies of the dynamics of this class of equations have been made in mathematical biology, where the main focus has been on systematic classifications of fixed points and bifurcations in low-dimensional cases BIBREF23, BIBREF24, BIBREF25, BIBREF26, BIBREF27, BIBREF28.From now on, we focus on the stationary state of the model for arbitrarily high values of the number $N$ of competing languages. The analysis of this goes as follows. The stationary values $x_i$ of the fractions of speakers are such that the right-hand sides of (DISPLAY_FORM20) vanish. For each language number $i$, there are two possibilities: either $x_i=0$, i.e., language $i$ gets extinct, or $x_i>0$, i.e., language $i$ survives forever. The non-zero fractions $x_i$ of speakers of surviving languages obey the coupled linear equationswhere the parameter $Z$ is determined by expressing that the sum rule (DISPLAY_FORM15) holds in the stationary state. For generic values of model parameters, there is a unique stationary state, and the system relaxes exponentially fast to the latter, irrespective of its initial conditions. The uniqueness of the attractor is characteristic of the specific form of the rate equations (DISPLAY_FORM20), (DISPLAY_FORM21), with skew-symmetric conversion rates $c_{ij}$ (see (DISPLAY_FORM18)). This has been demonstrated explicitly in the case of two competing languages, studied in detail in Section SECREF1. The problem is however more subtle than it seems at first sight, as the number $K$ of surviving languages depends on model parameters in a non-trivial way.It is useful to consider first the simple case where the (reduced) attractivenesses $a_i$ of the $N$ competing languages are equally spaced between 0 and some maximal value that we denote by $2g$. Numbering languages in order of decreasing attractivenesses, so that language 1 is the most attractive and language $N$ the least attractive, this readsWe haveThe parameter $g$ is therefore the mean attractiveness.The (reduced) conversion rates readso that the fixed-point equations (DISPLAY_FORM21) take the formAlready in this simple situation the number $K$ of surviving languages depends on the mean attractiveness $g$ in a non-trivial way.Consider first the situation where all languages survive ($K=N$). This is certainly true for $g=0$, where there are no conversions, so that the solution is simply $x_i=1/N$. There, all languages are indeed equally popular, as nothing distinguishes them. More generally, as long as all languages survive, the stationary solution obeying (DISPLAY_FORM26) readsfor $i=1,\dots ,N$. The above solution ceases to hold when the fraction of speakers of the least attractive language vanishes, i.e., $x_N=0$. This first extinction takes place for the threshold valueof the mean attractiveness $g$.Consider now the general case where only $K$ among the $N$ languages survive. These are necessarily the $K$ most attractive ones, shown as red symbols in Figure FIGREF29.In this situation, (DISPLAY_FORM26) yieldsfor $i=1,\dots ,K$. The linear relationship between the attractiveness $a_i$ of language $i$ and the stationary fraction $x_i$ of speakers of that language, observed in (DISPLAY_FORM27) and (DISPLAY_FORM30), is a general feature of the model (see Section SECREF34). The fraction $x_K$ of speakers of the least attractive of the surviving languages vanishes at the following threshold mean attractiveness:for $K=2,\dots ,N$.The following picture therefore emerges for the stationary state of $N$ competing languages with equally spaced attractivenesses. The number $K$ of surviving languages decreases as a function of the mean attractiveness $g$, from $K=N$ (all languages survive) near $g=0$ to $K=1$ (consensus) as very large $g$. Less attractive languages become extinct one by one as every single one of the thresholds (DISPLAY_FORM31) is traversed, so thatFigure FIGREF33 illustrates this picture for 5 competing languages. In each of the sectors defined in (DISPLAY_FORM32), the stationary fractions $x_i$ of speakers of the surviving languages are given by (DISPLAY_FORM30). They depend continuously on the mean attractiveness $g$, even though they are given by different expressions in different sectors. In particular, $x_i$ is flat, i.e., independent of $g$, in the sector where $K=2i-1$. The fraction $x_1$ of speakers of the most attractive language grows monotonically as a function of $g$, whereas all the other fractions of speakers eventually go to zero.When the number of languages $N$ is large, the range of values of $g$ where the successive transitions take place is very broad. The threshold at which a consensus is reached, $g_{N,2}=N/2$, is indeed much larger than the threshold at which the least attractive language disappears, $g_{N,N}=1/(N-1)$. The ratio between these two extreme thresholds reads $N(N-1)/2$.We now turn to the general case of $N$ competing languages with arbitrary reduced attractivenesses $a_i$. Throughout the following, languages are numbered in order of decreasing attractivenesses, i.e.,We shall be interested mostly in the stationary state of the model. As already mentioned above, the number $K$ of surviving languages depends on model parameters in a non-trivial way. The $K$ surviving languages are always the most attractive ones (see Figure FIGREF29). The fractions $x_i$ of speakers of those languages, obeying the fixed-point equations (DISPLAY_FORM21), can be written in full generality asfor $i=1,\dots ,K$, withThe existence of an explicit expression (DISPLAY_FORM36) for the solution of the fixed-point equations (DISPLAY_FORM21) in full generality is a consequence of their simple linear-minus-bilinear form, which also ensures the uniqueness of the attractor.The number $K$ of surviving languages is the largest such that the solution (DISPLAY_FORM36) obeys $x_i>0$ for $i=1,\dots ,K$. Equivalently, $K$ is the largest integer in $1,\dots ,N$ such thatEvery single one of the differences involved in the sum is positive, so that:From now on, we model attractivenesses as independent random variables. More precisely, we setwhere $w$ is the mean attractiveness, and the rescaled attractivenesses $\xi _i$ are positive random variables drawn from some continuous distribution $f(\xi )$ such that $\left\langle \xi \right\rangle =1$. For any given instance of the model, i.e., any draw of the $N$ random variables $\lbrace \xi _i\rbrace $, languages are renumbered in order of decreasing attractivenesses (see (DISPLAY_FORM35)).For concreteness we assume that $f(0)$ is non-vanishing and that $f(\xi )$ falls off more rapidly than $1/\xi ^3$ at large $\xi $. These hypotheses respectively imply that small values of $\xi $ are allowed with non-negligible probability and ensure the convergence of the second moment $\left\langle \xi ^2\right\rangle =1+\sigma ^2$, where $\sigma ^2$ is the variance of $\xi $.Some quantities of interest can be expressed in closed form for all language numbers $N$. One example is the consensus probability ${\cal P}$, defined as the probability of reaching consensus, i.e., of having $K=1$ (see (DISPLAY_FORM39)). This readsWe havefor all $N\ge 2$, whereis the cumulative distribution of $\xi $.In forthcoming numerical and analytical investigations we use the following distributions:We begin our exploration of the model by looking at the dynamics of a typical instance of the model with $N=10$ languages and a uniform distribution of attractivenesses with $w=0.3$. Figure FIGREF45 shows the time-dependent fractions of speakers of all languages, obtained by solving the rate equations (DISPLAY_FORM20) numerically, with the uniform initial condition $x_i(0)=1/10$ for all $i$. In this example there are $K=6$ surviving languages. The plotted quantities are observed to converge to their stationary values given by (DISPLAY_FORM36) for $i=1,\dots ,6$, and to zero for $i=7,\dots ,10$. They are ordered as the corresponding attractivenesses at all positive times, i.e., $x_1(t)>x_2(t)>\dots >x_N(t)$. Some of the fractions however exhibit a non-monotonic evolution. This is the case for $i=5$ in the present example.Figure FIGREF48 shows the distribution $p_K$ of the number $K$ of surviving languages, for $N=10$ (top) and $N=40$ (bottom), and a uniform distribution of attractivenesses for four values of the productThis choice is motivated by the analysis of Appendix SECREF5. Each dataset is the outcome of $10^7$ draws of the attractiveness profile. The widths of the distributions $p_K$ are observed to shrink as $N$ is increased, in agreement with the expected $1/\sqrt{N}$ behavior stemming from the law of large numbers. The corresponding mean fractions $\left\langle K\right\rangle /N$ of surviving languages are shown in Table TABREF49 to converge smoothly to the asymptotic prediction (DISPLAY_FORM126), i.e.,with $1/N$ corrections.An overall picture of the dependence of the statistics of surviving languages on the mean attractiveness $w$ is provided by Figure FIGREF50, showing the mean number $\left\langle K\right\rangle $ of surviving languages against $w$, for $N=10$ and uniform and exponential attractiveness distributions. The plotted quantity decreases monotonically, starting from the value $\left\langle K\right\rangle =N$ in the absence of conversions ($w=0$), and converging to its asymptotic value $\left\langle K\right\rangle =1$ in the $w\rightarrow \infty $ limit, where consensus is reached with certainty. Its dependence on $w$ is observed to be steeper for the exponential distribution. These observations are corroborated by the asymptotic analysis of Appendix SECREF5. For the uniform distribution, (DISPLAY_FORM126) yields the scaling law $\left\langle K\right\rangle \approx (N/w)^{1/2}$. Concomitantly, the consensus probability becomes sizeable for $w\sim N$ (see (DISPLAY_FORM124)). For the exponential distribution, (DISPLAY_FORM130) yields the decay law $\left\langle K\right\rangle \approx 1/w$, irrespective of $N$, and the consensus probability is strictly independent of $N$ (see (DISPLAY_FORM128)).As mentioned in the Introduction, different competing languages may coexist in distinct geographical areas, because they are more or less favoured locally, despite the homogenising effects of migration and language shift BIBREF8, BIBREF9, BIBREF10. The aim of this section is to provide a quantitative understanding of this scenario. We continue to use the approach and the formalism of Section SECREF2. We however take the liberty of adopting slightly different notations, as both sections are entirely independent.We consider the dynamics of two competing languages in a structured territory comprising several distinct geographic areas. For definiteness, we assume that the population of each area is homogeneous. We restrict ourselves to the geometry of an array of $M$ areas, where individuals can only migrate along the links joining neighbouring areas, as shown in Figure FIGREF51. We assume for simplicity that the migration rates $\gamma $ between neighbouring areas are uniform, so that in the very long run single individuals eventually perform random walks across the territory. The relative attractivenesses of both competing languages are distributed inhomogeneously among the various areas, so that the net conversion rate $C_m$ from language 2 to language 1 depends on the area number $m$. Finally, in order to emphasise the effects of spatial inhomogeneity on their own, we simplify the model by neglecting imbalance and thus set $q=1$.Let $X_m(t)$ and $Y_m(t)$ denote the respective numbers of speakers of language 1 and of language 2 in area number $m=1,\dots ,M$ at time $t$. The dynamics of the model is defined by the coupled rate equationsThe extremal sites $m=1$ and $m=M$ have only one neighbour. The corresponding equations have to be modified accordingly. The resulting boundary conditions can be advantageously recast asand similarly for other quantities. These are known as Neumann boundary conditions.The total populations $P_m(t)=X_m(t)+Y_m(t)$ of the various areas obeyirrespective of the conversion rates $C_m$. As a consequence, in the stationary state all areas have the same population, which reads $P_m=1$ in our reduced units. The corresponding stability matrix is given in (DISPLAY_FORM137). The population profile $P_m(t)$ therefore converges exponentially fast to its uniform stationary value, with unit relaxation time ($\omega =1$).From now on we assume, for simplicity, that the total population of each area is unity in the initial state. This property is preserved by the dynamics, i.e., we have $P_m(t)=1$ for all $m$ and $t$, so that the rate equations (DISPLAY_FORM52) simplify toThe rate equations (DISPLAY_FORM55) for the fractions $X_m(t)$ of speakers of language 1 in the various areas provide another example of the broad class of replicator equations (see e.g. BIBREF19, BIBREF20, BIBREF21). The above equations are the starting point of the subsequent analysis. In the situation where language 1 is uniformly favoured or disfavoured, so that the conversion rates are constant ($C_m=C$), the above rate equations boil down to the discrete Fisher-Kolmogorov-Petrovsky-Piscounov (FKPP) equation BIBREF29, BIBREF30, which is known to exhibit traveling fronts, just as the well-known FKPP equation in the continuum BIBREF31, BIBREF32. In the present context, the focus will however be on stationary solutions on finite arrays, obeyingWe begin with the case of two geographic areas connected by a single link. The problem is simple enough to allow for an explicit exposition of its full solution. The rate equations (DISPLAY_FORM55) becomeBecause of the migration fluxes, for any non-zero $\gamma $ it is impossible for any of the languages to become extinct in one area and survive in the other one. The only possibility is that of a uniform consensus, where one and the same language survives in all areas. The consensus state where language 1 survives is described by the stationary solution $X_1=X_2=1$. The corresponding stability matrix iswhere $\mathop {{\rm diag}}(\dots )$ denotes a diagonal matrix (whose entries are listed), whereas ${\Delta }_2$ is defined in (DISPLAY_FORM135). The stability condition amounts toSimilarly, the consensus state where language 2 survives is described by the stationary solution $X_1=X_2=0$. The corresponding stability matrix isThe conditions for the latter to be stable readFigure FIGREF66 shows the phase diagram of the model in the $C_1$–$C_2$ plane for $\gamma =1$. Region I1 is the consensus phase where language 1 survives. It is larger than the quadrant where this language is everywhere favoured (i.e., $C_1$ and $C_2$ are positive), as its boundary (red curve) reads $C_1C_2+\gamma (C_1+C_2)=0$. Similarly, region I2 is the consensus phase where language 2 survives. It is larger than the quadrant where this language is everywhere favoured (i.e., $C_1$ and $C_2$ are negative), as its boundary (blue curve) reads $C_1C_2-\gamma (C_1+C_2)=0$. The regions marked IIA and IIB are coexistence phases. These phases are located symmetrically around the line $C_1+C_2=0$ (black dashed line) where none of the languages is globally favoured. There, the fractions $X_1$ and $X_2$ of speakers of language 1 in both areas vary continuously between zero on the blue curve and unity on the red one, according towithWe have thereforeall over the coexistence phases IIA and IIB. The right-hand-side equals 0 on the blue curve, 1 on the black dashed line, and 2 on the red curve.From now on we consider the general situation of $M$ geographic areas, as shown in Figure FIGREF51. The basic properties of the model can be inferred from the case of two areas, studied in section SECREF57. In full generality, because of migration fluxes, it is impossible for any of the languages to become extinct in some areas and survive in some other ones. The only possibility is that of a uniform consensus, where one and the same language survives in all areas.The consensus state where language 1 survives is described by the uniform stationary solution where $X_m=1$ for all $m=1,\dots ,M$. The corresponding stability matrix isSimilarly, the consensus state where language 2 survives corresponds to the stationary solution where $X_m=0$ for all $m=1,\dots ,M$. The corresponding stability matrix isThese expressions respectively generalise (DISPLAY_FORM59) and (DISPLAY_FORM61).If all the conversion rates $C_m$ vanish, both the above matrices read $-\gamma {\Delta }_M$, whose spectrum comprises one vanishing eigenvalue (see (DISPLAY_FORM136)). In the regime where all the conversion rates $C_m$ are small with respect to $\gamma $, perturbation theory tells us that the largest eigenvalues of ${S}_M^{(0)}$ and ${S}_M^{(1)}$ respectively read $\overline{C}$ and $-\overline{C}$, to leading order, whereWe therefore predict that the average conversion rate $\overline{C}$ determines the fate of the system in the regime where conversion rates are small with respect to $\gamma $. If language 1 is globally favoured, i.e., $\overline{C}>0$, the system reaches the consensus where language 1 survives, and vice versa.In the generic situation where the conversion rates $C_m$ are comparable to $\gamma $, their dispersion around their spatial average $\overline{C}$ broadens the spectra of the matrices ${S}_M^{(1)}$ and ${S}_M^{(0)}$. As a consequence, the condition $\overline{C}>0$ (resp. $\overline{C}<0$) is necessary, albeit not sufficient, for the consensus where language 1 (resp. language 2) survives to be stable.In the following we shall successively consider ordered attractiveness profiles in Section SECREF71 and random ones in Section SECREF84.This section is devoted to a simple situation where the attractiveness profiles of both languages are ordered spatially. More specifically, we consider the case where language 1 is favoured in the $K$ first (i.e., leftmost) areas, whereas language 2 is favoured in the $L$ last (i.e., rightmost) areas, with $K\ge L$ and $K+L=M$. For the sake of simplicity, we choose to describe this situation by conversion rates that have unit magnitude, as shown in Figure FIGREF73:The symmetric situation where $M$ is even and $K=L=M/2$, so that $\overline{C}=0$, can be viewed as a generalisation of the case of two geographic areas, studied in Section SECREF57, for $C_1+C_2=0$, i.e., along the black dashed line of Figure FIGREF66. Both languages play symmetric roles, so that no language is globally preferred, and no consensus can be reached. As a consequence, both languages survive everywhere, albeit with non-trivial spatial profiles, which can be thought of as avatars of the FKPP traveling fronts mentioned above, rendered stationary by being pinned by boundary conditions. The upper panel of Figure FIGREF76 shows the stationary fraction $X_m$ of speakers of language 1 against area number, for $M=20$ (i.e., $K=L=10$) and several $\gamma $. The abscissa $m-1/2$ is chosen in order to have a symmetric plot. As one might expect, each language is preferred in the areas where it is favoured, i.e., we have $X_m>1/2$ for $m=1,\dots ,K$, whereas $X_m<1/2$ for $m=K+1,\dots ,M$. Profiles get smoother as the migration rate $\gamma $ is increased. The width $\xi $ of the transition region is indeed expected to grow asThis scaling law is nothing but the large $\gamma $ behaviour of the exact dispersion relation(see (DISPLAY_FORM148)) between $\gamma $ and the decay rate $\mu $ such that either $X_m$ or $1-X_m$ falls off as ${\rm e}^{\pm m\mu }$, with the natural identification $\xi =1/\mu $.The asymmetric situation where $K>L$, so that $\overline{C}=(K-L)/M>0$, implying that language 1 is globally favoured, is entirely different. The system indeed reaches a consensus state where the favoured language survives, whenever the migration rate $\gamma $ exceeds some threshold $\gamma _c$. This threshold, corresponding to the consensus state becoming marginally stable, only depends on the integers $K$ and $L$. It is derived in Appendix SECREF6 and given by the largest solution of (DISPLAY_FORM153).This is illustrated in the lower panel of Figure FIGREF76, showing $X_m$ against $m-1/2$ for $K=12$ and $L=8$, and the same values of $\gamma $ as on the upper panel. The corresponding threshold reads $\gamma _c=157.265$. The whole profile shifts upwards while it broadens as $\gamma $ is increased. It tends uniformly to unity as $\gamma $ tends to $\gamma _c$, demonstrating the continuous nature of the transition where consensus is formed.The threshold migration rate $\gamma _c$ assumes a scaling form in the regime where $K$ and $L$ are large and comparable. Settingso that the excess fraction $f$ identifies with the average conversion rate $\overline{C}$, the threshold rate $\gamma _c$ grows quadratically with the system size $M$, according towhere $g(f)$ is the smallest positive solution of the implicit equationwhich is a rescaled form of (DISPLAY_FORM153).The quadratic growth law (DISPLAY_FORM78) is a consequence of the diffusive nature of migrations. The following limiting cases deserve special mention.For $f\rightarrow 0$, i.e., $K$ and $L$ relatively close to each other ($K-L\ll M$), we haveyielding to leading orderFor $f\rightarrow 1$, i.e., $L\ll K$, we have $g(f)\approx \pi /(4(1-f))$, up to exponentially small corrections, so thatThe situation considered in the lower panel of Figure FIGREF76, i.e., $M=20$, $K=12$ and $L=8$, corresponds to $f=1/5$, hence $g=0.799622814\dots $, so thatThis scaling result predicts $\gamma _c\approx 156.397$ for $M=20$, a good approximation to the exact value $\gamma _c=157.265$.We now consider the situation of randomly disordered attractiveness profiles. The conversion rates $C_m$ are modelled as independent random variables drawn from some symmetric distribution $f(C)$, such that $\left\langle C_m\right\rangle =0$ and $\left\langle C_m^2\right\rangle =w^2$.The first quantity we will focus on is the consensus probability ${\cal P}$. It is clear from a dimensional analysis of the rate equations (DISPLAY_FORM56) that ${\cal P}$ depends on the ratio $\gamma /w$, the system size $M$, and the distribution $f(C)$. Furthermore, ${\cal P}$ is expected to increase with $\gamma /w$. It can be estimated as follows in the limiting situations where $\gamma /w$ is either very small or very large.In the regime where $\gamma \ll w$ (e.g. far from the center in Figure FIGREF66), conversion effects dominate migration effects. There, a consensus where language 1 (resp. language 2) survives can only be reached if all conversion rates $C_m$ are positive (resp. negative). The total consensus probability thus scales asConsensus is therefore highly improbable in this regime. In other words, coexistence of both languages is overwhelmingly the rule.In the opposite regime where $\gamma \gg w$ (e.g. in the vicinity of the center in Figure FIGREF66), migration effects dominate conversion effects. There, we have seen in Section SECREF67 that the average conversion rate defined in (DISPLAY_FORM70) essentially determines the fate of the system. If language 1 is globally favoured, i.e., $\overline{C}>0$, then the system reaches the uniform consensus where language 1 survives, and vice versa. Coexistence is therefore rare in this regime, as it requires $\overline{C}$ to be atypically small. The probability ${\cal Q}$ for this to occur, to be identified with $1-{\cal P}$, has been given a precise definition in Appendix SECREF6 by means of the expansion (DISPLAY_FORM143) of $D_M=\det {S}_M^{(1)}$ as a power series in the $C_m$, and estimated within a simplified Gaussian setting. In spite of the heuristic character of its derivation, the resulting estimate (DISPLAY_FORM147) demonstrates that the consensus probability scales asall over the regime where the ratio $\gamma /w$ and the system size $M$ are both large. Furthermore, taking (DISPLAY_FORM147) literally, we obtain the following heuristic prediction for the finite-size scaling function:The scaling result (DISPLAY_FORM86) shows that the scale of the migration rate $\gamma $ which is relevant to describe the consensus probability for a typical disordered profile of attractivenesses readsThis estimate grows less rapidly with $M$ than the corresponding threshold for ordered profiles, which obeys a quadratic growth law (see (DISPLAY_FORM78)). The exponent $3/2$ of the scaling law (DISPLAY_FORM88) can be put in perspective with the anomalous scaling of the localisation length in one-dimensional Anderson localisation near band edges. There is indeed a formal analogy between the stability matrices of the present problem and the Hamiltonian of a tight-binding electron in a disordered potential, with the random conversion rates $C_m$ replacing the disordered on-site energies. For the tight-binding problem, the localisation length is known to diverge as $\xi \sim 1/w^2$ in the bulk of the spectrum, albeit only as $\xi \sim 1/w^{2/3}$ in the vicinity of band edges BIBREF33, BIBREF34, BIBREF35, BIBREF36, BIBREF37. Replacing $\xi $ by the system size $M$ and remembering that $w$ stands for $w/\gamma $, we recover (DISPLAY_FORM88). The exponent $3/2$ is therefore nothing but the inverse of the exponent $2/3$ of anomalous band-edge localisation.Figure FIGREF89 shows a finite-size scaling plot of the consensus probability ${\cal P}$ against $x=\gamma /M^{3/2}$. Data correspond to arrays of length $M=20$ with uniform and Gaussian distributions of conversion rates with $w=1$. Each data point is the outcome of $10^6$ independent realisations. The thin black curve is a guide to the eye, suggesting that the finite-size scaling function $\Phi $ is universal, i.e., independent of details of the conversion rate distribution. It has indeed been checked that the weak residual dependence of data points on the latter distribution becomes even smaller as $M$ is further increased. The full green curve shows the heuristic prediction (DISPLAY_FORM87), providing a semi-quantitative picture of the finite-size scaling function. For instance, consensus is reached with probability ${\cal P}=1/2$ and ${\cal P}=2/3$ respectively for $x\approx 0.18$ and $x\approx 0.33$, according to actual data, whereas (DISPLAY_FORM87) respectively predicts $x=1/\sqrt{12}=0.288675\dots $ and $x=1/2$.Besides the value of the consensus probability ${\cal P}$, the next question is what determines whether or not the system reaches consensus. In Section SECREF67 it has been demonstrated that the average conversion rate $\overline{C}$ defined in (DISPLAY_FORM70) essentially determines the fate of the system in the regime where migration effects dominate conversion effects. It has also been shown that the consensus denoted by I1, where language 1 survives, can only be stable for $\overline{C}>0$, whereas the consensus denoted by I2, where language 2 survives, can only be stable for $\overline{C}<0$. The above statements are made quantitative in Figure FIGREF90, showing the probability distribution of the average conversion rate $\overline{C}$, for a Gaussian distribution of conversion rates with $w=1$. The total (i.e., unconditioned) distribution (black curves) is Gaussian. Red and blue curves show the distributions conditioned on consensus. They are indeed observed to live entirely on $\overline{C}>0$ for I1 and on $\overline{C}<0$ for I2. Finally, the distributions conditioned on coexistence (green curves, denoted by II) exhibit narrow symmetric shapes around the origin. Values of the migration rate $\gamma $ are chosen so as to have three partial histograms with equal weights, i.e., a consensus probability ${\cal P}=2/3$. This fixes $\gamma \approx 0.351$ for $M=2$ (top) and $\gamma \approx 10.22$ for $M=10$ (bottom).An area of interest that is common to both physicists and linguists concerns the evolution of competing languages. It was long assumed that such competition would result in the dominance of one language above all its competitors, until some recent work hinted that coexistence might be possible under specific circumstances. We argue here that coexistence of two or more competing languages can result from two symmetry-breaking mechanisms – due respectively to imbalanced internal dynamics and spatial heterogeneity – and engage in a quantitative exploration of the circumstances which lead to this coexistence. In this work, both symmetry-breaking scenarios are dealt with on an equal footing.In the first case of competing languages in a single geographical area, our introduction of an interpolation parameter $q$, which measures the amount of imbalance in the internal dynamics, turns out to be crucial for the investigation of language coexistence. It is conceptually somewhat subtle, since it appears only in the saturation terms in the coupled logistic equations used here to describe language competition; in contrast to the conversion terms (describing language shift from a less to a more favoured language), its appearance is symmetric with respect to both languages. For multiply many competing languages, the ensuing rate equations for the fractions of speakers are seen to bear a strong resemblance to a broad range of models used in theoretical ecology, including Lotka-Volterra or predator-prey systems.We first consider the case where the $N$ languages in competition in a single area have equally spaced attractivenesses. This simple situation allows for an exact characterisation of the stationary state. The range of attractivenesses is measured by the mean attractiveness $g$. As this parameter is increased, the number $K$ of surviving languages decreases progressively, as the least favoured languages successively become extinct at threshold values of $g$. Importantly, the range of values of $g$ between the start of the disappearances and the appearance of consensus grows proportionally to $N^2$. There is therefore a substantial amount of coexistence between languages that are significantly attractive.In the general situation, where the attractivenesses of the competing languages are modelled as random variables with an arbitrary distribution, the outcomes of numerical studies at finite $N$ are corroborated by a detailed asymptotic analysis in the regime of large $N$. One of the key results is that the quantity $W=Nw$ (the product of the number of languages $N$ with the mean attractiveness $w$) determines many quantities of interest, including the mean fraction $R=\left\langle K\right\rangle /N$ of surviving languages. The relation between $W$ and $R$ is however non-universal, as it depends on the full attractiveness distribution. This non-universality is most prominent in the regime where the mean attractiveness is large, so that only the few most favoured languages survive in the stationary state. The number of such survivors is found to obey a scaling law, whose non-universal critical exponent is dictated by the specific form of the attractiveness distribution near its upper edge.As far as symmetry breaking via spatial heterogeneity is concerned, we consider the paradigmatic case of two competing languages in a linear array of $M$ geographic areas, whose neighbours are linked via a uniform migration rate $\gamma $. In the simplest situation of two areas, we determine the full phase diagram of the model as a function of $\gamma $ as well as the conversion rates ruling language shift in each area. This allows us to associate different regions of phase space with either consensus or coexistence. Our analysis is then generalised to longer arrays of $M$ linked geographical regions. We first consider ordered attractiveness profiles, where language 1 is favoured in the $K$ leftmost areas, while language 2 is favoured in the $L$ rightmost ones. If the two blocks are of equal size so that no language is globally preferred, coexistence always results; however, the spatial profiles of the language speakers themselves are rather non-trivial. For blocks of unequal size, there is a transition from a situation of coexistence at low migration rates to a situation of uniform consensus at high migration rates, where the language favoured in the larger block is the only survivor in all areas. The critical migration rate at this transition grows as $M^2$. We next investigate disordered attractiveness profiles, where conversion rates are modelled as random variables. There, the probability of observing a uniform consensus is given by a universal scaling function of $x=\gamma /(M^{3/2}w)$, where $w$ is the width of the symmetric distribution of conversion rates.The ratio between migration and conversion rates beyond which there is consensus – either with certainty or with a sizeable probability – grows with the number of geographic areas as $M^2$ for ordered profiles of attractivenesses, and as $M^{3/2}$ for disordered ones. The first exponent is a consequence of the diffusive nature of migrations, whereas the second one has been derived in Appendix SECREF134 and related to anomalous band-edge scaling in one-dimensional Anderson localisation. If geographical areas were arranged according to a more complex geometric structure, these exponents would respectively read $2d/d_s$ and $(4-d_s)/(2d_s)$, with $d$ and $d_s$ being the fractal and spectral dimensions of the underlying structure (see BIBREF38, BIBREF39, and BIBREF40, BIBREF41 for reviews).Finally, we remark on another striking formal analogy – that between the rate equations (DISPLAY_FORM20) presented here, and those of a spatially extended model of competitive dynamics BIBREF42, itself inspired by a model of interacting black holes BIBREF43. In the latter, the non-trivial patterns of survivors on various networks and other geometrical structures were a particular focus of investigation, and led to the unearthing of universal behaviour. We believe that a network model of competing languages which combines both the symmetry-breaking scenarios discussed in this paper, so that every node corresponds to a geographical area with its own imbalanced internal dynamics, might lead to the discovery of similar universalities.AM warmly thanks the Leverhulme Trust for the Visiting Professorship that funded this research, as well as the Faculty of Linguistics, Philology and Phonetics at the University of Oxford, for their hospitality.Both authors contributed equally to the present work, were equally involved in the preparation of the manuscript, and have read and approved the final manuscript.This Appendix is devoted to an analytical investigation of the statistics of surviving languages in a single geographic area, in the regime where the numbers $N$ of competing languages is large.The properties of the attractiveness distribution of the languages are key to determining whether coexistence or consensus will prevail. In particular the transition to consensus depends critically, and non-universally, on the way in which the attractiveness distribution decays, as will be shown below.Statistical fluctuations between various instances of the model become negligible for large $N$, so that sharp (i.e., self-averaging) expressions can be obtained for many quantities of interest.Let us begin with the simplest situation where all languages survive. When the number $N$ of competing languages is large, the condition for this to occur assumes a simple form. Consider the expression (DISPLAY_FORM36) for $x_N$. The law of large numbers ensures that the sum $S$ converges towhereas $a_N$ is relatively negligible. The condition that all the $N$ competing languages survive therefore takes the form of a sharp inequality at large $N$, i.e.,All over this regime, the expression for $x_N$ simplifies toThe above analysis can be extended to the general situation where the numbers $N$ of competing languages and $K$ of surviving ones are large and comparable, with the fraction of surviving languages,taking any value in the range $0<R<1$.The rescaled attractiveness of the least favoured surviving language, namelyturns out to play a key role in the subsequent analysis. Let us introduce for further reference the truncated moments ($k=0,1,2$)First of all, the relationship between $R$ and $\eta $ becomes sharp in the large-$N$ regime. We have indeedThe limits of all quantities of interest can be similarly expressed in terms of $\eta $. We have for instancefor the sum introduced in (DISPLAY_FORM37). The marginal stability condition, namely that language number $K$ is at the verge of becoming extinct, translates toThe asymptotic dependence of the fraction $R$ of surviving languages on the rescaled mean attractiveness $W$ is therefore given in parametric form by (DISPLAY_FORM97) and (DISPLAY_FORM99). The identitydemonstrates that $R$ is a decreasing function of $W$, as it should be.When the parameter $W$ reaches unity from above, the model exhibits a continuous transition from the situation where all languages survive. The parameter $\eta $ vanishes linearly aswith unit prefactor, irrespective of the attractiveness distribution. The fraction of surviving languages departs linearly from unity, according toIn the regime where $W\gg 1$, the fraction $R$ of surviving languages is expected to fall off to zero. As a consequence of (DISPLAY_FORM97), $R\ll 1$ corresponds to the parameter $\eta $ being close to the upper edge of the attractiveness distribution $f(\xi )$. This is to be expected, as the last surviving languages are the most attractive ones. As a consequence, the form of the relationship between $W$ and $R$ for $W\gg 1$ is highly non-universal, as it depends on the behavior of the distribution $f(\xi )$ near its upper edge. It turns out that the following two main classes of attractiveness distributions have to be considered.Class 1: Power law at finite distance.Consider the situation where the distribution $f(\xi )$ has a finite upper edge $\xi _0$, and either vanishes or diverges as a power law near this edge, i.e.,The exponent $\alpha $ is positive. The density $f(\xi )$ diverges near its upper edge $\xi _0$ for $0<\alpha <1$, whereas it vanishes near $\xi _0$ for $\alpha >1$, and takes a constant value $f(\xi _0)=A$ for $\alpha =1$.In the relevant regime where $\eta $ is close to $\xi _0$, the expressions (DISPLAY_FORM97) and (DISPLAY_FORM99) simplify toEliminating $\eta $ between both above estimates, we obtain the following power-law relationship between $W$ and $R$:In terms of the original quantities $K$ and $w$, the above result readsSetting $K=1$ in this estimate, we predict that the consensus probability ${\cal P}$ becomes appreciable whenClass 2: Power law at infinity.Consider now the situation where the distribution extends up to infinity, and falls off as a power law, i.e.,The exponent $\beta $ is larger than 2, in order for the first two moments of $\xi $ to be convergent.In the relevant regime where $\eta $ is large, the expressions (DISPLAY_FORM97) and (DISPLAY_FORM99) simplify toEliminating $\eta $ between both above estimates, we obtain the following power-law relationship between $W$ and $R$:In terms of the original quantities $K$ and $w$, the above result readsSetting $K=1$ in this estimate, we predict that the consensus probability ${\cal P}$ becomes appreciable whenWe now summarise the above discussion. In the regime where $W\gg 1$, the fraction $R$ of surviving languages falls off as a power law of the formwhere the positive exponent $\lambda $ varies continuously, according to whether the distribution of attractivenesses extends up to a finite distance or infinity (see (DISPLAY_FORM106), (DISPLAY_FORM112)):In the marginal situation between both classes mentioned above, comprising e.g. the exponential distribution, the decay exponent sticks to its borderline valueThe decay law $R\sim 1/W$ might however be affected by logarithmic corrections.Another view of the above scaling laws goes as follows. When the number of languages $N$ is large, the number of surviving languages decreases from $K=N$ to $K=1$ over a very broad range of mean attractivenesses. The condition for all languages to survive (see (DISPLAY_FORM92)) sets the beginning of this range asThe occurrence of a sizeable consensus probability ${\cal P}$ sets the end of this range aswhere the exponent $\mu >-1/2$ varies continuously, according to (see (DISPLAY_FORM108), (DISPLAY_FORM114)):In the marginal situation between both classes, the above exponent sticks to its borderline valueThe extension of the dynamical range, defined as the ratio between both scales defined above, diverges asWe predict in particular a linear divergence for the exponential distribution ($\mu =0$) and a quadratic divergence for the uniform distribution ($\mu =1$). This explains the qualitative difference observed in Figure FIGREF50. The slowest growth of the dynamical range is the square-root law observed for distributions falling off as a power-law with $\beta \rightarrow 2$, so that $\mu =-1/2$.To close, let us underline that most of the quantities met above assume simple forms for the uniform and exponential distributions (see (DISPLAY_FORM44)).Uniform distribution.The consensus probability (see (DISPLAY_FORM42)) readsFor large $N$, this becomes ${\cal P}\approx \exp (-N/(2w))$, namely a function of the ratio $w/N$, in agreement with (DISPLAY_FORM119) and (DISPLAY_FORM120), with exponent $\mu =1$, since $\alpha =1$.The truncated moments readWe thus obtainwith exponent $\lambda =1/2$, in agreement with (DISPLAY_FORM106) and (DISPLAY_FORM116) for $\alpha =1$.Exponential distribution.The consensus probability readsirrespective of $N$, in agreement with (DISPLAY_FORM119), with exponent $\mu =0$ (see (DISPLAY_FORM121)).The truncated moments readWe thus obtainwith exponent $\lambda =1$, in agreement with (DISPLAY_FORM117).This Appendix is devoted to stability matrices and their spectra. Let us begin by reviewing some general background (see e.g. BIBREF44 for a comprehensive overview). Consider an autonomous dynamical system defined by a vector field ${E}({x})$ in $N$ dimensions, i.e., by $N$ coupled first-order equations of the formwith $m,n=1,\dots ,N$, where the right-hand sides depend on the dynamical variables $\lbrace x_n(t)\rbrace $ themselves, but not explicitly on time.Assume the above dynamical system has a fixed point $\lbrace x_m\rbrace $, such that $E_m\lbrace x_n\rbrace =0$ for all $m$. Small deviations $\lbrace \delta x_m(t)\rbrace $ around the fixed point $\lbrace x_m\rbrace $ obey the linearised dynamics given by the stability matrix ${S}$, i.e., the $N\times N$ matrix defined bywhere right-hand sides are evaluated at the fixed point. The fixed point is stable, in the strong sense that small deviations fall off exponentially fast to zero, if all eigenvalues $\lambda _a$ of ${S}$ have negative real parts. In this case, if all the $\lambda _a$ are real, their opposites $\omega _a=-\lambda _a>0$ are the inverse relaxation times of the linearised dynamics. In particular, the opposite of the smallest eigenvalue, simply denoted by $\omega $, characterises exponential convergence to the fixed point for a generic initial state. If some of the $\lambda _a$ have non-zero imaginary parts, convergence is oscillatory.The analysis of fixed points and bifurcations in low-dimensional Lotka-Volterra and replicator equations has been the subject of extensive investigations BIBREF23, BIBREF24, BIBREF25, BIBREF26, BIBREF27, BIBREF28 (see also BIBREF19, BIBREF20, BIBREF21).The remainder of this Appendix is devoted to the stability matrices involved in the array models considered in Section SECREF3, for an arbitrarily large number $M$ of geographical areas. All those stability matrices are related to the symmetric $M\times M$ matrixrepresenting (minus) the Laplacian operator on a linear array of $M$ sites, with Neumann boundary conditions. References BIBREF45, BIBREF46 provide reviews on the Laplacian and related operators on graphs.The eigenvalues $\lambda _a$ of ${\Delta }_M$ and the corresponding normalised eigenvectors ${\phi }_a$, such that ${\Delta }_M{\phi }_a=\lambda _a{\phi }_a$ and ${\phi }_a\cdot {\phi }_b=\delta _{ab}$, read($a=0,\dots ,M-1$). The vanishing eigenvalue $\lambda _0=0$ corresponds to the uniform eigenvector $\phi _{0,m}=1/\sqrt{M}$.Let us begin by briefly considering the simple example of the stability matrixof the rate equations (DISPLAY_FORM54) for the total populations $P_m(t)$. Its eigenvalues are $-1-\gamma \lambda _a$. The smallest of them is $-1$, so that the inverse relaxation time is given by $\omega =1$, as announced below (DISPLAY_FORM54).Let us now consider the stability matricesrespectively defined in (DISPLAY_FORM68) and (DISPLAY_FORM69), and corresponding to both uniform consensus states for an arbitrary profile of conversion rates $C_m$. The ensuing stability conditions have been written down explicitly in (DISPLAY_FORM60) and (DISPLAY_FORM62) for $M=2$. It will soon become clear that it is virtually impossible to write them down for an arbitrary size $M$. Some information can however be gained from the calculation of the determinants of the above matrices. They only differ by a global sign change of all the conversion rates $C_m$, so that it is sufficient to consider ${S}_M^{(1)}$. It is a simple matter to realise that its determinant readswhere $u_m$ is a generalised eigenvector solving the following Cauchy problem:with initial conditions $u_0=u_1=1$. We thus obtain recursivelyand so on. The expression (DISPLAY_FORM141) for $D_2$ agrees with the second of the conditions (DISPLAY_FORM60) and with the equation of the red curve in Figure FIGREF66, as should be. The expression () for $D_3$ demonstrates that the complexity of the stability conditions grows rapidly with the system size $M$.In the case of random arrays, considered in Section SECREF84, the conversion rates $C_m$ are independent random variables such that $\left\langle C_m\right\rangle =0$ and $\left\langle C_m^2\right\rangle =w^2$.The regime of most interest is where the conversion rates $C_n$ are small with respect to $\gamma $. In this regime, the determinant $D_M$ can be expanded as a power series in the conversion rates. The $u_m$ solving the Cauchy problem (DISPLAY_FORM140) are close to unity. Settingwhere the $u_m^{(1)}$ are linear and the $u_m^{(2)}$ quadratic in the $C_n$, we obtain after some algebrawhereare respectively linear and quadratic in the $C_n$. We haveIn Section SECREF84 we need an estimate of the probability ${\cal Q}$ that $\overline{C}=X/M$ is atypically small. Within the present setting, it is natural to define the latter event as $\left|X\right|<\left|Y\right|$. The corresponding probability can be worked out proviso we make the ad hoc simplifying assumptions – that definitely do not hold in the real world – that $X$ and $Y$ are Gaussian and independent. Within this framework, the complex Gaussian random variablehas an isotropic density in the complex plane. We thus obtainThe aim of this last section is to investigate the spectrum of the stability matrix ${S}_M^{(1)}$ associated with the ordered profile of conversion rates given by (DISPLAY_FORM72).In this case, the generalised eigenvector $u_m$ solving the Cauchy problem (DISPLAY_FORM140) can be worked out explicitly. We have $C_m=1$ for $m=1,\dots ,K$, and therefore $u_m=a{\rm e}^{m\mu }+b{\rm e}^{-m\mu }$, where $\mu >0$ obeys the dispersion relationThe initial conditions $u_0=u_1=1$ fix $a$ and $b$, and soSimilarly, we have $C_m=-1$ for $m=K+\ell $, with $\ell =1,\dots ,L$, and therefore $u_m=\alpha {\rm e}^{{\rm i}\ell q}+\beta {\rm e}^{-{\rm i}\ell q}$, where $0<q<\pi $ obeys the dispersion relationMatching both solutions for $m=K$ and $K+1$ fixes $\alpha $ and $\beta $, and soInserting the latter result into (DISPLAY_FORM139), we obtain the following expression for the determinant of ${S}_M^{(1)}$, with $M=K+L$:The vanishing of the above expression, i.e.,signals that one eigenvalue of the stability matrix ${S}^{(1)}$ vanishes. In particular, the consensus state where language 1 survives becomes marginally stable at the threshold migration rate $\gamma _c$, where the largest eigenvalue of ${S}^{(1)}$ vanishes. Equation (DISPLAY_FORM153) amounts to a polynomial equation of the form $P_{K,L}(\gamma )=0$, where the polynomial $P_{K,L}$ has degree $K+L-1=M-1$. All its zeros are real, and $\gamma _c$ is the largest of them. The first of these polynomials read

Q: What languages do they look at?