Two-community noisy Kuramoto model

TL;DR

This paper investigates a two-community noisy Kuramoto model with mixed positive and negative interactions, revealing non-symmetric synchronized solutions, bifurcation phenomena, and a detailed phase diagram through analytical and numerical methods.

Contribution

It introduces a detailed analysis of bifurcations and phase transitions in a two-community noisy Kuramoto model with mixed interactions, extending understanding beyond symmetric solutions.

Findings

Existence of non-symmetric synchronized solutions bifurcating from symmetric ones.

Derivation of bifurcation conditions and phase diagram classification.

Numerical simulations confirming stability and transitions between solutions.

Abstract

We study the noisy Kuramoto model for two interacting communities of oscillators, where we allow the interaction in and between communities to be positive or negative (but not zero). We find that, in the thermodynamic limit where the size of the two communities tends to infinity, this model exhibits non-symmetric synchronized solutions that bifurcate from the symmetric synchronized solution corresponding to the one-community noisy Kuramoto model, even in the case where the phase difference between the communities is zero and the interaction strengths are symmetric. The solutions are given by fixed points of a dynamical system. We find a critical condition for existence of a bifurcation line, as well as a pair of equations determining the bifurcation line as a function of the interaction strengths. Using the latter we are able to classify the types of solutions that are possible and thereby identify the phase diagram of the system. We also analyze properties of the bifurcation line in the phase diagram and its derivatives, calculate the asymptotics, and analyze the synchronization level on the bifurcation line. Part of the proofs are numerically assisted. Lastly, we present some simulations illustrating the stability of the various solutions as well as the possible transitions between these solutions.

Figures

Equations

$$\mathrm{d}{\theta }_{1,i}(t)=$$

$$\mathrm{d}{\theta }_{2,j}(t)=$$

$$H_N({\theta }_1,\dots ,{\theta }_N)=-\frac{1}{N}\sum _{i=1}^N\sum _{j=1}^N{A}_{i,j}\cos ({\theta }_j(t)-{\theta }_i(t))+\sum _{i=1}^N{\theta }_i(t){\omega }_i$$

$$A:=(A_{i,j}{)}_{i,j=1,\dots ,N}=\left[\begin{array}{cccccccc}0& K_1& \dots & K_1& L_1& L_1& \dots & L_1\\ K_1& 0& \dots & K_1& L_1& L_1& \dots & L_1\\ ⋮& ⋮& \ddots & ⋮& L_1& L_1& \dots & L_1\\ K_1& K_1& \dots & 0& L_1& L_1& \dots & L_1\\ L_2& L_2& \dots & L_2& 0& K_2& \dots & K_2\\ L_2& L_2& \dots & L_2& K_2& 0& \dots & K_2\\ L_2& L_2& \dots & L_2& ⋮& ⋮& \ddots & ⋮\\ L_2& L_2& \dots & L_2& K_2& K_2& \dots & 0\end{array}\right]=\left[\begin{array}{cc}{K}_1{\text{1}}_{\ast }& L_1\text{1}\\ L_2\text{1}& K_2{\text{1}}_{\ast }\end{array}\right],$$

$$\mathrm{d}{\theta }_i(t)={\partial }_{{\theta }_i}{H}_N({\theta }_1,\dots ,{\theta }_N)\mathrm{d}t+D\mathrm{d}{W}_i(t),i=1,\dots ,N,$$

$$r_{1,N_1}(t){\mathrm{e}}^{\mathrm{i}{\psi }_{1,N_1}(t)}=$$

$$r_{2,N_2}(t){\mathrm{e}}^{\mathrm{i}{\psi }_{2,N_2}(t)}=$$

$$\mathrm{d}{\theta }_{1,i}(t)=$$

$$\mathrm{d}{\theta }_{2,j}(t)=$$

$${\nu }_{N_1,t}(\mathrm{d}\theta ,\mathrm{d}\omega )$$

$${\nu }_{N_2,t}(\mathrm{d}\theta ,\mathrm{d}\omega )$$

$$\displaystyle\frac{\partial p_{1}(t;\theta,\omega)}{\partial t}=\frac{D}{2}\frac{\partial^{2}p_{1}(t;\theta,\omega)}{\partial\theta^{2}}-\frac{\partial}{\partial\theta}\big{[}v_{1}(t;\theta,\omega)p_{1}(t;\theta,\omega)\big{]}$$

$$v_1(t;\theta ,\omega )=\omega +{\alpha }_1{K}_1{r}_1(t)\sin ({\psi }_1(t)-\theta )+{\alpha }_2{L}_1{r}_2(t)\sin ({\psi }_2(t)-\theta ),$$

$$\displaystyle\frac{\partial p_{2}(t;\theta,\omega)}{\partial t}=\frac{D}{2}\frac{\partial^{2}p_{2}(t;\theta,\omega)}{\partial\theta^{2}}-\frac{\partial}{\partial\theta}\big{[}v_{2}(t;\theta,\omega)p_{2}(t;\theta,\omega)\big{]}$$

$$v_2(t;\theta ,\omega )=\omega +{\alpha }_2{K}_2{r}_2(t)\sin ({\psi }_2(t)-\theta )+{\alpha }_1{L}_2{r}_1(t)\sin ({\psi }_1(t)-\theta ).$$

$$r_1(t){\mathrm{e}}^{\mathrm{i}{\psi }_1(t)}:=$$

$$r_2(t){\mathrm{e}}^{\mathrm{i}{\psi }_2(t)}:=$$

$$\displaystyle 0=\frac{D}{2}\frac{\partial^{2}p_{1}(\theta,\omega)}{\partial\theta^{2}}-\frac{\partial}{\partial\theta}\big{[}v_{1}(\theta,\omega)p_{1}(\theta,\omega)\big{]},$$

$$p_1(\theta ,\omega )=$$

$$\displaystyle A_{1}(\theta,\omega)=B_{1}(\theta,\omega)\Big{(}\mathrm{e}^{\frac{4\pi\omega}{D}}\int_{\mathbb{S}}\frac{\mathrm{d}\phi}{B_{1}(\phi,\omega)}+(1-\mathrm{e}^{\frac{4\pi\omega}{D}})\int_{0}^{\theta}\frac{\mathrm{d}\phi}{B_{1}(\phi,\omega)}\Big{)}$$

$$\displaystyle B_{1}(\theta,\omega)=\exp\Big{[}\frac{2\omega\theta}{D}+\frac{2\alpha_{2}L_{1}r_{2}\cos(\psi_{2}-\theta)}{D}+\frac{2\alpha_{1}K_{1}r_{1}\cos(\psi_{1}-\theta)}{D}\Big{]}.$$

$$\displaystyle 0=\frac{D}{2}\frac{\partial^{2}p_{2}(\theta,\omega)}{\partial\theta^{2}}-\frac{\partial}{\partial\theta}\big{[}v_{2}(\theta,\omega)p_{2}(\theta,\omega)\big{]},$$

$$p_2(\theta ,\omega )=$$

$$\displaystyle A_{2}(\theta,\omega)=B_{2}(\theta,\omega)\Big{(}\mathrm{e}^{\frac{4\pi\omega}{D}}\int_{\mathbb{S}}\frac{\mathrm{d}\phi}{B_{2}(\phi,\omega)}+(1-\mathrm{e}^{\frac{4\pi\omega}{D}})\int_{0}^{\theta}\frac{\mathrm{d}\phi}{B_{2}(\phi,\omega)}\Big{)}$$

$$\displaystyle B_{2}(\theta,\omega)=\exp\Big{[}\frac{2\omega\theta}{D}+\frac{2\alpha_{1}L_{2}r_{1}\cos(\psi_{1}-\theta)}{D}+\frac{2\alpha_{2}K_{2}r_{2}\cos(\psi_{2}-\theta)}{D}\Big{]}.$$

$$r_1$$

$$r_2$$

$$0$$

$$\mathrm{d}{\psi }_m(t)=\sum _{j=1}^{N_m}\frac{\partial {\psi }_m}{\partial {\theta }_{m,j}}\mathrm{d}{\theta }_{m,j}+\frac{1}{2}\sum _{j=1}^{N_m}\frac{{\partial }^2{\psi }_m}{\partial {\theta }_{m,j}^2}(\mathrm{d}{\theta }_{m,j}{)}^2,m\in \{1,2\}.$$

$$\frac{\partial {\psi }_m}{\partial {\theta }_{m,j}}=\frac{1}{N_m{r}_m(t)}\cos ({\psi }_m(t)-{\theta }_{m,j}(t)),m\in \{1,2\},$$

$$\frac{{\partial }^2{\psi }_m}{\partial {\theta }_{m,j}^2}=$$

$$-\frac{2}{(N_m{r}_m(t){)}^2}\sin ({\psi }_m(t)-{\theta }_{m,j}(t))\cos ({\psi }_m(t)-{\theta }_{m,j}(t)),m\in \{1,2\}.$$

$$\mathrm{d}{\psi }_1(t)=$$

Full text

Key words and phrases. Two-community network, phase oscillators, noisy Kuramoto model, McKean-Vlasov equations, phase diagram, bifurcation line.

Acknowledgment. This work is supported by NWO Gravitation Grant 024.002.003–NETWORKS. The author is grateful to F. den Hollander for discussions and detailed comments, and to C. da Costa and E. Verbitskiy for fruitful input.

1 Background and motivation

The motivation for studying the two-community noisy Kuramoto model is two-fold. On the one hand, the suprachiasmatic nucleus (SCN) in the brain of mammals is responsible for biological time-keeping and consists of two communities of cells that exhibit synchronization [28]. On the other hand, there are recent studies of interacting particle systems with community structure, that reveal vast richness in behavior [12, 3, 9, 5]. The noisy Kuramoto model consists of a collection of oscillators with a mean-field interaction that favors alignment subject to external noisy [23].

The SCN is a cluster of neurons responsible for dictating the rhythm of bodily functions, most significantly the sleep-cycle. Malfunctioning of the SCN leads to a variety of health problems, ranging from epilepsy to narcolepsy. Remarkably, the network structure of the cluster is similar in all mammals, with the universal feature that it is split into two communities. In humans each cluster has a size of about $10^{4}$ neurons. It seems that this two-community structure is ideal, both for the robustness of the rhythm of the cluster not to be disturbed by unusual light inputs, as well as for the cluster to be adaptable enough to re-synchronize when there is a change in the light-dark cycle it is exposed to. As we will see below, this is reflected by the mathematical properties of the two-community noisy Kuramoto model, for which the interplay between positive and negative interactions introduces new features. The negative interaction, studied before in [15], [16], seems to play a key role in the appearance of a negative correlation between the neurons in the two communities in the SCN, resulting in new emergent behavior such as phase splitting [17].

In the mathematics literature there have been recent studies on bipartite mean-field spin systems [12], as well as on the Ising block model [3] and the asymmetric Curie-Weiss model [9], [5], where the splitting into two communities introduces interesting features, for example, the appearance of periodic orbits. These are discrete models which makes them hard to analyze. What makes the Kuramoto model considered here hard to analyze is that the interaction between phase oscillators in the Kuramoto model is non-linear.

Also in [25] the authors consider the two-community noisy Kuramoto model. They find an intricate phase diagram, with the system being able to take on a variety of different states. This confirms the observation that a simple modification in the network structure can greatly increase the complexity of the system. The results in [25], however, depend strongly on a Gaussian approximation for the phase distribution in each community (explained in [26]), which allows for a reduction of the dynamics to a low-dimensional setting. In this paper we do not rely on any such approximation.

We have recently studied the noisy Kuramoto model on the hierarchical lattice [13], finding conditions for synchronization either to propagate to all levels in the hierarchy or to vanish at a finite level. This analysis came about by writing down renormalized evolution equations for the average phases in a block-community at a given hierarchical level in the hierarchical mean-field limit. In the present paper we allow for negative interactions across the communities, a situation we did not consider in the hierarchical model.

In Section 2 we introduce the noisy Kuramoto model on the two-community network (see Fig. 1) and show that the empirical measures defined for each community evolve according to a McKean-Vlasov equation in the thermodynamic limit. We also give the steady-state solutions to these McKean-Vlasov equations and conjecture which values the phase difference between the average phases of the two communities can take in the steady state. In Section 3 we present results on the critical condition for synchronization in the case of symmetric interaction strengths and equal community sizes, first without disorder and then with disorder. By disorder we mean that the natural frequencies of the oscillators are taken from a distribution while without disorder means that all oscillators are assumed to have a natural frequency of zero. In Section 4 we prove the conjecture from the previous section for a simplified version of the model where we take the interaction strengths to be symmetric and prove the existence of non-symmetric solutions in this case. Here symmetric solutions are solutions in which the synchronization level is the same in both communities while non-symmetric solutions are solutions where the synchronization level are non-zero and not the same in both communities. We also characterize the bifurcation line at which the non-symmetric solutions split off from the symmetric solutions, and expound a collection of results on the (asymptotic) properties of the bifurcation line in the phase diagram. Furthermore we analyze the synchronization level along the bifurcation line. Some of the proofs in Section 4 are numerically assisted. Finally, in Section 5 we present some simulations illustrating the stability of the various solutions as well as the possible transitions between various steady-states.

2 Basic properties

In Section 2.1 we define the model, in Section 2.2 we take the McKean-Vlasov limit, and in Section 2.3 we identify the stationary solutions.

2.1 Model

We consider two communities of oscillators of size $N_{1}$ and $N_{2}$ with internal mean-field interactions of strength $\frac{K_{1}}{N_{1}}$ and $\frac{K_{2}}{N_{2}}$, respectively. In addition, the oscillators in community 1 experience a mean-field interaction with the oscillators in community 2 of strength $\frac{L_{1}}{N_{2}}$ and the oscillators in community 2 experience a mean-field interaction of strength $\frac{L_{2}}{N_{1}}$ with the oscillators in community 1. Here we will take $K_{1},K_{2}\in\mathbb{R}$ to be positive and $L_{1},L_{2}\in\mathbb{R}\setminus\{0\}$.

Definition 2.1** (Two-community noisy Kuramoto model).**

The phase angles of the oscillators in community 1 are denoted by $\theta_{1,i}$, $i=1,\cdots,N_{1}$, and their evolution on $\mathbb{S}=\mathbb{R}/2\pi$ is governed by the SDE

[TABLE]

The phase angles of the oscillators in community 2 are denoted by $\theta_{2,j}$, $j=1,\cdots,N_{2}$, and their evolution on $\mathbb{S}=\mathbb{R}/2\pi$ is governed by the SDE

[TABLE]

Here, the natural frequencies $\omega_{1,i}$, $i=1,\ldots,N_{1}$, of the oscillators in community 1 are drawn independently from a probability distribution $\mu_{1}(\mathrm{d}\omega)$ on $\mathbb{R}$ and the natural frequencies $\omega_{2,i}$, $i=1,\ldots,N_{2}$, of the oscillators in community 2 are drawn independently from a probability distribution $\mu_{2}(\mathrm{d}\omega)$ on $\mathbb{R}$, while $D>0$ is the noise strength, and \big{(}W_{1,i}(t)\big{)}_{t\geq 0},i=1,\ldots,N_{1}, and \big{(}W_{2,j}(t)\big{)}_{t\geq 0},j=1,\ldots,N_{2}, are independent standard Brownian motions. For simplicity we take $\mu_{1},\mu_{2}$ to be symmetric and have the same mean which we can assume to be zero without loss of generality.

The model can alternatively be defined in terms of an interaction Hamiltonian and a weighted adjacency matrix, given by

[TABLE]

with

[TABLE]

where $\textbf{1}=$ all 1’s and $\textbf{1}_{*}=$ all 1’s, except for 0’s on the diagonal. The model then reads

[TABLE]

where $N=N_{1}+N_{2}$. Here, we identify phase angle $\theta_{i}$ with the oscillators in community 1 when $i\in[1,N_{1}]$ and with the oscillators in community 2 when $i\in(N_{1},N_{1}+N_{2}]$. This representation of the model illustrates the network structure of the underlying interactions and in principle the adjacency matrix can be replaced by a matrix arising from a random graph model and has recently been addressed by a number of authors [2, 7, 10, 20, 21]. This however significantly complicates the calculations since the interactions are no longer expressible in terms of a closed function of the empirical measure. The representation via the Hamiltonian may also provide a method for studying the stability properties of the stationary states.

The following order parameters allow us to monitor the dynamics in each community:

[TABLE]

where $r_{1,N_{1}}(t)\in[0,1]$ and $r_{2,N_{2}}(t)\in[0,1]$ represent the synchronization levels, and $\psi_{1,N_{1}}(t)$ and $\psi_{2,N_{2}}(t)$ represent the average phases, in community 1 and 2, respectively. Using these order parameters, we can rewrite the evolution equations in (2.1) and (2.2) as

[TABLE]

and

[TABLE]

2.2 McKean-Vlasov limit

We assume that the sizes of the communities are related to one another by setting $N_{1}=\alpha_{1}N$ and $N_{2}=\alpha_{2}N$, $\alpha_{1}+\alpha_{2}=1$. In the limit as $N\to\infty$, we expect the angle density of oscillators in each community to follow a McKean-Vlasov equation. Define the empirical measure for each community $(\theta\in\mathbb{S},\omega\in\mathbb{R})$:

[TABLE]

Proposition 2.2** (McKean-Vlasov limit).**

In the limit as $N\to\infty$, the empirical measure $\nu_{N_{1},t}(\mathrm{d}\theta,\mathrm{d}\omega)$ converges to $\nu_{1,t}(\mathrm{d}\theta,\mathrm{d}\omega)=p_{1}(t;\theta,\omega)\,\mathrm{d}\theta\,\mathrm{d}\omega$, and the empirical measure $\nu_{N_{2},t}(\mathrm{d}\theta,\mathrm{d}\omega)$ converges to $\nu_{2,t}(\mathrm{d}\theta,\mathrm{d}\omega)=p_{2}(t;\theta,\omega)\,\mathrm{d}\theta\,\mathrm{d}\omega$, where $p_{1}(t;,\theta,\omega)$ evolves according to

[TABLE]

with

[TABLE]

and $p_{2}(t;\theta,\omega)$ evolves according to

[TABLE]

with

[TABLE]

Here, $r_{1}(t),r_{2}(t),\psi_{1}(t)$ and $\psi_{2}(t)$ are defined by

[TABLE]

The convergence is in $\mathcal{C}([0,T],\mathcal{M}_{1}(\mathbb{S}\times\mathbb{R}))$ and takes place for any $T>0$. Here we consider annealed convergence with respect to the natural frequencies.

Proof.

The proof is analogous to that in the case of the one-community noisy Kuramoto model in [8] with straightforward modifications. ∎

2.3 Stationary solutions

The stationary solutions of the McKean-Vlasov limit in Proposition 2.2 give the possible states the system can assume in the long time limit. These are presented in the next proposition.

Proposition 2.3** (Stationary solutions).**

In the cases $r_{1}=r_{2}=0$ and $r_{1},r_{2}>0$, the stationary density $p_{1}(\theta,\omega)$ solves the equation

[TABLE]

which has solution

[TABLE]

where

[TABLE]

with

[TABLE]

The stationary density $p_{2}(\theta,\omega)$, solves the equation

[TABLE]

which has solution

[TABLE]

where

[TABLE]

with

[TABLE]

In addition, the following self-consistency equations must be satisfied:

[TABLE]

Proof.

Note that in the case when $r_{1}=r_{2}=0$, both stationary densities are uniform on $\mathbb{S}$, i.e., $p_{1}(\theta,\omega)=p_{2}(\theta,\omega)=\frac{1}{2\pi}$, which satisfies (2.18) and (2.22). The proof in the case when $r_{1},r_{2}>0$ is analogous to the calculation given in [14, Solution to Exercise X.33]. ∎

Remark 2.4**.**

In the simplified version of the model we will consider below, we are able to prove that solutions of the type $r_{1}=0$ and $r_{2}>0$ (or vice versa) are not possible, but it is difficult to prove this in the general case considered above.

In order to understand the steady-state phase difference between the communities, we proceed heuristically as follows. For the stationary solutions we assume that $\hskip 4.26773ptr_{1}(t),\hskip 4.26773ptr_{2}(t),\\ \psi_{1}(t),\psi_{2}(t)$ reach their steady-state values $r_{1},r_{2},\psi_{1},\psi_{2}$ as $t\rightarrow\infty$ and assume that the parameters of the system are such that $r_{1},r_{2}>0$. For the synchronization levels the possible steady-state values are computed by solving the self-consistency equations in (2.26). For the average phases we use standard Itô-calculus to compute their evolution

[TABLE]

From the definition of the order parameters we have

[TABLE]

and

[TABLE]

Substituting (2.28)–(2.29) and (2.8)–(2.9) into (2.27), setting $N_{m}=\alpha_{m}N$ and taking the large-$N$ limit, we get the equations

[TABLE]

[TABLE]

Due to the last two self-consistency equations in (2.26) the last line of (2.30) and (2.31) is zero. For the steady-state average phases in the case when $\mu_{1}=\mu_{2}=\delta_{0}$, we must therefore simultaneously solve the equations

[TABLE]

[TABLE]

Since the system is invariant under rotations, we can set one of the two angles to zero. If we set $\psi_{1}=0$, then we see that the equation for $\psi_{2}$ is satisfied by taking $\psi_{2}=0$ or $\psi_{2}=\pi$. The above calculation is not rigorous, but does suggest the following conjecture.

Conjecture 2.5** (Steady-state phase difference).**

In the system without disorder, the phase difference $\psi=\psi_{2}-\psi_{1}$ between the two communities in the two-community noisy Kuramoto model with $K_{1}=K_{2}=K$ and $L_{1}=L_{2}=L\neq 0$ in the steady state can only be $\psi=0$ or $\psi=\pi$.

The intuition for this conjecture is that the system will try to maximize the interaction strength between oscillators in order to achieve the highest synchronization in each community. This will be achieved at $\psi=0$ when $L>0$ and at $\psi=\pi$ when $L<0$. The other combinations ($\psi=0$ with $L<0$ and $\psi=\pi$ with $L>0$) should also be possible, but should not be stable. For an illustration of stability properties obtained via simulations, we refer the reader to Section 5.

3 Symmetric interaction with fixed phase difference

In this section we pick $L_{1}=L_{2}=L$, $K_{1}=K_{2}=K$, $\alpha_{1}=\alpha_{2}$, $D=1$. In Section 3.1 we consider the case where the natural frequency of the oscillators is zero, and in Section 3.2 the case where the natural frequency of the oscillators is drawn from a symmetric distribution $\mu$ on $\mathbb{R}$.

3.1 Without disorder

Here we take $\mu_{1}=\mu_{2}=\delta_{0}$. This simplifies (2.19) and (2.23) to

[TABLE]

The self-consistency equations for $r_{1}$ and $r_{2}$ in (2.26) can be written in the form

[TABLE]

where $W(x)=\frac{2V(x)}{x},x\in(0,\infty),$ with

[TABLE]

The definitions of $a_{1},a_{2},b_{1}$ and $b_{2}$ will be given below. The function $V(x)$ is the same function that appears in the self-consistency equation of the one-community noisy Kuramoto model [11, Equation 2.2]. To see why the self-consistency equations can be written as in (3.3), note that

[TABLE]

with $\text{I}_{m}(x):=\frac{1}{2\pi}\int_{\mathbb{S}}\mathrm{d}\theta(\cos\theta)^{m}\exp(x\cos\theta)$ the modified Bessel functions of the first kind, so that

[TABLE]

Here we have used the identity $I_{0}(x)=I_{1}(x)$ given in [1, 9.6.27]. Using (3.6) and the trigonometric identity $\cos(a-b)=\cos a\cos b+\sin a\sin b$, $a,b\in\mathbb{R}$, we can rewrite the self-consistency equations for $r_{1}$ and $r_{2}$ as

[TABLE]

where

[TABLE]

Note that

[TABLE]

where we recall $\psi=\psi_{2}-\psi_{1}$. The most suggestive form of the self-consistency equations is in terms of $K,L$ and the phase difference $\psi$:

[TABLE]

and is obtained by substituting the expressions for $a_{1},a_{2},b_{1}$ and $b_{2}$ into (3.3).

Proposition 3.1** (Properties of $V$).**

$V(0)=0$. 2. 2.

$V^{\prime}(0)=\frac{1}{2}$. 3. 3.

$x\mapsto V(x)$* is strictly increasing on $[0,\infty)$.* 4. 4.

$x\mapsto V(x)$* is strictly concave on $[0,\infty)$.* 5. 5.

$V(x)<\frac{x}{2}$* for $x\in(0,\infty)$.* 6. 6.

$\lim_{x\to\infty}V(x)=1$. 7. 7.

$V(-x)=-V(x)$* for all $x\in(0,\infty)$.*

Proof.

Properties 1, 2, 3 and 6 are easily verified. Property 4 is proven by applying Lemma 4 in [22] (see Appendix A for a comprehensive proof). Property 5 is a direct consequence of properties 1, 2 and 4. For Property 7, use $-\cos(\theta)=\cos(\pi-\theta)$ to write

[TABLE]

By performing the change of variable $\phi=\pi-\theta$, we get $V(-x)=-V(x)$. ∎

Proposition 3.2** (Properties of $W$).**

$\lim_{x\downarrow 0}W(x)=1$. 2. 2.

$x\mapsto W(x)$* is continuous and strictly decreasing on $[0,\infty)$.* 3. 3.

$\lim_{x\to\infty}W(x)=0$.

Proof.

Properties 1 and 3 are easily verified. For property 2, note that

[TABLE]

so we need to verify that $V^{\prime}(x)<\frac{V(x)}{x}$. This is true by properties 1 and 4 in Proposition 3.1. ∎

In the case without disorder Conjecture 2.5 can be proven.

Proposition 3.3**.**

Fix $\psi_{1}=0$ and assume that $\mu_{1}=\mu_{2}=\delta_{0}$. Then the order parameters of the system are either $r_{1},r_{2}=0$ or $r_{1},r_{2}>0$ and $\psi\in\{0,\pi\}$.

Proof.

Here the set of self-consistency equations (2.26) simplify to

[TABLE]

Since the system is invariant under rotations we can set one of the average phase angles to zero. So take $\psi_{1}=0$ such that $\psi=\psi_{2}$. To determine which phase differences are possible we are left to solve

[TABLE]

Let us first consider the case when $r_{1}=0$. In this case (3.14) becomes

[TABLE]

and (3.15) becomes

[TABLE]

which is exactly the self-consistency equation for the one-community noisy Kuramoto model without disorder, and can be divided into two cases: Either $K\leq 2$, in which case $r_{2}=0$, making $(r_{1},r_{2})=(0,0)$ the only stationary solution, or $K>2$, in which case there is a unique $r_{2}>0$ solving (3.21). By making the change of variable $\vartheta=\psi_{2}-\theta$ in (3.20) and using the trigonometric identity $\cos(\psi_{2}-\vartheta)=\cos\psi_{2}\cos\vartheta+\sin\psi_{2}\sin\vartheta$ in (3.20), we see that (3.20), in this case, is only solved by $\psi_{2}=\frac{\pi}{2}$ or $\psi_{2}=\frac{3\pi}{2}$. In order to satisfy the self-consistency equations, these angles must satisfy (3.18) and (3.19) with $r_{1}=0$:

[TABLE]

The second equation is satisfied for all $\psi_{2}$, but the first equation is incompatible with $\psi_{2}=\frac{\pi}{2}$ as well as $\psi_{2}=\frac{3\pi}{2}$. so that the solution $r_{1}=0$ and $r_{2}>0$ is not possible, leaving only the solution $(r_{1},r_{2})=(0,0)$. Note that in this case the average angles are not well defined.

Let us next consider the case when $r_{1}>0$ (so that we must also have $r_{2}>0$). The allowed angles have to satisfy (3.18) and (3.19) simultaneously. These are satisfied only when $\sin\psi=0$, so that $\psi\in\{0,\pi\}$. ∎

Theorem 3.4** (Critical line without disorder).**

Fix $\psi=\psi_{2}-\psi_{1}\in\{0,\pi\}$. Then the parameter space $\{(K,L):K,L\in\mathbb{R}^{2}\}$ splits into two regions:

1. a)

In the region $K+L\cos\psi\leq 2$, there is precisely one solution: the unsynchronized solution $(r_{1},r_{2})=(0,0)$. 2. b)

In the region $K+L\cos\psi>2$, there are at least two solutions: the unsynchronized solution $(r_{1},r_{2})=(0,0)$ and the symmetric synchronized solution $(r_{1},r_{2})=(r,r)$ for some $r\in(0,1)$.

Proof.

For part a), note that $(0,0)$ always solves the self-consistency equations in (3.11), due to property 1 of Proposition 3.2 and the fact that $a_{1},a_{2},b_{1}$, $b_{2}$ are zero when $(r_{1},r_{2})=(0,0)$. The calculation given in the proof of Proposition 3.3 when $r_{1}=0$ shows that a solution of the form $r_{1}=0$ and $r_{2}>0$ is not possible, and due to symmetry the same is true for solutions with $r_{2}=0$ and $r_{1}>0$. To have strictly positive $r_{1},r_{2}$, we use property 5 in Proposition 3.1 to get

[TABLE]

Adding these equations, we get

[TABLE]

which is the condition to have positive synchronized solutions and defines the critical line. Let us next consider the case $\psi=0$ and $r_{1},r_{2}>0$. Then the self-consistency equations in (3.11) reduce to

[TABLE]

If we consider symmetric solutions so that $r_{1}=r_{2}=r$, then these two equations are identical and correspond to the self-consistency equation for the one-community noisy Kuramoto model with the replacement $2K\rightarrow K+L$, which has a positive solution when $K+L>2$. The same can be done when $\psi=\pi$ and yields $K-L>2$ as critical condition. ∎

It is tempting to conclude that the two-community model is the same as the one-community model with the replacement $2K\to K+L\cos\psi$. This is, however, not the case as we will see in Section 4.

3.2 With disorder

In this section we identify the critical line when we include disorder. We simplify the system by taking the distributions from which the natural frequencies are drawn in the two communities to be the same, i.e., $\mu_{1}=\mu_{2}=\mu$. Then the self-consistency equations in (2.26) read

[TABLE]

In light of Conjecture 2.5 we will restrict the following theorem to the two cases $\psi=0$ and $\psi=\pi$. Define

[TABLE]

Conjecture 3.5** (Critical line with disorder).**

Fix $\psi=\psi_{2}-\psi_{1}\in\{0,\pi\}$. If the disorder in the two communities is drawn from a symmetric unimodal distribution $\mu$, then the parameter space $\{(K,L):K,L\in\mathbb{R}^{2}\}$ splits into two regions:

1. a)

In the region $K+L\cos\psi\leq\chi^{-1}$, there is precisely one solution: the unsynchronized solution $(r_{1},r_{2})=(0,0)$. 2. b)

In the region $K+L\cos\psi>\chi^{-1}$, there are at least two solutions: the unsynchronized solution $(r_{1},r_{2})=(0,0)$ and the symmetric synchronized solution $(r_{1},r_{2})=(r,r)$ for some $r\in(0,1)$.

Heuristic Proof.

Following the method used in [23] for the one-community model, we Taylor expand the self-consistency equations for $r_{1}$ and $r_{2}$ in the two variables $r_{1}$ and $r_{2}$. The equations in (3.2) read, to first order,

[TABLE]

We can verify that $V_{1}^{\mu}(0,0)=V_{2}^{\mu}(0,0)=0$, and calculate the derivatives at zero. This leads to

[TABLE]

Adding these equations, we get

[TABLE]

Since we are considering the case where $\mu$ is symmetric, the last term vanishes and we obtain the critical line in Theorem 3.5. This shows that below the critical line the self-consistency equations are a contraction, making $(r_{1},r_{2})=(0,0)$ a fixed point. In order to show that solutions of the form $r_{1}=0$ and $r_{2}>0$ are not possible, we would have to repeat the calculation used in the proof of Proposition 3.3 for the general case. This turns out to be non-trivial, but we expect that it is possible to prove this for symmetric, unimodal $\mu$ by proving that $p_{1}(\theta+\psi_{2},\omega)=p_{1}(-\theta+\psi_{2},-\omega)$ (in the case that $\psi_{1}=0$) and using this symmetry to show that the first and third equation in (2.26) cannot be simultaneously satisfied when $r_{1}=0$ and $r_{2}>0$. If $\mu$ is symmetric and unimodal, then it is conjectured that the analog of $V_{1}^{\mu}(r_{1},r_{2})$ and $V_{2}^{\mu}(r_{1},r_{2})$ in the one-community noisy Kuramoto model is concave [see Conjecture [3.12], Chapter 3 in [19]]. We assume that this conjecture also holds in this case for both $V_{1}^{\mu}(r_{1},r_{2})$ and $V_{2}^{\mu}(r_{1},r_{2})$, at least for symmetric solutions. In the case $\psi=0$ the symmetric solution $r_{1}=r_{2}>0$ reduces the system of self-consistency equations in (3.2) to a single equation that is analogous to the one-community noisy Kuramoto model self-consistency equation [19, Proposition 3.10, Chapter 3] with the replacement $K\rightarrow K+L$. In the case $\psi=\pi$ we can perform a change of variable in the integral of the second line (3.2), namely, $\phi=\psi_{2}-\theta$, to see that the equations again reduce to the equation for the one-community case with the replacement $K\rightarrow K-L$. Thus, we see that in both cases we can apply the conjecture in [19, Conjecture 3.12] to ensure that the line $K+L\cos\psi=\chi^{-1}$ is the critical condition for symmetrically synchronized solutions, which settles the conditions in a) and b).

For part b), we must still show that the symmetric solution is possible above the critical line. Due to the reduction of the system to the one-community noisy Kuramoto model, both for $\psi=0$ and for $\psi=\pi$, we see that the symmetric solution indeed exists above the critical line. ∎

4 Bifurcation of non-symmetric solutions

In this section we consider the system with the same parameter specifications and simplifications as in Section 3, but without disorder and with $\psi=0$. The analysis with $\psi=\pi$ carries over after the replacement $L\to-L$ in the self-consistency equations in (3.1) (the resulting modified phase diagram is shown in the right panel of Fig. 8). The proofs in this section rely on numerics.

The self-consistency equations can be visualized as a vector field, in which the solutions to the equations appear as fixed points, by plotting

[TABLE]

For a certain range of parameters non-symmetric solutions appear, as seen in Fig. 3. The non-symmetric solutions appear to be saddle-points, having a stable and an unstable manifold under the vector field representing the self-consistency equations. Note that this vector field does not represent the dynamics of the system, since the self-consistency equations contain only the stationary densities. By plotting the possible solutions as functions of $K$ while keeping $L$ fixed, we see that the non-symmetric solutions bifurcate from the symmetric solutions, as is seen in Fig. 4 for the case where $L=-2$. The symmetric solutions correspond to equal amounts of synchronization in the two communities. This is also the only solution possible between $K=4$ and $K=4.9953\ldots$. At $K=4.9953\ldots$, the non-symmetric solutions appear, corresponding to one community having a larger synchronization level than the other community. Due to the symmetry of the system, both communities can have a higher level of synchronization in the non-symmetric solution.

In Section 4.1 we prove a necessary and sufficient condition for the existence of non-symmetric solutions. In Section 4.2 we show that the non-symmetric solutions are ordered and are such that the symmetric solution is wedged in between the two non-symmetric solutions. In Section 4.3 we analyze the (asymptotic) properties of the bifurcation line as well as the synchronization level along the bifurcation line.

4.1 Existence and characterization of non-symmetric solutions

Theorem 4.1** (Characterization of the bifurcation line).**

The existence of non-symmetric solutions requires $L<0$, in which case the bifurcation point $K^{*}=K^{*}(L)$ is the unique solution to the equation

[TABLE]

and the synchronization level at the bifurcation point is given by

[TABLE]

Proof.

We assume that a non-zero symmetric solution exists, so that $r_{1}=r_{2}=r$ and $r=V((K+L)r)$, which is the case when $K+L>2$. Let $(K^{*},r^{*})$ be a bifurcation point for fixed $L$. We will show via a perturbation argument that this bifurcation point exists and is unique. At the bifurcation point the non-symmetric solutions split off from the symmetric solution since $V$ is continuous. This allows us to perform a perturbation around $r^{*}$, namely,

[TABLE]

where $\epsilon$ and $\delta$ are small, either positive or negative, and are related, as will be shown shortly. We Taylor expand around the point $(K+L)r^{*}$ and use $r^{*}=V((K+L)r^{*})$, to get

[TABLE]

where by $\sim$ (here and in the rest of the paper) we mean that the ratio tends to 1 asymptotically. Abbreviate $C^{*}=V^{\prime}((K+L)r^{*})$. Then the equations in (4.6) combine to give

[TABLE]

which implies

[TABLE]

and $\epsilon\sim\pm\delta$. Using the negative sign would require the following two equations to be satisfied:

[TABLE]

However, these equations cannot be satisfied simultaneously with $r^{*}>0$. Indeed, the first finds the intersection point of $V$ with the line of slope $\frac{1}{K+L}$ passing through zero. But due to properties 1 and 4 we know that $V$ has slope $\frac{1}{K+L}$ before this intersection point. Thus, the two equations that must be satisfied at the bifurcation point are

[TABLE]

For fixed $L$, these equations determine both the value $r^{*}=r^{*}(L)$ of the synchronization level at the bifurcation point and the internal coupling strength $K^{*}=K^{*}(L)$ at which the bifurcation occurs. The first equation finds the intersection point of $V$ and the line with slope $\frac{1}{K+L}$ passing through zero. The second equation requires the derivative of $V$ at this point to be $\frac{1}{K-L}$. Due to the concavity of $V$ (Property 4 of Proposition 3.1), this gives the relation

[TABLE]

which implies that $L<0$, as claimed. To visualize the procedure for determining the bifurcation point, we plot the appropriate lines in Fig. 5. It is clear that the slope of the thickly dashed line must be less than that of the solid line, which gives $L<0$.

We can find an expression for the derivative of $V$ in (3.4) by writing

[TABLE]

For the first term in the right-hand side we can use the identity from [4, Eq. (2.21)], so that in our case

[TABLE]

where we have used (4.11) for the second term. This reduces (4.12) to

[TABLE]

To find $r^{*}=r^{*}(K^{*},L)$, we must find $K^{*}=K^{*}(L)$ that solves (4.11). Substituting (4.16) into (4.11), we obtain (4.2).

We will first prove that, given $r$, there is a unique $K^{*}$. In order to do this, we solve the equation

[TABLE]

for $L$ to find

[TABLE]

where we have taken the negative since we are dealing with the case $L<0$. In order to have a real solution, we require

[TABLE]

The equation for the bifurcation point in (4.2) reads

[TABLE]

where

[TABLE]

Clearly, $K\mapsto f_{r}(K)$ is strictly decreasing on $(\frac{2}{1-r^{2}},\infty)$ for $r\in(0,1)$. Since $x\mapsto V(x)$ is strictly increasing, $K\mapsto V(f_{r}(K)r)$ is strictly decreasing on $(\frac{2}{1-r^{2}},\infty)$. However, in order to satisfy (4.20) with $r\in(0,1)$, by property 5 in Proposition 3.1, we must have

[TABLE]

i.e.,

[TABLE]

Moreover,

[TABLE]

For fixed $r\in(0,1/\sqrt{2})$, $V(f_{r}(K)r)$ decreases from $V(\frac{2r}{1-r^{2}})$ to $V(2r)$ as $K$ increases from $\frac{2}{1-r^{2}}$ to $\frac{2(1-r^{2})}{1-2r^{2}}$ while for $r\in[1/\sqrt{2},1)$, $V(f_{r}(K)r)$ decreases from $V(\frac{2r}{1-r^{2}})$ to $V(\frac{r}{1-r^{2}})$ as $K$ increases from $\frac{2}{1-r^{2}}$ to $\infty$. In order to prove uniqueness, we need to show that

[TABLE]

Uniqueness follows because is be a unique $K^{*}$ satisfying (4.16), due to $V$ decreasing continuously from the upper to the lower bounds in (4.26) and (4.27) and the line $r$ being wedged between the bounds (note that $r\mapsto V(f_{r}(K)r)$ intersects $r$ exactly once). The curves $V(\frac{2r}{1-r^{2}})$, $r$, $V(\frac{r}{1-r^{2}})$ are plotted numerically in Fig. 6, which shows that the bounds in (4.27) and the upper bound in (4.26) hold for all $r\in(0,1)$. The lower bound in (4.26) is immediate from property 5 in Proposition 3.1. Indeed, we see that the bifurcation point exists and that $K^{*}$ is unique given $r$. We will show later that $r^{*}$ is also unique given $K$ by showing that $\frac{\partial r^{*}}{\partial K}>0$ in Theorem  4.4. ∎

The uniqueness of the bifurcation point corroborates the picture in Fig. 4.

Remark 4.2**.**

Note that (4.2) can also be solved for $L^{*}=L^{*}(K)$. The way this should be understood is that, after one of the variables $K$ and $L$ is fixed, the bifurcation point for the other variable is determined. A plot of the bifurcation point as a function of $K$ and $L$ is shown in Fig. 7.

4.2 Ordering of non-symmetric solutions

Due to the symmetry of the system, if $(r_{1},r_{2})$ is a solution to (3.11) with $\psi=0$, then so is $(r_{2},r_{1})$. When non-symmetric solutions exist, we have the following ordering of the synchronization levels in the two communities.

Theorem 4.3** (Ordered solutions).**

Fix $L$ and take $K>K^{*}$ where $K^{*}$ is the bifurcation point obtained by solving (4.2). Furthermore take only positive solutions so that $r_{1},r_{2},r>0$. Without loss of generality, consider a non-symmetric solution with $r_{1}>r_{2}$. Then

[TABLE]

Proof.

The symmetric solution $r$ solves the equation

[TABLE]

To prove that $r<r_{1}$, we consider the self-consistency equation (3.1) for $r_{1}$,

[TABLE]

and recall that we must have $L<0$ for non-symmetric solutions to exist. Since $\frac{r_{2}}{r_{1}}<1$, we know that $K+L\frac{r_{2}}{r_{1}}>K+L$ and, due to the fact that $x\mapsto V(x)$ is strictly increasing, also $r<r_{1}$. Note that we are not quantifying the difference $r_{1}-r_{2}$. The strict inequality follows purely from the fact that $\frac{r_{2}}{r_{1}}<1$, making it impossible to match the solutions for $r$ and $r_{1}$. Similarly, we can show that $r_{2}<r$. ∎

4.3 Properties of the bifurcation line

We cannot solve (4.11) analytically for $K^{*}$. We can, however, plot (4.11) numerically, which refines the phase diagram in Fig. 2 for $\psi=0$, as shown in Fig. 8. In this section we first list some basic properties of $r^{*}(K)$ and its derivatives, defined as the solution of (4.11) when we eliminate $L$ with the help of (4.18). After that we state a theorem on the asymptotic properties of the bifurcation line $L^{*}(K)$ defined implicitly by (4.2).

Theorem 4.4** (Properties of $K\mapsto r^{*}(K)$).**

$\lim_{K\downarrow 2}r^{*}(K)=0$. 2. 2.

$\lim_{K\rightarrow\infty}r^{*}(K)=1$. 3. 3.

$r^{*}(K)\sim\sqrt{\frac{K-2}{2}}\quad\text{as }K\downarrow 2$. 4. 4.

$1-r^{*}(K)\sim\frac{1}{2\sqrt{K}}\quad\text{as }K\rightarrow\infty$. 5. 5.

$\frac{\partial r^{*}(K)}{\partial K}>0\quad\text{for all }K>2$. 6. 6.

$\frac{\partial^{2}r^{*}(K)}{\partial K^{2}}<0\quad\text{for all }K>2$.

Proof.

We use (4.19) to get

[TABLE]

from which property 1 follows by taking the limit $K\downarrow 2$. The inequality in (4.19) also implies

[TABLE]

If both $K\mapsto r^{*}(K)$ and $r\mapsto K^{*}(r)$ are continuous, then property 2 follows. To show that $K\mapsto r^{*}(K)$ is continuous, we apply the implicit function theorem (IFT) to calculate the derivative

[TABLE]

which we find by rewriting (4.20). From the conditions for the IFT [18], we have that, in order for $K\mapsto r^{*}(K)$ to be continuous, we need that

[TABLE]

which we obtain by differentiating $h(K,r)$ with respect to $r$ and setting the derivative to zero. From this we obtain the following bound on $r^{*}(K)$:

[TABLE]

In order to rigorously show that this bound is satisfied, we can use the sequence of (iteratively defined) upper bounds $xu_{1}^{(k)}(x),k\in\mathbb{N}_{0}$, for $V(x)$ given in [24, Theorem 4], which converge to $V(x)$ as $k\rightarrow\infty$. Here we will use

[TABLE]

as suggested in[24, Equation (22)]. If substitution of the right-hand side of (4.35) for $r$ into $f_{r}(K)ru_{1}^{(k)}(f_{r}(K)r)-r$ makes it less than 0, then we know that the bound in (4.35) is satisfied. To see why, note that then

[TABLE]

Now, if $f_{r_{-}^{*}}(K)r_{-}^{*}u_{1}^{(k)}(f_{r_{-}^{*}}(K)r_{-}^{*})-r_{-}^{*}<0$ for all $K$, then so is $h(K,r_{-}^{*})$, so that $r_{-}^{*}$ does not satisfy $h(K,r^{*})=0$ and the solution satisfies $r^{*}>r_{-}^{*}$. Using $xu_{1}^{(k)}(x)$ with $k=2$,

[TABLE]

as an upper bound, we get that the bound in (4.35) is at least satisfied for $K\in(2,K_{k=2})$, where $K_{k=2}=15.8684$. By increasing $k$, we see that the upper bound of this interval increases and we expect that in the limit as $k\rightarrow\infty$, (4.35) is satisfied on $K\in(2,\infty)$. Numerically, we indeed see that this bound is satisfied, as shown in Fig. 9 (this figure shows that $K^{*}(r)<\frac{2-r^{2}}{(1-r^{2})}$, which is the same as (4.35)). For the continuity of $r\mapsto K^{*}(r)$ we require, again by the condition for the IFT, that

[TABLE]

which is satisfied by all $K^{*}>0$ and $r^{*}\in(0,1)$, so that property 2 is proved.

We know that $\lim_{K\downarrow 2}r(K)=0$ (by property 1), so that we can expand $V$ around 0 in the self-consistency equation (4.20). This leads to

[TABLE]

The corresponding asymptotic equation can be solved for $r^{*}(K)$ to obtain property 3. Property 4 follows from a similar calculation, by using the expansion of $V$ around infinity, and gives

[TABLE]

This equation gives rise to a cubic polynomial in $r^{*}(K)$, which can be solved and gives

[TABLE]

where

[TABLE]

The complex parts in the right-hand side of (4.42) compensate one another, making it real. Taking only the leading order terms in $K$, we obtain the asymptotics in property 4. We can calculate $\partial_{K}r^{*}(K)$ by differentiating (4.20), i.e.,

[TABLE]

where in the right-hand side we have abbreviated $r^{*}=r^{*}(K)$, and

[TABLE]

It follows from (4.12) that

[TABLE]

which simplifies (4.44) to

[TABLE]

Due to the inequality in [4, Equation (2.4)], we have that

[TABLE]

which makes the numerator positive. The denominator becomes zero when

[TABLE]

Rewriting the lower bound for $r^{*}$ in (4.35), we get

[TABLE]

which ensures that the denominator of (4.44) is never zero. For values of $K$ satisfying (4.50) the derivative is positive. This we find by substituting a pair of values $r^{*}(K),K$, calculated numerically, into (4.44), and proves property 5 because the derivative does not change sign in the range of $K$. To prove property 6, we take the derivative with respect to $K$ of (4.47) and substitute the expression for the first derivative. This leads to a lengthy equation with denominator

[TABLE]

which is positive by the same argument as for the first derivative. Setting the numerator to zero and solving for $K$, we find that there are no solutions when $r$ is between zero and the appropriate root of a 9th order polynomial in $r$, which numerically is 0.946819. Between this value and 1 there are two solutions, for which the numerator is zero, given by the solutions to the two roots of a quartic polynomial in $K$. We can plot these solutions together with the upper and lower bounds for $K^{*}(r)$ and compare them with the true $K^{*}(r)$, calculated numerically, as shown in the right panel of Fig. 9.

The right panel of Fig. 9 suggests that the second derivative also does not change sign. Numerically solving for a pair $(K,r^{*}(K))$, and substituting this into the numerator, we see that the second derivative is negative. This is confirmed by the left panel of Fig. 9. ∎

To confirm the asymptotic solutions for $r^{*}(K)$ in properties 3 and 4, we plot them and compare them to the numerical solutions in Fig. 10.

The next theorem gives the asymptotics of $L^{*}(K)$ implicitly defined by (4.2) in the limit as $K\rightarrow\infty$ and close to $(K,L)=(2,0)$.

Theorem 4.5** (Asymptotic properties of the bifurcation line).**

The derivative of $L^{*}(K)$, defined implicitly by (4.2), has the following properties:

$\lim_{K\rightarrow\infty}\frac{\partial L^{*}(K)}{\partial K}=-1$. 2. 2.

$\lim_{K\downarrow 2}\frac{\partial L^{*}(K)}{\partial K}=-\frac{1}{2}$.

Proof.

We begin by proving the existence of the limits, for which we need the following lemma.

Lemma 4.6** (Derivatives of $K\mapsto L^{*}(K)$).**

For all $K>2$,

$\frac{\partial L^{*}(K)}{\partial K}<0$. 2. 2.

$\frac{\partial^{2}L^{*}(K)}{\partial K^{2}}<0$.

Proof.

In order for $L^{*}(K)$ to be continuous by the IFT (in a similar way as in the proof of Theorem 4.4), we require that

[TABLE]

We will see, numerically, that this bound is satisfied because it lies below another lower bound of $L^{*}(K)$. To rigorously show that this bound is satisfied, we expect that it is possible to use the same procedure as outlined for the bound on $r^{*}(K)$ in (4.35). Now we start by differentiating (4.2) with respect to $K$ and solving for $\partial_{K}L^{*}(K)$. This leads to

[TABLE]

Setting the numerator, which is a quartic polynomial in $L^{*}(K)$, equal to zero and solving for $L^{*}(K)$, we find one solution that lies above the critical condition for $L$ when fixing $K$, $-K+2$. The expression is too lengthy to present here and does not lead to any useful insight. Taking the derivative with respect to $K$ of (4.53), substituting the expression for the first derivative (4.53) and setting the resulting numerator to zero, we are left with solving a 7th order polynomial for $L^{*}(K)$. Again the expression is lengthy and does not lead to any insight. Only one of the solutions to the 7th order polynomial lies above the critical line. Comparing these two solutions, one coming from the quartic polynomial and the other from the 7th order polynomial, we see numerically that the first is a lower bound for $L^{*}(K)$ and the second is an upper bound for $L^{*}(K)$, as seen in the right panel of Fig. 11. This lower bound is an upper bound for the right-hand side of (4.52), so that the conditions for the IFT are satisfied. The expression determining when the denominator of both the first and the second derivative is zero, obtained by setting their respective denominators to zero (which makes the derivatives diverge), is the same, and the only solution falling above the critical condition is upper bounded by the lower bound for $L^{*}(K)$ found above (as the solution to the quartic polynomial), so that the derivatives do not diverge.

The right panel of Fig. 11 suggests that both the first derivative and the second derivative of $L^{*}(K)$ do not change sign as a function of $K$. Substituting a pair of values $K,L^{*}(K)$, solved for numerically, we confirm the statements in Lemma 4.6. This is also corroborated by the left panel of Fig. 11. ∎

Remark 4.7**.**

For the mathematical reader the numerical assistance in the argument above might not be satisfying. We suspect that the proof can be made rigorous by using the sequences of upper and lower bounds in [24, Theorem 4] on $V$ in (4.2), in order to get upper and lower bounds for $L^{*}(K)$ that give a tighter wedge than the one in the right panel of Fig. 11.

Due to Lemma 4.6 and the fact that $L^{*}(K)$ is bounded below by $-K+2$, we have that the limits exists.

We now turn to the proof of Theorem 4.5. Abbreviate

[TABLE]

By the implicit function theorem, we have

[TABLE]

Compute

[TABLE]

and

[TABLE]

where $S(x)=\frac{\text{I}_{2}(x)}{\text{I}_{0}(x)}$.

For property 1, we make the Ansatz $L^{*}(K)=-aK+c,K\rightarrow\infty$ where $c=c(K)=o(K)$ (which is confirmed in Fig. 8). Taking the limit $K\to\infty$, we get zero for the first terms in the right-hand sides of (4.3)–(4.3), i.e.,

[TABLE]

where we have used the expression for $r^{*}(K,L)$ from (4.3). The multiplication factors in the last line of the right-hand sides of (4.3)–(4.3) are the same, so we are left with calculating the limit as $K\rightarrow\infty$ of the quotient

[TABLE]

A straightforward but tedious calculation (with the help of MATHEMATICA) shows that this limit is $-1$.

For property 2, we must find the limit of $-\frac{\partial_{K}g(K,L)}{\partial_{L}g(K,L)}$ as we approach the point $(K,L)=(2,0)$ along the line $L^{*}(K)$. We make the Ansatz $L^{*}(K)=(K-2)b+o(1)$, $K\downarrow 2$. Making this replacement in the expression for the derivative and doing a Taylor expansion around $K=2$, we obtain after a tedious calculation (with the help of MATHEMATICA),

[TABLE]

for the terms in the numerator and

[TABLE]

for the terms in the denominator. Combining (4.61)–(4.62) we obtain

[TABLE]

so that $b=-\frac{1}{2}$. ∎

Properties 1 and 2 are confirmed by the left panel of Fig 11. It seems that $K\mapsto L^{*}(K)$ for large $K$ does not have an asymptote, since when we take the limit after the replacement $L^{*}(K)=-K+c$ we get an equation for the bifurcation point that reads

[TABLE]

The only solution to this equation is $c=1$, which is not possible because it would place the asymptote below the critical line. This suggests that $c=c(K)$ grows as a function of $K$, but that this growth is sublinear.

5 Simulation

Fixing the phase difference is not physical, since the system will relax into a steady state and will choose the angles that are the least costly energetically. Studying the dynamics of the transitions between states or the stability properties of the possible states are both difficult tasks. However, we expect that the non-symmetric state is either unstable or metastable and using simulations we can observe what type of transitions one might expect between the possible states. To see this, we take the initial distribution for both populations to have mean $\pi$, but choose the second community to have a slightly larger variance initially, meaning that the synchronization level starts lower. The outcome of the simulation can be seen in Fig 12. It seems that the community with less synchronization initially is suppressed by the community with more synchronization, until the ‘push’ from the latter becomes too strong. This is reflected in the angles, which stay relatively close together for a while, before moving apart. This type of transition seems only to occur when the parameters are chosen such that the non-symmetric solutions discussed above exist.

We expect that the most stable state is the symmetric solution with the largest synchronization level (i.e., the largest effective interaction strength). For example, if $K=5$ and $L=2$, then the symmetric solution with phase difference $\psi=\pi$ has $r_{1}=r_{2}=0.724\ldots$, while the symmetric solution with phase difference $\psi=0$ has $r_{1}=r_{2}=0.918\ldots$. The first state is unstable/metastable, the second state is stable. The transition from the one to the other is shown in Fig 13.

To be clear, these simulations are not meant to prove any stability properties or transitions, but are useful for determining what types of stability properties or transitions we may expect. They also suggest that much interesting work remains to be done.

Appendix A Concavity of V

Recall that

[TABLE]

The first derivative of (A.1) is

[TABLE]

We can rewrite

[TABLE]

where we have used Holder’s inequality in the second line. Taking the square on both sides, we obtain

[TABLE]

which proves that (A.2) is non-negative. We evaluate (A.2) at $x=0$, to get

[TABLE]

For the second derivative we rewrite

[TABLE]

where $f(x,\theta)=\mathrm{e}^{x\cos\theta}$ and the prime refers to the derivative with respect to $x$. The integrals are always from [math] to $2\pi$. Taking the derivative of the first term, we find

[TABLE]

while for the second we find

[TABLE]

Using a common denominator, we can write the difference as

[TABLE]

To continue, we first let

[TABLE]

making the desired expression into

[TABLE]

With the functions

[TABLE]

we can perform the change of variable $u=\cos\theta$, i.e., $\theta=\arccos u$ and

[TABLE]

Here we get

[TABLE]

where we have used $\arccos_{1}$ for the first integral and $\arccos_{2}$ for the second. Note that, when $k=0$,

[TABLE]

With the change of measure

[TABLE]

we obtain

[TABLE]

and defining $m=m(x)=\int_{-1}^{1}ue^{xu}\,\mathrm{d}\nu(u)$ we get

[TABLE]

Note that, due to (A.13), and since $\int_{-1}^{1}3m^{2}\,u\,e^{xu}\,\mathrm{d}\nu(u)=3m^{3}$, (A.16) equals

[TABLE]

We can check this by writing

[TABLE]

To complete the proof we state [22, Lemma 4], suitably adapted.

Lemma A.1**.**

Let $\nu$ be an even probability measure with support on $[-1,1]$, and suppose that $\nu$ is absolutely continuous, i.e., $\mathrm{d}\nu(\sigma)=f(\sigma)\mathrm{d}\sigma$, with $f$ non-decreasing on $[0,1]$. Then $\nu\in\mathscr{P}$, the class of all probability measures on $\mathbb{R}$ with compact support, is such that

[TABLE]

where

[TABLE]

and $k=Jm+h$, where $J$ is the mean-field interaction strength and $h$ is the magnetic field strength of the spin system.

With the identification

[TABLE]

and taking out a negative, we complete the proof. To get the strict inequality we note that the equality in the lemma does not hold for our choice of $\nu(u)$.