OPINION DYNAMICS ON A GENERAL COMPACT RIEMANNIAN MANIFOLD

. This work formulates the problem of deﬁning a model for opinion dynamics on a general compact Riemannian manifold. Two approaches to modeling opinions on a manifold are explored. The ﬁrst deﬁnes the distance between two points using the projection in the ambient Euclidean space. The second approach deﬁnes the distance as the length of the geodesic between two agents. Our analysis focuses on features such as equilibria, the long term be- havior, and the energy of the system, as well as the interactions between agents that lead to these features. Simulations for speciﬁc manifolds, S 1 , S 2 , and T 2 , accompany the analysis. Trajectories given by opinion dynamics may resem- ble n − body Choreography and are called “social choreography”. Conditions leading to various types of social choreography are investigated in R 2 .

. Similarly, applications to satellite or ground vehicle coordination have motivated the development of models on special orthogonal groups [16,17]: satellites evolve on SO(3) while ground vehicles evolve on SE (2) or SE (3). A nonlinear model of opinion formation on the sphere was also developed in [3]. A discussion of the so called "flocking realizability problem" for a sphere is given in [5], which focuses on the particular equilibrium that we refer to as consensus for a holonomic dynamical system on a sphere.
The present work defines a general model of opinion dynamics on a Riemannian manifold. We investigate how the manifold on which the model is defined affects the global configurations resulting from opinion dynamics. These are the first steps to build a robust theory of opinion dynamics on general Riemannian manifolds.
There is an inherent difficulty in defining opinion dynamics on a general Riemannian manifold. Using the Riemannian distance, an agent will move towards a point by following the manifold's geodesics, which are well defined only locally. On a larger scale, there might not exist a unique geodesic. Another challenge is the extreme complexity of computing geodesics, even on a relatively simple manifold such as the torus [8]. One way around this issue is to consider the embedding of the manifold into a Euclidean space. Each agent's velocity is defined by projection of the other agents' influence onto the tangent space at that point. This is the choice made in [3].
Other than the mentioned practical aspect, there is an intrinsic rationale for choosing one approach over the other. When evolving along the geodesics of the manifold, one assumes that each agent has a global understanding of the manifold's geometry and is able to choose the shortest path among all possible ones. On the other hand, the approach based on the projection of the desired destination onto the tangent space implies that each agent only holds local information about the space in which it evolves. It chooses to move in the direction which locally seems to bring it closer to the target.
We explore these two specific approaches for our generalized model. The first method, Approach A, uses projections in the Euclidean space in which the manifold is embedded. The second method, Approach B, uses only geodesics defined on the manifold to define strength and direction of interaction. We exhibit properties of the interaction matrix that lead to specific kinds of equilibria. Simulations and examples compare the two methods. Dancing equilibria for Approach B are shown (dancing equilibria were studied for Approach A in [3]).
We use the sphere and torus as sample manifolds to evaluate these approaches. Specifically, we simulate dynamics on the following manifolds: S 1 , S 2 and T 2 . These examples allow us to directly compare the two approaches, and see if one is more appropriate for a given manifold. We show the influence of the manifold's geometry on the dynamics by examining the dynamics resulting from the same interaction matrix in S 2 , T 2 and R 2 .
Opinion dynamics trajectories can resemble n−body choreography, that is, solutions to the well known n−body problem. These dynamics drive agents along orbits which either are periodic, or have a periodic feature, and that may be shared by multiple agents. We refer to opinion dynamics trajectories along such orbits as "Social Choreography". We show that a simple example of Social Choreography in R 2 does not hold on S 2 or T 2 , see Figures 10 and 11. We then focus on the case of R 2 and investigate initial conditions and properties of the interaction matrix which give rise to Social Choreography. Future investigations can consist of exploring similar properties of trajectories on general manifolds. 2. Choice of the model. This work will primarily discuss two approaches to define opinion dynamics on a Riemannian manifold. Let M be a Riemannian manifold. Let N ∈ N represent the number of agents with opinions evolving on M . We denote by x := (x i ) i∈{1,...,N } ∈ M N the set of opinions. For each i ∈ {1, ..., N }, x i ∈ T xi M . The opinions x i evolve according to the following general dynamics: where • a ij ∈ R is the interaction coefficient of the pair of agents i and j, • Ψ : R → R is the interaction potential, • d(·, ·) : M × M → R + represents the distance between opinions, • ν ij ∈ T xi M is a unit vector giving the direction of the influence of j over i. Each of these terms is further specified in the following.
2.1. Approaches. The evolution of each agent's opinion depends on the opinions of all other agents, with influences weighted by the interaction coefficients a ij . More specifically, an agent x j 's influence on x i is determined by two elements: the direction of influence ν ij ∈ T xi M and the magnitude of influence Ψ(d(x i , x j )) ∈ R + . We propose and study two different approaches for the choices of d and ν ij . Approach A uses the embedding of M in R n to define d(x i , x j ), whereas Approach B is intrinsic to M , with distance and direction of influence based on geodesics. Approach A. Assume that M of dimension m is embedded in a Euclidean space R n , with n ≥ m. Agent x j acts on agent x i via a projection onto T xi M ⊂ R m . Now considering points (x i , x j ) ∈ M 2 as points of R n , the difference x j − x i is a vector of R n . Given a vector subspace Y of R n , we denote by Π Y y the projection of y ∈ R n onto Y ⊂ R n and define d P (·, ·) as follows: where · denotes the Euclidean norm on R n . The same projection also defines the direction of influence of x j on x i : With the specific choice Ψ ≡ Id, system (1)-(2)-(3) becomes: This is the approach used in [3], applied to the sphere S 2 .
Notice that the magnitude of influence, d P (x i , x j ), is symmetric for the sphere in the sense that d P (x i , x j ) = d P (x j , x i ), but not symmetric for a general Riemannian manifold (see Figure 1). However, it is a continuous function defined for all pairs of points (x i , x j ) ∈ M 2 . The originality of this approach is that the influence of x j on x i is not related to a notion of distance between the points. The use of the projection of x j − x i onto T xi M reflects the concept of "local visibility." For the situation of 492 AYLIN AYDOGDU, SEAN T. MCQUADE AND NASTASSIA POURADIER DUTEIL (4), an agent is subject to "local visibility", and movement of x i along T xi M (dashed line through x i ) will not bring x i closer to x j in this local sense.
two agents evolving on a one dimensional manifold, if x j − x i ⊥ T xi M , then a local displacement of x i does not affect the distance between the points x i −x j . Indeed, a first order Taylor expansion gives: Supposing that x j is fixed, we have: Hence if x i only has local visibility, all directions of displacement seem equivalent (at first order), which justifies the influence of x j over x i to be zero if their difference is orthogonal to the tangent space of M at x i . This is illustrated in Figure 1. Approach B. This second approach defines d and ν ij using the manifold M itself, and does not require any reference to the space in which M is immersed. This would make Approach B a natural way to define system dynamics, however the complete knowledge of the geodesics between any two points on the manifold may be unrealistic. Furthermore, the geometry of the manifold may introduce difficulties to the uniqueness of ν ij , particularly at the cut-locus of a point.
Definition 2.1. The cut locus of a point q ∈ M is the set of points CL(q) ⊂ M for which there are multiple geodesics between q and p ∈ CL(q) (see also [4]).
Let γ ij : [0, 1] → M denote a geodesic connecting x i to x j , γ ij (0) = x i and γ ij (1) = x j . We then define the distance between x j and x i as the length of a geodesic, i.e. denoting by g y : T y M × T y M → R + the Riemannian metric at point y ∈ M , The direction of influence is determined by the same geodesic: otherwise.
Unlike in Approach A, the magnitude of influence is a symmetric function: Furthermore, this approach ensures that the magnitude of influence of one agent on another is a function of the exact Riemannian distance between the agents. Interaction networks. In finite-dimensional systems such as system (1), the set of interacting agents can be described by vertices of a graph. A directed edge exists from a vertex i to a vertex j if and only if a ij = 0. The system depends on the interaction network, and likewise, if the coefficients a ij are chosen to be functions of the state, the interaction network may change as a result of the dynamics. Two main types of interaction networks have been proposed in the literature: metric interactions and topological interactions. If interactions between agents occur only locally, only the neighbors of agent i influence agent i. Metric interactions define the set of neighbors of agent i, given a radius r > 0, as where d(·, ·) can represent either the projection or the geodesic distance, as specified in each of the two approaches described above (see equations (2) and (6)). The other main type of interactions specifies that an agent is influenced by only its k closest neighbors. We call these topological interactions [1]. We define the relative separation between two agents as α ij = card{k : The set of neighbors of agent i is then defined as the set of its k closest neighbors, i.e. for a given k ∈ N, Figures 2 and 3 illustrate differences between the metric and topological networks for the specific example of S 1 , with each of the approaches A and B.
x 3 x 4 x 1 x 2 x 5 x 6 (a) Approach A metric, r = π 4 Tx 1 S x 3 x 4 x 1 x 2 x 5 x 6 (b) Approach A topological, k = 2 x 3 x 4 x 1 x 3 x 4 x 1 x 2 x 5 x 6 π 4 (d) Approach B topological, k = 2  x 1 x 2 x 3 x 4 x 5 x 6 (a) Approach A metric, r = π 4 x 1 x 2 x 3 x 4 x 5 x 6 (b) Approach A topological, k = 2 x 1 x 2 x 3 x 4 x 5 x 6 (c) Approach B metric, r = π 4 x 1 x 2 x 3 x 4 x 5 x 6 (d) Approach B topological, k = 2 Figure 3. The agent x 1 is influenced by different agents depending on how the interaction network is defined. These networks may change as the dynamics move the agents on S 1 . Each agent x j , j ∈ {1, . . . , 6} will have a network describing which other agents influence x j . The interaction networks corresponding to systems from Figure 2.
Resolution of discontinuities. The definitions of ν ij for approaches A and B (given by equations (3) and (7)) allow discontinuities of ν ij at certain points. Thus, one must impose conditions on the interaction potential Ψ ∈ C 0 (R + , R + ), in order to ensure the continuity of the right-hand side of the system (1), and hence the existence and uniqueness of a solution. Table 1 lists the discontinuities of ν ij and gives necessary conditions on Ψ to ensure the continuity of Ψ(d(x i , x j ))ν ij . Firstly, notice that in both approaches, ν ij is discontinuous at the point x i = x j . Indeed, if x i = x j , ν ij = 0, whereas almost everywhere else, ν ij = 1. To ensure the continuity of Ψ(d(x i , x j ))ν ij at this point, we impose the following condition: In Approach A, we created a discontinuity of ν ij at the points x j ∈ N (x i ), where we denote by N (q) the set N (q) := {q ∈ M | Π TpM (q − p) = 0}. For convenience of notation, we will use interchangeably the notations N (x i ) and N i . More specifically, we have lim xj →Ni ν ij = 1 but ν ij = 0 if x j ∈ N i (see also Table 1). However, from the definition of d P (see equation (2)), we have lim xj →Ni d P (x i , x j ) = 0 and d(x i , x j ) = 0 for x j ∈ N i . Hence a sufficient condition for Ψ(d(x i , x j ))ν ij to be continuous is again: In Approach B, there is a discontinuity for x j ∈ CL(x i ). Denoting by B geo (p, ρ) the geodesic ball of center p and radius ρ, we require the following condition on the influence function Ψ: where := inf{ρ > 0 | ∀p ∈ M, B geo (p, ρ) ∩ CL(p) = ∅}. This distance , also known as injectivity radius, is known to exist and be greater than 0 for any compact Riemannian manifold (see [4]).

Approach A B
A and B Critical points Notice that in the case of the geodesics approach (Approach B), the condition Ψ(d) = 0 for all d ≥ is incompatible with the use of the topological network (9). Indeed, if agent j is among the k closest neighbors of agent i, the topological network would require: a ij = 0. However, the interaction between i and j would be canceled if d G (x i , x j ) > . On the other hand, the metric interaction network as defined by (8) is compatible with Approach A, and with Approach B if the interaction radius is smaller than the injectivity radius: r ≤ . For simplicity purposes, in the rest of this paper, we will consider that the interaction coefficients a ij are constant, thus not requiring the need to differentiate between metric and topological networks. While models with constant interaction coefficients are our focus here, these models are quite restrictive, and exclude all models with dynamic interactions.
This is a system of at most N d equations in the N 2 − N unknowns a ij , i = j, notice that Ψ(d(x i , x i )ν ii = 0, and diagonal values of A do not change the system. So if N > d + 1 there exists a nontrivial choice of the interaction coefficients for whichx is an equilibrium.

AYLIN AYDOGDU, SEAN T. MCQUADE AND NASTASSIA POURADIER DUTEIL
The kinetic energy of System (1)-(6)- (7) is the quantity Proposition 3. Let M be a general Riemannian compact manifold. Consider the dynamics given by projection onto the tangent space (Approach A) given by (4).
Indeed, notice that Then we compute where the third equality uses the property: a ij = a ji for all i, j.
This shows that d dt F (t) → 0 when t → ∞, which implies that lim t→∞ E P (t) = 0. Remark 2. Proposition 3 assumes that Ψ ≡ Id which creates a discontinuity for Approach B (see Table 1). A result is shown for the more restricted case of M = S 2 and Ψ(·) ≡ sin(·), (Corollary 1). Definition 2.4. Let x solve the differential equation (1). A dancing equilibrium is a configuration in which for all pairs of agents (i, j), In the context of a system of oscillators, this equilibrium is also known as a phase-locked state or entrainment state [10].
Remark 3. This definition is a generalization of the concept of dancing equilibrium described in [3].
Remark 4. It follows immediately from definition 2.4 that the kinetic energy of a system in dancing equilibrium is constant.
3. Analysis and simulations on S 1 .

Models. We study both approaches A and B in the case
Approach A. The projection onto an agent's tangent space can be rewritten as: So System (1)-(2)-(3) becomes: for alli ∈ {1, . . . , N }, where sgn(·) is the sign function defined by: We can then specify: This gives the system of scalar equations: In particular, in the case Ψ ≡ Id, the system becomes the Kuramoto model [12].
Approach B. For M = S 1 , the geodesics distance d G and the vector ν G ij are given by: System (1)-(6)-(7) is written: In order for the system to be well defined, the interaction function Ψ must satisfy the conditions given in Table 1. Notice that the injectivity radius is constant over S 1 , with = π. Possible choices involve choosing Ψ from a family of function defined as follows: where a ∈ (0, π) (see Figure 8). Another possible choice is: Ψ : x → sin(x). Notice that for the specific choices Ψ = Id for Approach A and Ψ : x → sin(x) for Approach B, the two approaches A and B are equivalent.

3.2.
Analysis. We first examine the different equilibria for both approaches.
Proof. Using the hypotheses from Theorem 3.1, we can easily compute: Interestingly, Theorem 3.1 is not applicable to Approach B. We illustrate the different behaviors of the two systems by studying the specific example of four agents initially in a rectangular configuration. According to Theorem 3.1, this configuration is an equilibrium for Approach A, independently of the choice of interaction function Ψ. However, one can easily prove that in the geodesics-based Approach B, with N = 4 and the choice Ψ := Ψ a with a = 3π 4 , the only equilibrium for which all agents have pairwise distinct positions is obtained by a regular polygon, i.e. all agents are evenly spaced out on the circle. This is illustrated by numerical simulations shown in Figure 4.
This highlights the fundamentally different behaviors of the systems (1)-(2)- (3) and (1) (26)). Notice that with Approach A, initial and final positions are identical since any rectangle configuration is an equilibrium. However, with Approach B, the system reaches a square configuration, the only possible equilibrium with pairwise distinct positions.
In both approaches A and B, conditions on the interaction matrix A can be found such that the system forms a dancing equilibrium (see Definiton 2.4).
Theorem 3.2. Consider the dynamics on S 1 given by: where d(·, ·) and ν are given either by Approach A (21) or Approach B (24). Let C ∈ R and suppose that for all i ∈ {1, . . . , N }, Then the system is in a dancing equilibrium.
Proof. If the interaction matrix satisfies (28), then at t = 0,   4.1. Models. We study both approaches A and B for M = S 2 , i.e. for a two dimensional sphere embedded in R 3 . We use spherical coordinates: let (θ i ) i∈{1,...,N } ∈ [0, 2π] N , and (φ i ) i∈{1,...,N } ∈ [0, π] N such that for all i ∈ {1, ..., N }, Choice of influence function. We choose an influence function Ψ(d) between two agents x i and x j so that the right-hand side of the system is continuous, the discontinuities are shown in Table 1. For Approach B, the only point in CL(x i ) for a given x i is the antipodal point (this is an end point of a diameter for which x i is the other end point.) As in the case of S 1 , for Approach B, we choose a function Ψ from a family of functions of the form Ψ a , see equation (26) (choices of Ψ are shown in Figure 8). Approach A. On S 2 , the derivative for system (1)-(2)-(3) with Ψ ≡ Id reduces to the sum of all projections onto the tangent space of agent x i , weighted by the corresponding interaction term a ij . This is rewritten as: (6) between two points x i , and x j on S 2 is given by: where · is the standard norm in R 3 . Noticing that for S 2 , Approach A with Ψ ≡ Id is equivalent to Approach B with Ψ ≡ sin, we extend the results of Proposition 3. Corollary 1. Consider the dynamics given by geodesic distance (Approach B) on S 2 , system (1)-(6)- (7), and let Ψ ≡ sin. If the interaction matrix A = (a ij ) i,j∈{1,...,N } 2 is symmetric, then Proof. System (1)-(6)- (7), with Ψ(·) ≡ sin(·) reads: Considering the embedding of the system in R 3 , we notice that for all i, j ∈ {1, . . . , N }, thus the system is (4), and by proposition 3, lim t→∞ E P (t) = 0. Finally,

4.2.
Simulations. We use a fourth order Runge-Kutta scheme to approximate the trajectories. The derivative is calculated as a vector in R 3 , and then we express this vector in spherical coordinates,θt θ andφt φ where t θ and t φ are the unit vectors in the direction of the azimuth angle(θ) and the polar angle(φ) respectively for the ith agent. {t θ , t φ } form an orthonormal basis for the tangent space of x i . Using the angular derivatives avoids having to calculate the agent's trajectory in R 3 and then project onto the sphere for every iteration which would cause significant numerical errors.
For an agent x i = (θ i , φ i ) we can write t θi and t φi as We express the derivative of an agent aṡ By direct computation, we get: We can also express the derivative of x i as the projection of the derivative in R 3 with t θi and t φiẋ It follows from (33) and (34) thaṫ Singularities: In (35), the factor 1 sin φ causes a singularity around φ = kπ for a nonnegative integer k. To avoid this practical problem, before each iteration of the RK4 scheme, we identify critical agents that have a polar angle close to 0 or π (φ = kπ), for non-negative integer k; for these critical agents, we rotate all agents π 2 around the x-axis to calculate the derivative. This is a concern for both Approach A and Approach B.
An additional concern is singularities for agents forming consensus. In equations (3) and (7), ν ij is normalized, and when agent i and agent j are very close together, dividing by d P and d G causes a singularity. We avoid this problem in our simulations by defining a minimum distance d min between agents for the sake of the normalization term in the denominator. When d(x i , x j ) < d min then Examples. We ran simulations using different choices of Ψ to see how this choice can impact the system. We show two simulations, the first uses Approach A ( Figure 7); the second uses Approach B (Figure 9). In both pictures the same interaction matrix is used (given in the appendix), and we see that our choice of Ψ ( Figure 8) may dramatically change the system behavior. The second example shows the effect of curvature on the system ( Figure 10) for comparison to T 2 and R 2 (Figure 11). Another example in the appendix shows unexpected behavior using Approach A.     To assess the influence of the curvature of S 2 on the dynamics, observe a simple case involving 3 agents evolving according to the interaction matrix: In Section 6.2, we prove that those dynamics in R 2 lead to periodic trajectories on a single orbit shared by all three agents, the orbit's parameters being fully determined by the initial conditions (see Theorem 6.4). However, the same dynamics on the sphere do not give rise to periodic trajectories. In sections 5.3, we also discuss the dynamics with this interactions matrix on T 2 , to assess the effect of curvature of the manifold.

5.
Analysis and simulations on T 2 . We now study how the general dynamics given by equation (1) apply to the specific case of the torus T 2 ⊂ R 3 . Let (e x , e y , e z ) denote the Euclidean basis of R 3 . Let (R, r) ∈ (R + ) 2 , with R > r. We define the manifold T 2 as the torus obtained by rotating the circle (x − R) 2 + z 2 = r 2 around the z-axis. Hence T 2 is defined by the equation (R − x 2 + y 2 ) 2 + z 2 = r 2 . The Figure 10. Dynamics with Approach A on S 2 , using the interactions matrix (36). If the agents' initial positions are close enough to each other, the agents with will form trajectories that remain in a neighborhood of their initial position.
parametric equations for such a torus are: The angles φ and θ are respectively referred to as the toroidal and poloidal angles. A set of points with the same toroidal angle is called a meridian.

5.1.
Model. We first investigate the behavior of system (1) with Approach B (using the geodesic distance) in the case of T 2 . Unlike in the cases of S 1 and S 2 presented in the sections 3 and 4, there exists no simple expression for the geodesic distance between two points on the torus. In 1903, Bliss studied and classified the different kinds of geodesic lines on the standard torus [2], using elliptic functions. Gravesen et al. determined the structure of the cut loci of a torus of revolution [8]. Several challenges arise when defining Approach B on T 2 . Firstly, computing the Riemannian distance between two points is highly non-trivial. One could consider approximating it numerically, but in the numerical discretization of equations (1)-(6)- (7), N (N − 1)/2 geodesics would have to be computed per time-step. That would require tremendous computing power.
Secondly, assuming that one is able to efficiently compute the geodesics on T 2 , one must take into account the cut-loci of each point to ensure that the dynamics (1)-(6)-(7) are well-defined. A method to guarantee well-defined dynamics would be to use a bounded confidence model [11], where the neighborhood of influence for an agent x i at point p is of smaller radius than the closest element in the cut locus of p. See section 2 for conditions on Ψ to make the right hand side of equation (1) continuous.
For simplicity, we thus focus on Approach A, where the dynamics are a function of the projection of each vector x j − x i onto the tangent space at x i . We will show that some restrictions still apply to the interaction function Ψ, but they are less restrictive and more easily determined than in Approach B.
Equations (1)-(2) reads: The vector ν ij depends on the influence that x j has over x i . It is zero if Π Tx i T 2 (x j − x i ) = 0, and it is a unit vector otherwise. Let N i be the set of points that have no influence on x i (see Table 1). Then, given i, j ∈ {1, . . . , N }, ν ij has the following expression: Let x i ∈ T 2 . We start by determining the set N i . For all i, we define the vectors u φi = cos φ i e x + sin φ i e y and u θi = cos θ i u φi + sin θ i e z , so that each agent's position vector reads: x i = Ru φi + ru θi . With these notations, u θi is the normal to the tangent space at the point x i . A basis for the tangent space at a point x i (φ i , θ i ) is given by the two tangent vectors t φi = (− sin φ i , cos φ i , 0) and t θi = (− sin θ i cos φ i , − sin θ i sin φ i , cos θ i ). Notice that x i , t φi = 0. Hence the condition Π Tx i T 2 (x j − x i ) = 0 reads: After computations, we get: If φ j = φ i , the second condition becomes: If φ j = φ i ± π, the second condition becomes: Notice that this last equation only has a solution if | sin θ i | ≤ r 2R . The set of positions that have no influence on x i thus comprises up to four points on the torus, depending on the values of r, R and sin θ i . We then have: r sin θ i |)), (−φ i , π − θ i + sgn(sin θ i ) arcsin(| 2R r sin θ i |))}. To ensure the continuity of the right-hand side of equation (37), one must impose the conditions of table 1.
We now go back to equation (37). We study the specific case where Ψ ≡ Id, which indeed satisfies (1). Then the system becomes: Hence the velocity reads: a ij x j is the sum of the influences of all agents on agent i. Notice that with the same notation, the system does not reduce to the simple formẋ i = α i − α i , x i x i for the same dynamics on the sphere (see [3]). This is due to the fact that on the torus, the position vector x i does not define the normal to the tangent space at x i , unlike in the cases of S 1 and S 2 .
The velocity of each agent is given by: (40) From (39) and (40) we get the angular velocities: Notice that unlike in the case of S 2 , see equation (35), here the derivativesφ i andθ i are not singular. This makes numerical simulations straightforward, not requiring the approximations described in Section 4.2.

5.2.
Properties. We now analyze the dynamics (1)-(2)-(3) on T 2 . We identify families of initial conditions that trivialize the dynamics. , then the torus dynamics simplify to the dynamics on S 1 given by (22) or (23).

Simulations.
To assess the influence of the curvature of the manifold on the dynamics, we compare a simple case involving 3 agents evolving according to the interaction matrix given in equation (36). As in the case of S 2 , the dynamics on the torus do not give rise to periodic trajectories (as opposed to the dynamics in R 2 , see Theorem 6.4). Instead, since T 2 can locally be identified with R 2 , if the initial mutual distances are small enough, the dynamics resemble those in R 2 . More specifically, the trajectories are quasi-periodic with a gradual shift of the center of mass (see Figure 11). However, if the initial distances between agents are large, the geometry and curvature of the torus changes radically the behavior of the system. 6. Social choreographies. As seen in Sections 3 and 5.3, when the interaction matrix A satisfies certain properties, for instance given by (28) on S 1 or by (36) in R 2 , then the trajectories exhibit special properties of symmetry or periodicity.
In [3], configurations on S 2 in which all mutual distances between agents remain constant were named dancing equilibrium. In this section, we investigate systems with similar properties of periodicity or symmetry. We use the term social choreography, drawing a parallel with the wellknown "n-body choreographies" discovered by Moore [13,14] in the context of point masses subject to gravitational forces. In the n-body problem, the interaction potentials between masses are predetermined, as they depend exclusively on the masses and distances between agents. Hence the conditions for a n-body choreography to occur only depend on the initial state of the system. In the case of social choreography, there are more degrees of freedom, as we design the interaction matrix as well as to set the initial conditions.
We study sufficient conditions on the interaction matrices for the trajectories of the system to be periodic or symmetric by focusing on the Euclidean space R 2 , with the specific choice of interaction potential Ψ ≡ Id. A future direction of this paper can consist in extending these results to general Riemannian manifolds. In R 2 and with Ψ ≡ Id, both approaches A and B are equivalent and the system simply reads as: We define the kinetic energy E = E G = E P as in Definition 2.3: A simple case of social choreography is that of a system with periodic trajectories, which we define as follows: ..N be a solution of (42). We refer to the system as having periodic trajectories if there exists τ > 0 such that for all i ∈ {1, ..., N }, for all t > 0, x i (t + τ ) = x i (t).
We will examine possible periodic behaviors of the system in sections 6.2, 6.3 and 6.4. 6.1. Rotationally invariant system. We now give sufficient conditions on the interaction matrix and on the initial conditions for the system to be invariant by rotation. denote the rotation matrix in R 2 for the angle θ ∈ [0, 2π). Suppose that initially, the system is invariant by rotation of angle 2kπ N , that is: Suppose that the interaction matrix A is invariant by change of basis, i.e. P −1 k AP k = A. Then the system remains invariant by rotation of angle 2kπ N at all time:

Proof. Let
can be rewritten as: From the theorem's hypotheses, Y (0) = Z(0). Let us show that Y and Z have the same evolution. One can easily prove that P −1 kÃ P k if and only if P −1 k AP k . Then notice thatẊ =ÃX = P −1 kÃ P k X. From that we compute: Since Y and Z satisfy the same differential equation and Y (0) = Z(0), then Y (t) = Z(t) for all t ≥ 0. This implies that at all time, 6.2. Unique orbit. Another example of social choreography is that of a system in which all agents share one unique orbit. Such choreographies have been discovered in the context of the n-body problem, for instance the "figure 8" orbit for three equal masses [13].  To illustrate Theorem 6.2, we study the evolution of N agents initially positioned at regular intervals on a circle, with an interaction matrix and initial conditions given by: Notice thatÃ = A, and the system satisfies the conditions of Theorem 6.2 with k = 1. Hence for all i ∈ {1, . . . , The 2N -dimensional system then reduces to a 2-dimensional one for the two coordinates x 11 and x 12 of x 1 , and all the other variables can be recovered by rotation of x 1 :ẋ This can be written as: Solving this linear system yields: x 11 (t) = x 11 (0) cos(2 sin( 2π N )t) − x 12 (0) sin(2 sin( 2π N )t) = cos(2 sin( 2π N )t) x 12 (t) = x 11 (0) sin(2 sin( 2π N )t) + x 12 (0) cos(2 sin( 2π N )t) = sin(2 sin( 2π N )t) .

OPINION DYNAMICS ON A GENERAL COMPACT RIEMANNIAN MANIFOLD 513
This proves that all agents share one common circular orbit, and their trajectories are periodic of period 2π(2 sin( 2π N )) −1 . Figure 15 provides a numerical illustration of this behavior, with 10 agents initially positioned at regular intervals on the unit circle. Another interesting example is that of 3 agents interacting according to the interaction matrix given previously, which, reduced to N = 3, gives: Theorem 6.4. Let N = 3. Consider the system (42) with interaction matrix given by (45). Then there exists a unique orbit shared by all agents, and all three trajectories are periodic.
Proof. The x and y-coordinates of the systems are decoupled, so that the 6-dimensional system can be reduced to two 3-dimensional ones. Notice thatÃ = A. Then for each coordinate j ∈ {1, 2}, the system reads: Due to the special structure of e tA , this can be rewritten as:

AYLIN AYDOGDU, SEAN T. MCQUADE AND NASTASSIA POURADIER DUTEIL
This shows that all three trajectories are periodic, or period 2π √ 3 . One can compute the positions of each agent after a third of a period and notice that: This shows that there is one unique shared orbit.

Coupled periodic trajectories. Other conditions on the interaction matrix
A give rise to different kinds of periodic behaviors. Here we provide sufficient conditions for the system to exhibit periodic trajectories, such that each orbit is shared by two agents.
Theorem 6.5 (Coupled periodic trajectories). Let N be even. Suppose that initially, the system is invariant by rotation of angle 4π N , that is: Let a, b > 0 and let (46) Then the system is periodic of period τ = π √ ab sin(2π/N ) . Furthermore, if N is divisible by 4, opposite agents share orbits two by two, i.e.: and the kinetic energy is periodic with period τ /2.
Proof. First remark that the system satisfies the hypotheses of Theorem 6.2, so Hence the system is entirely known from the positions of the first two agents, since all others can be obtained by simple rotations. We show that this 2N -dimensional problem can be rewritten as a 4-dimensional one. Indeed, using the fact that x N = R(−4π/N )x 2 and x 3 = R(4π/N )x 1 , the system . This can be rewritten in matrix form as:    ẋ 11 x 12 where One can easily show that this reduced interaction matrix A 4 has two purely imaginary conjugate eigenvalues, iλ and −iλ, each of multiplicity 2, where λ = 2 √ ab sin( 2π N ). Hence the solution of the system (47) can be written as a weighted sum of the functions t → cos(λt) and t → sin(λt). This implies that the system is periodic, of period .
Furthermore, if N is divisible by 4, according to Theorem 6.2, x N 2 +1 = −x 1 and x N 2 +2 = −x 2 . This implies that for all t > 0, x 1 (t + τ ) = −x 1 (t) = x N 2 +1 (t) and x 2 (t + τ ) = −x 2 (t) = x N 2 +2 (t), so the agents x 1 and x N 2 +1 share an orbit, as well as all pairs of agents x i and x N 2 +i for i ∈ {1, . . . , N 2 }. As a consequence, the kinetic energy is periodic, of period τ E = τ = π/( √ ab sin( 2π N )). If N is divisible by 4, every half period, the system is rotated by an angle π, so the kinetic energy is periodic with period τ E = τ /2.

Remark 7.
Notice that the agents sharing orbits do not interact with one another, as shown in Figure 16.
An example of such a choreography is given in Figure 17.
Remark 8. As a slight generalization, we provide numerical simulations illustrating a similar behavior, but with slightly different conditions: the periodic evolution of 9 agents on three distinct orbits shared three by three, see figures 18 and 19.
6.4. Helical trajectories. In sections 6.2 and 6.3, we provided conditions for the trajectories of the system to be periodic. Here, we explore further the notion of periodicity by studying systems with drift, displaying helical trajectories but periodic kinetic energy. Definition 6.6. Let (x i ) i=1...N be a solution of (42). We call the corresponding trajectories helical trajectories if there exists v ∈ R 2 and τ ∈ R * such that for all i ∈ {1, ..., N }, for all t > 0, x i (t + τ ) = x i (t) + τ v.
x 1 x 2 x 3 x 4 x 5 x 6 x 7 x 8 x 9 x 10 x 1 x 2 x 3 x 4 x 5 x 6 x 7 x 8 Figure 16. Left: Directed graph corresponding to the matrix A given in (44). Full arrows represent positive coefficients (a ij > 0) while dashed ones represent negative coefficients (a ij < 0). Right: Weighted directed graph corresponding to the matrix A given in (46). Thin arrows represent the weighted edges |a ij | = a while bold ones represent the weight |a ij | = b. Nodes with the same color and symbol share orbits but are not directly connected in the graph.
Notice that this definition generalizes the notion of periodic trajectories recalled in Definition 6.1, which corresponds to the case v = 0. When v = 0, the system has a drift term, meaning that the relative positions between agents remain periodic but their absolute positions evolve in space. Theorem 6.7. Sufficient conditions for helical trajectories. Let N = 4. Let (a, b, c, d) ∈ (R + ) 4 such that the interaction matrix reads Then the system exhibits helical trajectories. Kinetic Energy Figure 17. Left: Periodic trajectories of 8 agents sharing orbits two by two, in the situation of Theorem 6.5. Matrix A from (46) was constructed with (a, b) = (1, 3). The initial positions x 1 (0) and x 2 (0) were randomly generated and the other 6 were obtained by rotation. The period is τ = 2π/ √ 6. Right: Corresponding kinetic energy, of period τ /2.
Proof. First notice that the first and second components x i1 and x i2 of the i-th agent's position are decoupled, so that the system in matrix form readṡ (49) Hence the projections of x on the first and second axes solve the same differential equation. The matrixÃ has three distinct eigenvalues: There is one eigenvector associated with λ 1 : Let v 2 denote the eigenvector associated with λ 2 and let v R 2 and v I 2 denote respectively its real and imaginary components, i.e. v 2 := v R 2 + iv I 2 . Then the solution of System (49) can be written as: x j (t) = C j 1 v 1 + C j 2 (v 1 t + ν) + C j 3 v R 2 cos(λ 2 t) − v I 2 sin(λ 2 t) + C j 4 v R 2 sin(λ 2 t) + v I 2 cos(λ 2 t) where (C 1 , C 2 , C 3 , C 4 ) ∈ R 4 are constants depending on the initial conditions. Let τ = 2π λ2 . Then for all t > 0, for all i ∈ {1, . . . , 4}, for all j ∈ {1, 2}, x ij (t + τ ) = x ij (t) + C j 2 τ . This can be rewritten as: for all i ∈ {1, . . . , 4}, for all t > 0, x i (t + τ ) = x i (t) +  Theorem 6.8. A system with helical trajectories has periodic kinetic energy.

7.
Appendix. An example using Approach A which shows unexpected behavior in the first example (A.1), as well as the interactions matrix and initial positions used for simulations shown in Figure 7 and 9.
Example 7.1. 15 agents with a general interaction matrix A. We use the same notion of a general interactions matrix as used in [3]. The interactions matrix A is composed of integers a ij that are uniformly chosen between -5 and 5 inclusive. Ψ ≡ Id. Generally, a system with this kind of interaction matrix will exhibit simple oscillating kinetic energy, as in [3]. We show a rare simulation using this general interactions matrix A in the appendix (Figures 22 and 23).