Phase transitions and macroscopic limits in a BGK model of body-attitude coordination

TL;DR

This paper studies phase transitions in a kinetic model of body-attitude coordination, revealing how self-organization depends on local density and deriving macroscopic limits through stability and equilibrium analysis.

Contribution

It introduces a novel analysis of phase transitions in a BGK model of body-attitude coordination, linking it to nematic polymer models and deriving associated macroscopic equations.

Findings

Self-organization depends on local density.

Complete equilibrium description via nematic alignment models.

Convergence and stability proven for different regimes.

Abstract

In this article we investigate the phase transition phenomena that occur in a model of self-organisation through body-attitude coordination. Here, the body-attitude of an agent is modelled by a rotation matrix in $\mathbb{R}^3$ as in [Degond, Frouvelle, Merino-Aceituno, 2017]. The starting point of this study is a BGK equation modelling the evolution of the distribution function of the system at a kinetic level. The main novelty of this work is to show that in the spatially homogeneous case, self-organisation may appear or not depending on the local density of agents involved. We first exhibit a connection between body-orientation models and models of nematic alignment of polymers in higher dimensional space from which we deduce the complete description of the possible equilibria Then, thanks to a gradient-flow structure specific to this BGK model, we are able to prove the stability and the convergence towards the equilibria in the different regimes. We then derive the macroscopic models associated to the stable equilibria in the spirit of [Degond, Frouvelle, Merino-Aceituno, 2017] and [Degond, Frouvelle, Liu, 2015].

Figures

Equations

$$M_{\kappa \Omega }(\omega ):=\frac{e^{\kappa \Omega \cdot \omega }}{{\int }_{{\mathbb{S}}^{n-1}}{e}^{\kappa \Omega \cdot {\omega }^{\prime }}d{\omega }^{\prime }},$$

$$M_J(A):=\frac{e^{J\cdot A}}{{\int }_{S{O}_3(\mathbb{R})}{e}^{J\cdot A^{\prime }}d{A}^{\prime }},$$

$${\partial }_{t}f+(A{e}_1\cdot {\nabla }_x)f={\rho }_f{M}_{J_f}-f,$$

$${\rho }_f(t,x)={\int }_{S{O}_3(\mathbb{R})}f(t,x,A)dA\text{and}{J}_f={\int }_{S{O}_3(\mathbb{R})}f(t,x,A)AdA$$

$${\partial }_{t}f=\rho M_{J_f}-f,$$

$$\partial_{t}f^{\varepsilon}+(Ae_{1}\cdot\nabla_{x})f^{\varepsilon}=\frac{1}{\varepsilon}\Big{(}\rho_{f^{\varepsilon}}M_{J_{f^{\varepsilon}}}-f^{\varepsilon}\Big{)},$$

$${\rho }_f(t,x)={\int }_{S{O}_3(\mathbb{R})}f(t,x,A)dA\text{and}{J}_f(t,x)={\int }_{S{O}_3(\mathbb{R})}f(t,x,A)AdA.$$

$${\partial }_t{\rho }^{\epsilon }=\epsilon {\nabla }_x\cdot \left(\frac{\frac{1}{3}{\nabla }_x{\rho }^{\epsilon }}{1-\frac{{\rho }^{\epsilon }}{{\rho }_c}}\right),{\rho }_c=6.$$

$${\partial }_t\rho +{\nabla }_x\cdot (\rho c_1(\alpha (\rho ))\Lambda e_1)=0,$$

$$\rho ({\partial }_t\Lambda +{\tilde{c}}_2((\Lambda e_1)\cdot {\nabla }_x)\Lambda )+\widetilde{c_3}[(\Lambda e_1)\times {\nabla }_x\rho {]}_{\times }\Lambda$$

$$+c_4\rho [-{\mathbf{r}}_x(\Lambda )\times (\Lambda e_1)+{\delta }_x(\Lambda )\Lambda e_1{]}_{\times }\Lambda =0.$$

$$\mathrm{Orb}(M):=\{PMQ,P,Q\in S{O}_3(\mathbb{R})\}.$$

$$\forall i\in \{1,\dots ,N\},d{X}_t^{i,N}=(A_t^{i,N}{e}_1)dt,d{A}_t^{i,N}=0.$$

$$J^{i}(Z^{N}):=\frac{1}{N}\sum_{j=1}^{N}K\Big{(}|X^{i,N}-X^{j,N}|\Big{)}A^{j,N}\in\mathscr{M}_{3}(\mathbb{R}),$$

$${\partial }_{t}f+(A{e}_1\cdot {\nabla }_x)f={\rho }_f{M}_{J_{K\ast f}}-f$$

$$J_{K\ast f}(x):={\iint }_{{\mathbb{R}}^3\times S{O}_3(\mathbb{R})}K(x-y)Af(y,A)dydA\in {\mathcal{M}}_3(\mathbb{R}).$$

$$d{X}_k=A_k(t)e_{1}dt,$$

$$\displaystyle{dA_{k}}=P_{T_{A_{k}}}\circ\left(\Big{(}\frac{1}{N}\sum_{i=1}^{N}K(|X_{i}-X_{k}|)A_{i}\Big{)}dt+2\sqrt{D}dB_{t}\right).$$

$$J_t^N=\frac{1}{N}\sum _{i=1}^N{A}_t^{i,N},$$

$${\partial }_{t}f=M_{J_f}-f.$$

$${\partial }_{t}f+(A{e}_1\cdot {\nabla }_x)f={\rho }_f{M}_{J_f}-f,$$

$${\rho }_f(t,x)={\int }_{S{O}_3(\mathbb{R})}f(t,x,A)dA\text{and}{J}_f={\int }_{S{O}_3(\mathbb{R})}f(t,x,A)AdA.$$

$${\partial }_{t}f+(A{e}_1\cdot {\nabla }_x)f={\nabla }_A\cdot \left[M_{J_f}{\nabla }_A\left(\frac{f}{M_{J_f}}\right)\right],$$

$$\mathcal{F}[f]:={\int }_{S{O}_3(\mathbb{R})}f\log f-\frac{1}{2}|J_f{|}^2.$$

$$\frac{d}{dt}\mathcal{F}[f]=-\mathcal{D}[f]\le 0,$$

$$\mathcal{D}[f]={\int }_{S{O}_3(\mathbb{R})}(f-\rho M_{J_f})(\log f-\log (\rho M_{J_f}))\ge 0.$$

$${\partial }_{t}f+(A{e}_1\cdot {\nabla }_x)f={\nabla }_A\cdot \left[M_{{\Lambda }_f}{\nabla }_A\left(\frac{f}{M_{{\Lambda }_f}}\right)\right].$$

$$J_f={\int }_{{\mathbb{S}}^{n-1}}f(t,x,v)vdv,$$

$${\Omega }_f=\frac{{\int }_{{\mathbb{S}}^{n-1}}f(t,x,v)vdv}{\left|{\int }_{{\mathbb{S}}^{n-1}}f(t,x,v)vdv\right|}\in {\mathbb{S}}^{n-1},$$

$${\partial }_{t}f=\rho M_{J_f}-f,$$

$$\forall t\in {\mathbb{R}}_{+},{\rho }_f(t)={\rho }_{f_0},$$

$$f(t)=e^{-t}{f}_0+\rho {\int }_0^t{e}^{-(t-s)}{M}_{J_{f(s)}}ds,$$

Full text

Keywords: Collective motion; Vicsek model; generalized collision invariant; rotation group;

AMS Subject Classification: 34D05, 35Q92, 58J60, 82B26, 82C22, 92D50

1. 1 Introduction

2. 2 The BGK equation and other related models of self-organisation

3. 2.1 A review of the different IBM

4. 2.2 A review of the different kinetic equations

5. 3 Preliminaries: structure and calculus in $SO_{n}(\mathbb{R})$

6. 3.1 Structure and Haar measure on $SO_{n}(\mathbb{R})$

7. 3.2 Volume forms in $SO_{3}(\mathbb{R})$

8. 3.3 Singular Value Decomposition (SVD)

9. 4 Equilibria of the BGK operator

10. 4.1 Characterisation of the equilibria and compatibility equations

11. 4.2 Proof of Proposition 4.2

12. 4.3 Determination of the equilibria for each density $\rho$

13. 5 Convergence to equilibria

14. 5.1 A gradient-flow structure in $\mathbb{R}^{3}$

15. 5.1.1 Reduction to a nonlinear ODE in $\mathbb{R}^{3}$ and equilibria

16. 5.1.2 A gradient-flow structure

17. 5.2 Stability of the equilibria and conclusion

18. 5.3 Final remark: characterisation of the stable manifold of the unstable equilibria and critical cases

19. 6 Macroscopic limit for the stable equilibria

20. 6.1 Diffusion model in a disordered region

21. 6.2 Self-organised hydrodynamics in an ordered region

22. 6.2.1 Generalised collision invariants

23. 6.2.2 Hydrodynamic limit (formal proof of Theorem 9)

24. 7 Conclusion

25. A Quaternions and rotations

26. B More about $SO_{3}(\mathbb{R})$ and $SO_{n}(\mathbb{R})$

27. Acknowledgments

1 Introduction

The model studied in the present work is a new elaboration of the work initiated in [14] to model collective behaviour of agents described by their position and body-attitude. New results about emergence of phenomena of body-attitude coordination are presented in the context of a Bhatnagar-Gross-Krook (BGK) model. Such models can be applied to many biological systems such as flocking birds [37], fish school [36, 35] or sperm motion [20]. These systems are constituted by a large number of self-propelled agents which move at a constant speed and try to imitate their neighbours by moving in the same direction and trying to coordinate their body attitude. The agents are modelled by a moving frame in dimension 3, i.e. three orthogonal axes, one of which gives the direction of the motion and the two others the body orientation. In this work, as in [14], the body attitude is modelled by a rotation matrix in dimension 3, i.e. an element of the special orthogonal group $SO_{3}(\mathbb{R})$. In [17] agents are modelled by quaternions.

Collective behaviour in many-agent systems has been a thoroughly studied subject in the mathematical literature, from both theoretical and applied points of view. Among the models which have received the most attention, one can cite the Cucker-Smale model [10, 31, 47], attractive-repulsive models [7] or the Vicsek model for self-propelled particles [54]. The present work belongs to the class of Vicsek-inspired models. Such models have two main distinctive features, first the assumption that the particles are self-propelled and secondly a geometrical constraint: in the original work of Vicsek, the velocities of the particles have constant norm and the dynamics therefore takes place on the sphere $\mathbb{S}^{n-1}$ in dimension $n$ ($n=2$ in [54]). Here the dynamics takes place on the Riemannian manifold $SO_{3}(\mathbb{R})$.

The tools used to study models of collective behaviour are generally borrowed from the mathematical kinetic theory of gases which gives a mathematical framework to study many-particle systems. At a microscopic scale, the motion of each particle is detailed (Individual Based Model, IBM) through Ordinary Differential Equations (ODE) coming from Newton’s laws or through stochastic processes. When the number of particles is large, the whole system is described at a mesoscopic scale by a kinetic partial differential equation such as the Boltzmann, Fokker-Planck or BGK equation. Finally, large-scale dynamics is described by macroscopic equations (Euler, Navier-Stokes…). A review of the main results of kinetic theory of gases can be found in [11]. In particular the BGK equation (for Bhatnagar-Gross-Krook) was introduced in [3] as a substitute for the Boltzmann equation in the context of gas dynamics. The BGK operator is a relaxation operator towards a Maxwellian having the same moments as the distribution function of the system. Its mathematical properties and relevance in the mathematical kinetic theory of gases have been studied in particular in [49] and [51]. The BGK operator has been used in a model of collective dynamics of self-propelled particles in [22]. However, together with [16], it is the first time that it is rigorously studied in a body-attitude coordination model.

The main mathematical challenge in classical kinetic theory is the rigorous derivation of the kinetic equations from the IBM and of the macroscopic models from the kinetic equations. These questions are at the core of Hilbert’s sixth problem and have received much attention in the last decades. Many different techniques have been developed to derive kinetic equations from hard-sphere gases (Boltzmann-Grad limit [44, 25]), from systems of interacting particles (mean-field limit and propagation of chaos [40, 34, 53]) or from stochastic processes (and in particular jump processes [46, 43]). Some of these techniques have been adapted to problems arising in the study of collective behaviour [9, 4, 5]. The passage from kinetic equations to macroscopic models generally depends on physical constraints and in particular on conservation laws (hydrodynamic limits, Hilbert and Chapman-Enskog methods, see [8, 11] for a review) and is still an active research field [30, 23, 6, 29]. In the context of self-propelled particles, due to the lack of conservation laws which normally hold in the classical kinetic theory of gases, specific tools are needed. In [18], a methodological breakthrough has been achieved by introducing the so called Generalised Collisional Invariants (GCI) to rigorously link kinetic and macroscopic equations in the context of collective behaviour of self-propelled particles. This technique is now rigorously justified [42] and has already been successfully applied to a wide range of problems [41, 56]. It will be the key here to derive the macroscopic model in Section 6. This will lead to a system of partial differential equations on the mean density and body attitude, referred as the Self-Organised Hydrodynamics for Body-attitude coordination (SOHB) in [14].

The aim of this work is to show the emergence of collective behaviour and self-organisation which give rise to macroscopic scale patterns such as clusters, travelling bands etc. These patterns emerge from the collective interactions and are not directly encoded in the behaviour of the individual particles as described by the IBM. The continuum version of the Vicsek model [18] named the Self-Organised Hydrodynamics (SOH) model is an exemple of a model able to describe such emergence of self-organised dynamics. The Vicsek model describes a system where agents try to imitate their neighbours by adapting their direction of motion to the average direction of their neighbours. It has been shown that, in a certain scaling and when the equilibrium of the system is reached, the directions of motion of the agents are not uniformly distributed but follow a von Mises distribution. For $\kappa\in\mathbb{R}_{+}$ and $\Omega\in\mathbb{S}^{n-1}$ the von Mises distribution of parameters $\kappa$ and $\Omega$ is the Probability Density Function (PDF) on $\mathbb{S}^{n-1}$ defined by:

[TABLE]

where the dot product is the usual dot product in $\mathbb{R}^{n}$. This model [18] has been the starting point of many other models of self-organised dynamics, including [14] for the body-attitude coordination. In this context, we define the von Mises distribution of parameter $J\in\mathscr{M}_{3}(\mathbb{R})$ (a $3\times 3$ real matrix) as the following PDF on $SO_{3}(\mathbb{R})$:

[TABLE]

where the dot product and the measure on $SO_{3}(\mathbb{R})$ come from the Riemannian structure of $SO_{3}(\mathbb{R})$ detailed in Section 3.

In the present work, we focus on phase transition phenomena between non-organised and organised dynamics (collective motion). We will prove that the spatial density of agents is the key parameter which encodes the main features of phase transitions: in low density regions, no self-organised dynamics appears but when the density crosses a critical value, self-organised dynamics, given by a von Mises distribution for the body-attitude, becomes a stable equilibria of the system. This phase transition in the dynamics is purely an emergent phenomena, in the sense that at the macroscopic scale, different equations are required to describe the dynamics for different values of the density of agents, whereas for the IBM and at a mesoscopic level, the dynamics is described by one unique (system of) equation(s).

The starting point of this study is the BGK equation

[TABLE]

where $f(t,x,A)$ is a probability measure which gives the distribution of agents at position $x\in\mathbb{R}^{3}$ with body-orientation $A\in SO_{3}(\mathbb{R})$ at time $t\in\mathbb{R}_{+}$ and where:

[TABLE]

are the respective local density and flux. The measure on $SO_{3}(\mathbb{R})$ is the normalised Haar measure, the main properties of which are summarised in Section 3.

The left-hand side of the equation models the transport phenomenon: an agent with body orientation $A\in SO_{3}(\mathbb{R})$ moves in the direction $Ae_{1}$ where $e_{1}$ is the first vector of the canonical basis of $\mathbb{R}^{3}$. The right-hand side of the equation is the BGK operator which models the interactions between the agents: here we assume that $f$ relaxes towards a “moving equilibrium” which takes the form of a von Mises distribution. In particular, the von Mises distribution appears as the analog of the Maxwellian distribution of the classical gas dynamics. The flux $J_{f}$ plays the same role as the momentum density for gas dynamics or the average flux for the Vicsek model. The term $\rho_{f}M_{J_{f}}$ can therefore be seen as the analog of the “Maxwellian distribution with same moments as $f$” in the context of the BGK equation for gas dynamics.

The main results of this work are (informally) summarised in the two following theorems.

Theorem 1**.**

Let us consider the spatially homogeneous BGK equation:

[TABLE]

where $\rho\in\mathbb{R}_{+}$ is a given density of agents.

The equilibria $f^{\mathrm{eq}}$ of the spatially homogeneous BGK model are either the uniform equilibrium $f^{\mathrm{eq}}=\rho$ or of the form $f^{\mathrm{eq}}=\rho M_{\alpha\Lambda}$ or $f^{\mathrm{eq}}=\rho M_{\alpha\,p\otimes q}$ where $\Lambda\in SO_{3}(\mathbb{R})$ and $p,q\in\mathbb{S}^{2}$ and where $\alpha\in\mathbb{R}$ and $\rho$ are linked by a compatibility equation to be defined later (see Section 4 and equations (20) and (21)). 2. 2.

Depending on the density of agents $\rho\in\mathbb{R}_{+}$, the only stable equilibria are either the uniform equilibrium $f^{\mathrm{eq}}=\rho$ or the equilibria of the form $f^{\mathrm{eq}}=\rho M_{\alpha\Lambda}$ where $\Lambda\leavevmode\nobreak\ \in\leavevmode\nobreak\ SO_{3}(\mathbb{R})$ and where $\alpha\in\mathbb{R}_{+}$ is linked to $\rho$ by a compatibility equation to be defined later.

The first point of this theorem is detailed in Section 4 (see in particular Theorem 5 and Corollary 4.2). The second point is detailed in Section 5 (see in particular Theorem 7). We will then prove the following result.

Theorem 2** (Formal).**

Let us consider the rescaled spatially inhomogeneous problem

[TABLE]

where

[TABLE]

We assume that in a disordered region, $f^{\varepsilon}$ converges as $\varepsilon\to 0$ towards a density $\rho=\rho(t,x)$ uniform in the body-attitude variable. Then the density $\rho^{\varepsilon}\equiv\rho_{f^{\varepsilon}}$ satisfies at first order the following diffusion equation:

[TABLE] 2. 2.

We assume that in an ordered region, $f^{\varepsilon}$ converges as $\varepsilon\to 0$ towards an equilibrium of the form $\rho M_{\alpha\Lambda}$ with $\rho\in\mathbb{R}_{+}$, $\alpha\in\mathbb{R}_{+}$ and $\Lambda\in SO_{3}(\mathbb{R})$ defined above. Then the density $\rho=\rho(t,x)$ and mean body attitude $\Lambda=\Lambda(t,x)\in SO_{3}(\mathbb{R})$ satisfy the SOHB model given by the following system of partial differential equations:

[TABLE]

where $\alpha=\rho c_{1}(\alpha)$ and $\tilde{c}_{2}$, $\tilde{c_{3}}$, $c_{4}$ are functions of $\rho$ to be defined later and $\delta$ and $\mathbf{r}$ are the “divergence” and “rotational” operators defined in **[14]** (see Section 6)

This theorem is detailed in Section 6 (see in particular Proposition 6.1 and Theorem 9).

The phase transition problem has been completely treated in the space-homogeneous case for the Vicsek model in [13] but the geometrical structure inherent to body-orientation models requires specific tools and techniques. In particular, the rotation group $SO_{3}(\mathbb{R})$ is a compact Lie group, endowed with a Haar measure. The links between this topological structure and the Riemannian structure (detailed in Section 3 and Appendix B) will be the key to reduce the problem to a form that shares structural properties with the models of nematic alignment of polymers, studied in a completely different context to model liquid crystals [32, 55, 57, 2, 1]. These two worlds will be formally linked through the isomorphism between $SO_{3}(\mathbb{R})$ and the group of unit quaternions detailed in Section 4.2 and Appendix A. It will lead to the first point of Theorem 1 (the complete description of the equilibria, Section 4). As in [55] we will see that there exist a class of equilibria which cannot be interpreted as equilibria around a mean-body orientation. These equilibria were not studied in [14, 17]. A key point of the proof will be the reduction to a problem for diagonal matrices which will be a consequence of the left and right invariance of the Haar measure together with an adapted version of the Singular Value Decomposition of a matrix (Definition 3.2).

The stability of the different equilibria, are studied in Section 5.2. We will show that our model has an underlying gradient-flow structure which will allow us to determine the asymptotic behaviour of the system after a reduction to an ODE in $\mathbb{R}^{3}$. This is a specificity of the BGK model which doesn’t hold for the other models of body-attitude coordination [14, 17] and allows us to use different and simpler techniques. In particular, we will prove that the equilibria which cannot be interpreted as equilibria around a mean body-orientation are always unstable, which tends to justify the analysis carried out in [14] for a model where only equilibria around a mean body-orientation were considered.

Finally, the SOHB model (Section 6) will be obtained as in [14] by using the GCI. However, compared to [14], additional terms appear which require a specific treatment and in particular the coefficient $\tilde{c_{3}}$ that appears in Theorem 2 is different from the one that appears in [14]. The SOHB model (42) raises many questions, most of which are still open, and its mathematical and numerical analyses are still in progress. In particular, the hyperbolicity of the model is currently under study [15] and has been shown when $\tilde{c_{3}}$ is constant.

The organisation of the work is the following: in Section 2 we will give a review of the existing models at a microscopic and mesoscopic scales and motivate the study of the BGK equation among them. In Section 3, we gather the main technical results we will constantly use throughout this work. In Section 4, we will describe, depending on the density, all the possible equilibria of the system. We will use the tools developed to mathematically study the alignment of polymers [55, 57]. In Section 5 we will describe the asymptotic behaviour of the system and in particular which equilibria are attained, leading to a self-organised dynamics or not. This will be based on a specific underlying gradient-flow structure of the BGK equation. Finally in Section 6 we will write the macroscopic models for the stable equilibria.

Notations. For the convenience of the reader, we collect here the main notations we will use in the following.

1. $\bullet$

$\mathscr{M}_{n}(\mathbb{R})$ is the set of $n\times n$ real matrices. 2. $\bullet$

$\mathscr{D}_{n}(\mathbb{R})\subset\mathscr{M}_{n}(\mathbb{R})$ is the subspace of $n\times n$ diagonal real matrices. 3. $\bullet$

$\operatorname{Tr}(M)$ denotes the trace of the matrix $M\in\mathscr{M}_{n}(\mathbb{R})$ and $M^{T}$ its transpose. 4. $\bullet$

$I_{n}$ denotes the identity matrix in dimension $n$. 5. $\bullet$

$\operatorname{diag}:\mathbb{R}^{n}\to\mathscr{D}_{n}(\mathbb{R})$ is the vector space isomorphism such that for $(d_{1},\dots,d_{n})\in\mathbb{R}^{n}$, $D=\operatorname{diag}(d_{1},\dots,d_{n})$ is the diagonal matrix, the $(i,i)$-th coefficient of which is equal to $d_{i}$ for $i\in\{1,\ldots,n\}$. 6. $\bullet$

$\mathscr{S}_{n}(\mathbb{R})$ and $\mathscr{A}_{n}(\mathbb{R})$ denote respectively the sets of symmetric and skew-symmetric matrices of dimension $n$. 7. $\bullet$

$SO_{n}(\mathbb{R})$ is the special orthogonal group in dimension $n$, i.e. the group of matrices $P\in\mathscr{M}_{n}(\mathbb{R})$ such that $PP^{T}=I_{n}$ and $\det P>0$. 8. $\bullet$

$\mathbb{S}^{n}\subset\mathbb{R}^{n+1}$ is the sphere of dimension $n$. 9. $\bullet$

$\mathbb{H}$ is the group of unitary quaternions. 10. $\bullet$

$\langle\cdot\rangle_{g}$ denotes the mean for the probability density $g$ on $SO_{n}(\mathbb{R})$. We will simply write $\langle\cdot\rangle$ when $g$ is the uniform probability ($g\equiv 1$). 11. $\bullet$

$A$ will generically be a rotation matrix in $SO_{3}(\mathbb{R})$ and $a_{ij}$ its $(i,j)$ coefficient. 12. $\bullet$

$PD(M)$ is the orthogonal part of the polar decomposition of $M\in\mathscr{M}_{n}(\mathbb{R})$ when $\det M\neq 0$ : there exists a unique couple $(PD(M),S)\in\mathcal{O}_{n}(\mathbb{R})\times\SS_{n}(\mathbb{R})$ such that $M=PD(M)S$. The matrix $PD(M)$ is given by $PD(M)=M\left(\sqrt{M^{T}M}\right)^{-1}$. 13. $\bullet$

For a matrix $M\in\mathscr{M}_{3}(\mathbb{R})$, the orbit $\operatorname{Orb}(M)\subset\mathscr{M}_{3}(\mathbb{R})$ is defined by:

[TABLE] 14. $\bullet$

$\mathbb{R}_{+}:=[0,+\infty)$, $\mathbb{R}_{+}^{*}:=(0,+\infty)$

2 The BGK equation and other related models of self-organisation

In this section we give a review of the different existing models of collective dynamics at both microscopic and mesoscopic levels and emphasise the singularity of the BGK model among them.

2.1 A review of the different IBM

The rigorous proofs of the two following theorems (Theorems 3 and 4) can be found in [21] in a more general framework.

At a microscopic level, we fix a reference frame given by the canonical basis $(e_{1},e_{2},e_{3})$ of $\mathbb{R}^{3}$. The agents are described by their position $X\in\mathbb{R}^{3}$ and their body-attitude $A\in SO_{3}(\mathbb{R})$ which can be seen as a moving frame. We assume that an agent with body attitude $A\in SO_{3}(\mathbb{R})$ moves at a constant speed in the direction of the first vector of $A$ : the instantaneous velocity of the agent is $Ae_{1}$.

In the following we consider an increasing sequence of jump times $(T_{n})_{n}$ such that the increments between two jumps are independent and follow an exponential law of parameter $N\in\mathbb{N}^{*}$ (their expectation is $1/N$). The $N$ agents are described at time $t\in\mathbb{R}_{+}$ by their positions and body-attitudes Z^{N}_{t}=\big{\{}(X^{i,N}_{t},A^{i,N}_{t})\big{\}}_{i\in\{1,\ldots,N\}}\in\big{(}\mathbb{R}^{3}\times SO_{3}(\mathbb{R})\big{)}^{N}. The interactions between the agents can be modelled by the following Piecewise Deterministic Markov Process (PDMP) which has already been described heuristically in [16, 22]:

Between two jump times $({T}_{n},{T}_{n+1})$, the systems evolves in a deterministic way:

[TABLE] 2. 2.

At time ${T}_{n+1}$, a particle $i\in\{1,\ldots,N\}$ is chosen uniformly among the $N$ particles. At time ${T}_{n+1}^{+}$, the new body-orientation of particle $i$ is sampled from the PDF M_{J^{i}\big{(}Z_{{T}_{n+1}^{-}}^{N}\big{)}} where for $Z^{N}\in(\mathbb{R}^{3}\times SO_{3}(\mathbb{R}))^{N}$ we define the flux:

[TABLE]

and where $K$ is a smooth observation kernel.

The following theorem describes the limiting behaviour of the laws of the particles when $N\to+\infty$ under the assumption that the empirical measure of the $N$ processes converges (weakly) towards a smooth function $f$ (propagation of chaos property). The equation on $f$ can be derived as in [16, Section 4.2].

Theorem 3**.**

Let $f_{0}$ be a probability measure on the space $\mathbb{R}^{3}\times SO_{3}(\mathbb{R})$ and let $Z_{0}^{N}\in(\mathbb{R}^{3}\times SO_{3}(\mathbb{R}))^{N}$ be an initial state given by $N$ independent random variables, identically distributed with law $f_{0}$. Then for any $t\in\mathbb{R}_{+}$, the law $f^{N}_{t}$ of any of one of the processes $(X^{i,N}_{t},A^{i,N}_{t})_{t}$ at time $t$ converges weakly towards the solution $f_{t}$ of the following BGK equation with initial condition $f_{0}$:

[TABLE]

where $J_{K*f}$ is a matrix-valued function of the space variable $x\in\mathbb{R}^{3}$ defined by

[TABLE]

In the previous works, the IBM were typically given as in [14] by a system of stochastic differential equations such as the following:

[TABLE]

where $P_{T_{A_{k}}}$ denotes the projection on the tangent space of $SO_{3}(\mathbb{R})$ at $A_{k}\in SO_{3}(\mathbb{R})$ (see Section 3). In this case, the resulting equation when $N\to+\infty$ is a non-linear Fokker-Planck equation (see [5] for a rigorous proof in the Vicsek case).

In the spatially homogeneous case, we can take the observation kernel $K$ to be constantly equal to 1 to prove the mean-field limit. The agents are described at time $t\in\mathbb{R}_{+}$ only by their body-attitudes \big{\{}A^{i,N}\big{\}}_{i\in\{1,\ldots,N\}}\in SO_{3}(\mathbb{R})^{N} and they follow the following jump process: at each jump time $T_{n}$, compute the flux

[TABLE]

choose a particle $i\in\{1,\ldots N\}$ uniformly among the $N$ particles and draw the new body-orientation $A^{i,N}_{T_{n}^{+}}$ after the jump according to the law given by the PDF $M_{J^{N}_{T_{n}^{-}}}$. The following theorem describes analogously the limiting behaviour of the laws of the particles as $N\to+\infty$.

Theorem 4**.**

Let \big{\{}A^{i,N}_{0}\big{\}}_{i\in\{1,\ldots,N\}}\in SO_{3}(\mathbb{R})^{N} be an initial state given by $N$ independent random variables, identically distributed according to a law $f_{0}$ on $SO_{3}(\mathbb{R})$. Then for any $t\in\mathbb{R}_{+}$, the law $f^{N}_{t}$ of any of one of the processes $(A^{i,N}_{t})_{t}$ at time $t$ converges weakly towards the solution $f_{t}$ of the following spatially homogeneous BGK equation with initial condition $f_{0}$:

[TABLE]

2.2 A review of the different kinetic equations

The model studied in the present article belongs to a class of models, the study of which has been initiated in [18] as a continuum version of the Vicsek model [54]. These models can be classified in two types. First, in the Vicsek-type models, the agents are described by their orientation defined as a unit vector in $\mathbb{S}^{n-1}$. In the second type of models, we take into account their body-orientation, defined as a rotation matrix in $SO_{3}(\mathbb{R})$. Our study enters into this second framework.

The kinetic version of the Vicsek-type or Body-Orientation-type models is given either by a Fokker-Planck equation or by a BGK equation. In this work we will focus on the BGK equation

[TABLE]

where

[TABLE]

The Fokker-Planck version of our model corresponds to:

[TABLE]

where $\nabla_{A}$ and $\nabla_{A}\cdot$ are respectively the gradient and the divergence in $SO_{3}(\mathbb{R})$ for the Riemannian structure detailed in Section 3. Apart from the fact that the underlying interaction process [22, 16] which leads to the BGK model is different from the one that leads to the Fokker-Planck model, the BGK model is structurally different and can be treated independently by using specific and simpler mathematical techniques presented in the next sections. Nevertheless, the BGK and Fokker-Planck models share important properties. For instance, the following functional is a free-energy for both the spatially homogeneous BGK equation and the spatially homogeneous Fokker-Planck equation (though with a different dissipation term):

[TABLE]

It satisfies in both cases:

[TABLE]

where $\mathcal{D}[f]$ is the dissipation term which is equal for the BGK model to:

[TABLE]

In the context of the Vicsek model, this free energy was the key to study the phase transition phenomena [13] and we believe that the same kind of study can be made in the body-attitude coordination dynamics modelled by a Fokker-Planck equation (5). Moreover, in the Fokker-Planck case this dissipation inequality implies a gradient flow structure in the Wasserstein-2 distance which has been studied (in the Vicsek case) in [24]. However, the BGK model has another underlying gradient-flow dynamics (studied in Section 5) on which the present study will be based, and we will therefore not use this free-energy in the present work.

Both models (BGK and Fokker-Planck) have a normalised and a non-normalised version. The model (4) will be referred as the non-normalised BGK model. A normalised model is a model where the flux $J_{f}$ is replaced by the orthogonal part of its polar decomposition $\Lambda_{f}:=PD(J_{f})$ as defined in the introduction and under the assumption that $\det J_{f}>0$. The normalised Fokker-Planck model is the model studied in [14]:

[TABLE]

This terminology comes from the continuum version of the Vicsek model [18] where either the total flux

[TABLE]

or its normalisation

[TABLE]

is considered. A mathematical analysis of the normalised Vicsek model can be found in [24, 26]. The importance of this distinction in the context of phase transitions has been shown in [13] and [12]: phase transitions appear only in non-normalised models.

The following chart (Figure 1) shows the different models and gives references where they are studied (when such references exist).

Finally, in Sections 4 and 5, we will focus on the spatially homogeneous version of the BGK model (4) given by:

[TABLE]

where the probability distribution $f(t,A)$ only depends on the body-orientation variable and time. In the spatially homogeneous case, the local density of agents previously denoted by $\rho_{f}$ does not depend on $f$ in the sense that an initial density $\rho_{f_{0}}\in\mathbb{R}_{+}$ associated to the initial distribution $f_{0}$ is preserved by the dynamics:

[TABLE]

as it can be seen by integrating the equation over $SO_{3}(\mathbb{R})$. We therefore take $\rho\in\mathbb{R}_{+}$ as a fixed parameter of the problem. Note also that the well-posedness of (7) directly follows from Duhamel’s formula:

[TABLE]

since $J_{f}$ is given as the solution of the following differential equation on $\mathscr{M}_{3}(\mathbb{R})$:

[TABLE]

as it can be seen by multiplying (7) by $A$ and integrating over $SO_{3}(\mathbb{R})$. Note that it contrasts with the Fokker-Planck case where even the well-posedness of the spatially-homogeneous equation would require further investigations. This will be part of future work.

3 Preliminaries: structure and calculus in $SO_{n}(\mathbb{R})$

This paragraph collects the main properties of the Riemannian manifold $SO_{n}(\mathbb{R})$ and other technical results. In this paragraph $n\geq 3$ denotes the dimension, we will mainly consider the case $n=3$ in the next sections.

3.1 Structure and Haar measure on $SO_{n}(\mathbb{R})$

Lemma 3.1**.**

The following is an inner product on $\mathscr{M}_{n}(\mathbb{R})$ :

[TABLE]

and the following properties hold:

1. $\bullet$

Endowed with this metric, $SO_{n}(\mathbb{R})$ is a topological group and a Riemannian manifold. 2. $\bullet$

The sets $\mathscr{S}_{n}(\mathbb{R})$ and $\mathscr{A}_{n}(\mathbb{R})$ of symmetric and skew-symmetric matrices are orthogonal and $\mathscr{M}_{n}(\mathbb{R})=\mathscr{S}_{n}(\mathbb{R})\oplus\mathscr{A}_{n}(\mathbb{R})$. 3. $\bullet$

For $A\in SO_{n}(\mathbb{R})$, the tangent space to $SO_{n}(\mathbb{R})$ at $A$ is denoted by $T_{A}$ and

[TABLE]

The norm on $\mathscr{M}_{n}(\mathbb{R})$ associated to the inner product (8) will be denoted by $\|\cdot\|$.

The general theory of locally compact topological groups ensures the existence of a Haar measure $\mu$ on $SO_{n}(\mathbb{R})$ which satisfies for all $P\in SO_{n}(\mathbb{R})$ and all Borel set $\mathcal{E}$ of the Borel $\sigma$-algebra of $SO_{n}(\mathbb{R})$:

[TABLE]

where $P\mathcal{E}=\{PA,\,\,A\in\mathcal{E}\}$ and $\mathcal{E}P=\{AP,\,\,A\in\mathcal{E}\}$. We will assume that $\mu$ is the unique Haar measure which is a probability measure and simply write

[TABLE]

As a consequence if $P\in SO_{n}(\mathbb{R})$, $A\mapsto PA$ and $A\mapsto AP$ are two changes of variable with unit Jacobian. We will constantly use the following changes of variable :

Definition 3.1** (Useful changes of variable).**

Let us define the following matrices:

1. $\bullet$

For $i\neq j\in\{1,\ldots,n\}$, $D^{ij}\in SO_{n}(\mathbb{R})$ is the diagonal matrix such that all its coefficients are equal to $1$ except at positions $i$ and $j$ where they are equal to $-1$. 2. $\bullet$

For $i\neq j\in\{1,\ldots,n\}$, $P^{ij}\in SO_{n}(\mathbb{R})$ is the matrix such that $P^{ij}_{ii}=P^{ij}_{jj}=0$, $P^{ij}_{kk}=1$ for $k\neq i,j$, $P^{ij}_{ij}=1$ and $P^{ij}_{ji}=-1$. The other coefficients are equal to 0.

Then we define the following changes of variable with unit Jacobian:

1. $\bullet$

$A^{\prime}=D^{ij}A$ multiplies the rows $i$ and $j$ by $-1$. Everything else remains unchanged. 2. $\bullet$

$A^{\prime}=AD^{ij}$ multiplies the columns $i$ and $j$ by $-1$. Everything else remains unchanged. 3. $\bullet$

$A^{\prime}=D^{ij}AD^{ij}$ multiplies the elements $(k,i)$, $(k,j)$ and $(i,k)$, $(j,k)$ by $-1$ for $k\neq i,j$. Everything else remains unchanged 4. $\bullet$

$A^{\prime}=P^{ij}A$ multiplies row $i$ by $-1$ and permutes the rows $i$ and $j$. 5. $\bullet$

$A^{\prime}=P^{ij}A(P^{ij})^{T}$ exchanges the diagonal coefficients $(i,i)$ and $(j,j)$ (and involves other changes).

The two following lemmas are important applications of these results.

Lemma 3.2**.**

Let $D\in\mathscr{M}_{n}(\mathbb{R})$ be a diagonal matrix and $M_{D}$ the von Mises distribution with parameter $D$, then

[TABLE]

is diagonal.

Proof.

Let $k\neq\ell$ and $m\neq k,\ell$. The change of variable $A\mapsto D^{km}AD^{km}$ gives:

[TABLE]

where we have used that $D^{k,m}DD^{k,m}=D$. ∎

Lemma 3.3**.**

For any $n\geq 3$ and any $J\in\mathscr{M}_{n}(\mathbb{R})$,

[TABLE]

Lemma 3.4**.**

Let $n\geq 3$, $n\neq 4$. Let $g:SO_{n}(\mathbb{R})\to\mathbb{R}$ such that for all $A,P\in SO_{n}(\mathbb{R})$, $g(A)=g(A^{T})\leavevmode\nobreak\ =\leavevmode\nobreak\ g(PAP^{T})$. For all $J\in\mathscr{M}_{n}(\mathbb{R})$ we have:

[TABLE]

for given $a,b,c\in\mathbb{R}$ depending on $g$ and on the dimension, the expressions of which can be found in the proof.

The proof of these lemmas and other technical results about $SO_{3}(\mathbb{R})$ and $SO_{n}(\mathbb{R})$ are postponed to Appendix B.

3.2 Volume forms in $SO_{3}(\mathbb{R})$

When an explicit calculation will be needed, we will use one of the two following parametrisations of $SO_{3}(\mathbb{R})$ which give two explicit expressions of the normalised Haar measure in dimension 3.

1. $\bullet$

To a matrix $A\in SO_{3}(\mathbb{R})$ there is an associated angle $\theta\in[0,\pi]$ and a vector $\mathbf{n}\in\mathbb{S}^{2}$ such that $A$ is the rotation of angle $\theta$ around the axis $\mathbf{n}$. Rodrigues’ formula gives a representation of $A$ knowing $\theta$ and $\mathbf{n}=(n_{1},n_{2},n_{3})$ :

[TABLE]

where

[TABLE]

and we have:

[TABLE]

If $f(A(\theta,\mathbf{n}))=\bar{f}(\theta,\mathbf{n})$ the volume form of $SO_{3}(\mathbb{R})$ is given by:

[TABLE]

With the usual parametrisation of the sphere $\mathbb{S}^{2}$ we can take $\mathbf{n}=(n_{1},n_{2},n_{3})^{T}$ with

[TABLE]

where $\psi\in[0,\pi]$ and $\varphi\in[0,2\pi]$. The volume form for the sphere is given by:

[TABLE] 2. $\bullet$

We have the following one to one map :

[TABLE]

where

[TABLE]

and for $p=(\sin\phi_{1}\,\sin\phi_{2},\cos\phi_{1}\,\sin\phi_{2},\cos{\phi_{2}})^{T}$ in spherical coordinates $\phi_{1}\in[0,2\pi]$ and $\phi_{2}\in[0,\pi]$, we define:

[TABLE]

The matrix $A^{a}$ performs an arbitrary rotation of the first $2$ coordinates and the matrix $M(p)\in SO_{3}(\mathbb{R})$ maps the vector $e_{3}$ to $p\in\mathbb{S}^{2}$. A matrix $A\in SO_{3}(\mathbb{R})$ can thus be written as the product:

[TABLE]

where $\phi_{1},\theta\in[0,2\pi]$ and $\phi_{2}\in[0,\pi]$. With this parametrisation:

[TABLE]

and the volume form on the sphere is given by:

[TABLE]

This parametrisation can be extended in any dimension and comes from the Lie groups quotient:

[TABLE]

3.3 Singular Value Decomposition (SVD)

We recall the following classical result proved in [50, Section 1.9].

Proposition 3.1** (Singular Value Decomposition, SVD).**

Any square matrix $M\in\mathscr{M}_{n}(\mathbb{R})$ can be written:

[TABLE]

where $P,Q\in\mathcal{O}_{n}(\mathbb{R})$ and $D$ diagonal with nonnegative coefficients listed in decreasing order.

In order to use the properties of the Haar measure, we will need the matrices $P$ and $Q$ to belong to $SO_{3}(\mathbb{R})$ (not only $\mathcal{O}_{3}(\mathbb{R}))$ and we define therefore another decomposition, called the Special Singular Value Decomposition (SSVD) in the following.

Definition 3.2** (SSVD in $SO_{3}(\mathbb{R})$).**

Let $M\in\mathscr{M}_{3}(\mathbb{R})$. A Special Singular Value Decomposition (SSVD) of $M$ is a decomposition of the form

[TABLE]

where $P,Q\in SO_{3}(\mathbb{R})$ and $D=\operatorname{diag}(d_{1},d_{2},d_{3})$ with

[TABLE]

The existence of a SSVD follows from Proposition 3.1. Let us start from a SVD

[TABLE]

1. $\bullet$

If $\det M>0$, either $P^{\prime},Q^{\prime}\in SO_{3}(\mathbb{R})$ and the SVD is a SSVD or $P^{\prime},Q^{\prime}$ have both negative determinant and in this case we can take

[TABLE]

where $\tilde{D}=\operatorname{diag}(1,1,-1)$. 2. $\bullet$

If $\det M<0$, either $P^{\prime}\in SO_{3}(\mathbb{R})$ or $Q^{\prime}\in SO_{3}(\mathbb{R})$ (only one of them). Assume without loss of generality that $Q^{\prime}\in SO_{3}(\mathbb{R})$. Then we can take:

[TABLE] 3. $\bullet$

If $\det M=0$, then the last coefficient of $D^{\prime}$ is equal to 0 so $\tilde{D}D^{\prime}=D^{\prime}$ and $D^{\prime}\tilde{D}=D^{\prime}$. We can take $D=D^{\prime}$. If $P^{\prime}\notin SO_{3}(\mathbb{R})$ we can take $P=P^{\prime}\tilde{D}$ and if $Q^{\prime}\notin SO_{3}(\mathbb{R})$ we can take $Q=\tilde{D}Q^{\prime}$.

Remark 3.1**.**

As for the polar decomposition and the standard SVD, the matrix $D$ is always unique. However the matrices $P$ and $Q$ may not be unique.

The subset $\mathscr{D}\subset\mathscr{M}_{3}(\mathbb{R})$ of the diagonal matrices which are the diagonal part of a SSVD is the cone delimited by the image by the isomorphism $\operatorname{diag}$ of the three planes $\{d_{1}=d_{2}\}$, $\{d_{2}=d_{3}\}$ and $\{d_{2}=-d_{3}\}$ in $\mathbb{R}^{3}$ and depicted in Figure 3 :

[TABLE]

4 Equilibria of the BGK operator

In this section we determine the equilibria for the BGK operator:

[TABLE]

that is to say the distributions $f$ such that $Q_{BGK}(f)=0$. In Section 4.1 we characterise these equilibria (Theorem 5) and show that for them to exist, compatibility equations must be fulfilled. These compatibility equations depend on the density $\rho$. Therefore, for different values of the density $\rho$, there exists different equilibria. These will be determined in Section 4.2 by studying the compatibility equations. A full description of the equilibria of the BGK operator is finally given in Corollary 4.2.

4.1 Characterisation of the equilibria and compatibility equations

The main result of this section is Theorem 5 which gives all the equilibria of the BGK operator (13). Before stating and proving it we will need the following lemma which is the analog of lemma 4.4 in [14]. The proof of this lemma is an application of the results presented in Section 3.

Lemma 4.1** (Consistency relations).**

The following holds:

1. (i)

There exists a function $c_{1}=c_{1}(\alpha)$ defined for all $\alpha\in\mathbb{R}$ such that for all $\Lambda\in SO_{3}(\mathbb{R})$,

[TABLE]

The function $c_{1}$ can be explicitly written c_{1}(\alpha)=\frac{1}{3}\big{\{}(2\cos\theta+1)\big{\}}_{\alpha} where $\{\cdot\}_{\alpha}$ denotes the mean with respect to the probability density

[TABLE] 2. (ii)

Consider the set $\mathscr{B}\subset\mathscr{M}_{3}(\mathbb{R})$ defined by:

[TABLE]

There exists a function $c_{2}=c_{2}(\alpha)$ defined for all $\alpha\in\mathbb{R}$ such that for all $B\in\mathscr{B}$,

[TABLE]

The function $c_{2}$ can be explicitly written: $c_{2}(\alpha)=[\cos\phi]_{\alpha}$, where $[\cdot]_{\alpha}$ denotes the mean with respect to the probability density

[TABLE]

Remark 4.1**.**

The relevance of the set $\mathscr{B}$ will become apparent in Proposition 4.2.

Proof.

1. (i)

Using the left invariance of the Haar measure, it is enough to prove the result for $\Lambda=I_{3}$, since

[TABLE]

When $\Lambda=I_{3}$, Lemma 3.2 first ensures that $\langle A\rangle_{M_{\alpha I_{3}}}$ is diagonal, then the change of variable $A^{\prime}=P^{12}A(P^{12})^{T}$ (see Definition 3.1) shows that:

[TABLE]

Proceeding analogously with the other coefficients we have that $\langle A\rangle_{M_{\alpha I_{3}}}$ is proportional to $I_{3}$, i.e. there exists $c_{1}=c_{1}(\alpha)\in\mathbb{R}$ such that

[TABLE]

The parametrisation of $SO_{3}(\mathbb{R})$ using Rodrigues’ formula (9) then gives the explicit expression of $c_{1}$ by taking the trace in Equation (18) and using that for $A=A(\theta,\mathbf{n})$, $\operatorname{Tr}(A)=2\cos\theta+1$. 2. (ii)

As before, using the left and right invariance of the Haar measure it is enough to prove the result for $B=\operatorname{diag}(1,0,0)$. Now if $D=\operatorname{diag}(a,b,-b)$ for $a,b\in\mathbb{R}$, then the change of variable $A\mapsto P^{23}A(P^{23})^{T}$ followed by the change of variable $A\mapsto D^{23}A$ (see Definition 3.1) show that

[TABLE]

which proves with lemma 3.2 that $\langle A\rangle_{M_{D}}$ is diagonal of the form $\operatorname{diag}(\tilde{a},\tilde{b},-\tilde{b})$ for $\tilde{a},\tilde{b}\in\mathbb{R}$. Similarly, if $D=\operatorname{diag}(a,b,b)$ then $\langle A\rangle_{M_{D}}$ is of the form $\operatorname{diag}(\tilde{a},\tilde{b},\tilde{b})$. These two results prove that $\langle A\rangle_{M_{\alpha B}}$ is proportional to $B$, i.e. there exists $c_{2}=c_{2}(\alpha)\in\mathbb{R}$ such that (16) holds. The parametrisation of $SO_{3}(\mathbb{R})$ coming from the isomorphism (10) then gives the explicit expression of $c_{2}$ by taking $B=\operatorname{diag}(1,0,0)$ in Equation (16). First, using the change of variable $A\mapsto P^{13}A(P^{13})^{T}$ it holds that,

[TABLE]

where

[TABLE]

Then, using the parametrisation (11), it follows that:

[TABLE]

∎

Remark 4.2**.**

We could alternatively use one of the two parametrisations of $SO_{3}(\mathbb{R})$ given in Section 3.2 or the quaternion formulation to prove that $\langle A\rangle_{\alpha I_{3}}$ and $\langle A\rangle_{\alpha B}$ are proportional to $I_{3}$ and $B$. However, the proof that we have just presented here holds in any dimension (the value of the constants $c_{1}(\alpha)$ and $c_{2}(\alpha)$ depends on the dimension but not the form of the matrices) whereas the volume forms and the quaternion formulation strongly depend on the dimension $n=3$.

We can now state the main result of this section:

Theorem 5** (Equilibria for the homogeneous Body-Orientation BGK equation).**

Let $\rho\leavevmode\nobreak\ \in\leavevmode\nobreak\ \mathbb{R}_{+}$ be a given density. The equilibria of the spatially homogeneous BGK equation (7) are the distributions of the form $f=\rho M_{J}$ where $J\in\mathscr{M}_{3}(\mathbb{R})$ is a solution of the matrix compatibility equation:

[TABLE]

The solutions of the compatibility equation (19) are:

the matrix $J=0$, 2. 2.

the matrices of the form $J=\alpha\Lambda$ with $\Lambda\in SO_{3}(\mathbb{R})$ and where $\alpha\in\mathbb{R}$ satisfies the scalar compatibility equation

[TABLE] 3. 3.

the matrices of the form $J=\alpha B$ where $B\in\mathscr{B}$ and where $\alpha\in\mathbb{R}$ satisfies the scalar compatibility equation

[TABLE]

where the set $\mathscr{B}$ and the functions $c_{1}$ and $c_{2}$ are defined in Lemma 4.1.

Remark 4.3**.**

Notice that the existence of a non-zero solution for the scalar compatibility equations (20) and (21) is not guaranteed for all values of $\rho>0$ . The existence of non-zero solutions for these equations will be explored in Section 4.2. They will determine the existence of equilibria for Equation (7) for a given value of $\rho$ (Corollary 4.2).

Remark 4.4**.**

The fact that these matrices are solutions of the matrix compatibility equation (19) follows directly from the consistency relations (14) and (16) as it will be shown in the proof of Theorem 5. The main difficulty of the proof is therefore the necessary condition: we will prove that a solution of the matrix compatibility equation (19) is necessarily of one of the forms listed in Theorem 5.

The proof of this theorem will use the two following propositions. The first one and its corollary (Proposition 4.1 and Corollary 4.1) show that the compatibility equation (19) can be reduced to a compatibility equation on diagonal matrices (equation (22)). The second one (Proposition 4.2) provides a necessary condition for a diagonal matrix to be a solution of (22). The proof of Proposition 4.2 is deferred to the next section.

Proposition 4.1** (Orbital reduction).**

The following equivalence holds: $J\in\mathscr{M}_{3}(\mathbb{R})$ is a solution of the matrix compatibility equation (19) if and only if for all $J^{\prime}\in\operatorname{Orb}(J)$, $J^{\prime}$ is a solution of the matrix compatibility equation (19).

Proof.

This is a consequence of the left and right invariance of the Haar measure which ensures that for any $J\in\mathscr{M}_{3}(\mathbb{R})$ and any $P,Q\in SO_{3}(\mathbb{R})$ :

[TABLE]

∎

Since the diagonal part of the SSVD of a matrix $J$ is in the orbit of $J$, we obtain the following corollary:

Corollary 4.1** (Reduction to diagonal matrices).**

Let $J\in\mathscr{M}_{3}(\mathbb{R})$ with SSVD given by $J\leavevmode\nobreak\ =\leavevmode\nobreak\ PDQ$. The following equivalence holds: $J$ is a solution of (19) if and only if $D$ is a solution of (19).

We will therefore consider only the following problem in dimension 3: find all the diagonal matrices $D\in\mathscr{M}_{3}(\mathbb{R})$ such that

[TABLE]

where the set $\mathscr{D}$ is the subset of diagonal matrices which are the diagonal part of a SSVD and is defined by (12). Notice that Equation (22) is just Equation (19) restricted to the set $\mathscr{D}$.

Remark 4.5**.**

The diagonal part $D\in\mathscr{M}_{3}(\mathbb{R})$ of a SSVD of a matrix $J\in\mathscr{M}_{3}(\mathbb{R})$ is unique so the problems (19) and (22) are equivalent. Notice that there might be other diagonal matrices in $\operatorname{Orb}(J)$ (take for example $J$ diagonal which does not satisfy the conditions (12)). However the diagonal part of any SSVD of these matrices is $D$ : the diagonal part of the SSVD characterises the orbit of a matrix. In the following, we will find all the diagonal solutions of (19) (i.e. the solutions of (22) without the restriction $D\in\mathscr{D}$) and then only consider the ones which belong to $\mathscr{D}$. For instance we will see that there are solutions of (19) of the form $\operatorname{diag}(0,-\alpha,0)$ where $\alpha>0$. The diagonal part of their SSVD is $\operatorname{diag}(\alpha,0,0)$ and is a solution of (22).

Remark 4.6**.**

A diagonal solution $D$ of the matrix compatibility equation (19) verifies that $D/\rho$ belongs to the set:

[TABLE]

where $\mathcal{P}(SO_{3}(\mathbb{R}))$ is the set of probability measures on $SO_{3}(\mathbb{R})$. The set $\operatorname{diag}^{-1}(\Omega)\subset\mathbb{R}^{3}$ is exactly the tetrahedron $\mathscr{T}$ defined as the convex hull of the points $(\pm 1,\pm 1,\pm 1)$ with an even number of minuses (which we will call Horn’s tetrahedron). It is a consequence of Horn’s theorem [39, Theorem 8] which states that $\mathscr{T}$ is exactly the set of vectors which are the diagonal of an element of $SO_{3}(\mathbb{R})$. It ensures that if $f$ is a probability measure, we have by convexity of $\mathscr{T}$ :

[TABLE]

and therefore $\operatorname{diag}^{-1}(\Omega)\subset\mathscr{T}$. Conversely, taking the Dirac deltas $\delta_{I_{3}}$ and similarly for the other vertices of $\mathscr{T}$, we see that the four vertices of Horn’s tetrahedron belong to $\operatorname{diag}^{-1}(\Omega)$. Since $\Omega$ is convex, we conclude that $\mathscr{T}\subset\operatorname{diag}^{-1}(\Omega)$.

The diagonal solutions of the matrix compatibility equation (19) satisfy the following necessary condition.

Proposition 4.2**.**

The diagonal solutions of the compatibility equation (19) are necessarily of one of the following the types :

1. (a)

$D=0$. 2. (b)

$D=\alpha\operatorname{diag}(\pm 1,\pm 1,\pm 1)$* with an even number of minus signs and where $\alpha\in\mathbb{R}\setminus\{0\}$.*

If $\alpha\in(0,+\infty)$, the diagonal part of the SSVD of these diagonal matrices is equal to $D=\alpha I_{3}$.

If $\alpha\in(-\infty,0)$, the diagonal part of the SSVD of these diagonal matrices is equal to $D=\alpha\operatorname{diag}(-1,-1,1)=|\alpha|\operatorname{diag}(1,1,-1)$. 3. (c)

$D=\alpha\operatorname{diag}(\pm 1,0,0)$* and the matrices obtained by permutation of the diagonal coefficients and where $\alpha\in\mathbb{R}\setminus\{0\}$.*

The diagonal part of the SSVD of these diagonal matrices is equal to $D=\operatorname{diag}(|\alpha|,0,0)$.

Section 4.2 will be devoted to the proof of this proposition. We are now ready to prove Theorem 5.

Proof (of Theorem 5).

An equilibria of the BGK equation is of the form

[TABLE]

where

[TABLE]

It is straightforward to check that $J=0$ is a solution of (19). Now, let $D$ a matrix of one the form described in Proposition 4.2 with a parameter $\alpha\in\mathbb{R}$. For instance, for a matrix of type (c) like $D=\alpha\operatorname{diag}(0,-1,0)$, thanks to Lemma 4.1 we have:

[TABLE]

Similarly for the other diagonal matrices of type (c), we prove that they are solution of the matrix compatibility equation (19) if and only if their parameter $\alpha\in\mathbb{R}$ is solution of the scalar compatibility equation (21). Analogously one can check that the diagonal matrices of type (b) are solutions of the matrix compatibility equation (19) if and only if their parameters $\alpha\in\mathbb{R}$ are solutions of the scalar compatibility equation (20). This yields all the diagonal solutions of (19). Now, the solutions of (19) are exactly the matrices $J\in\operatorname{Orb}(D)$ where $D$ is a diagonal solution of (19) and the set $\operatorname{Orb}(D)\subset\mathscr{M}_{3}(\mathbb{R})$ is the orbit of $D$ defined in the introduction. We conclude by noticing that if $D$ is of type $(b)$ then $\operatorname{Orb}(D)=SO_{3}(\mathbb{R})$ and if $D$ is of type $(c)$ then $\operatorname{Orb}(D)=\mathscr{B}$. ∎

Remark 4.7**.**

When applied to diagonal matrices, the last part of Theorem 5 states that the diagonal solutions of (19) are necessarily of one of the types (a), (b) or (c) defined in Proposition 4.2 and that, it holds that

the matrix [math] is always a solution of (19), 2. 2.

a matrix of type (b) is a solution of (19) iff its parameter $\alpha\in\mathbb{R}\setminus\{0\}$ satisfies (20), 3. 3.

a matrix of type (c) is a solution of (19) iff its parameter $\alpha\in\mathbb{R}\setminus\{0\}$ satisfies (21).

4.2 Proof of Proposition 4.2

The proof of Proposition 4.2 is based on two results. The first one has been proved in [55, Section 4] to study the nematic alignment of polymers in higher dimensional spaces:

Theorem 6** ([55]).**

Let $n\geq 3$, $b\in\mathbb{R}_{+}$ and $\mathbf{s}=(s_{1},s_{2},\ldots,s_{n})\in\mathbb{R}^{n}$ a solution of the nonlinear system

[TABLE]

where the average is taken with respect to the PDF on the sphere $\mathbb{S}^{n-1}$ :

[TABLE]

where $Z$ is the normalisation constant which ensures that $g_{\mathbf{s},b}$ is a PDF on the sphere $\mathbb{S}^{n-1}$. Then $\operatorname{Card}\{s_{1},s_{2},\ldots,s_{n}\}\leq 2$.

The second tool that we will use to prove Proposition 4.2 is an isomorphism between $SO_{3}(\mathbb{R})$ and the space of unitary quaternions which transforms the compatibility equation (22) into the compatibility equation (23) studied in Theorem 6.

Proposition 4.3**.**

There is an isomorphism between the group $SO_{3}(\mathbb{R})$ and the quotient group $\mathbb{H}/\pm 1$, where $\mathbb{H}$ is the group of unit quaternions. Since $\mathbb{H}$ is homeomorphic to $\mathbb{S}^{3}$, there is an isomorphism $\Phi$ :

[TABLE]

Moreover $\Phi$ is an isometry in the sense that it maps the volume form of $\mathbb{S}^{3}/\pm 1$ (defined as the image measure of the usual measure on $\mathbb{S}^{3}$ by the projection on the quotient space) to the volume form on $SO_{3}(\mathbb{R})$: for all measurable function $f$ on $SO_{3}(\mathbb{R})$,

[TABLE] 2. 2.

There is a linear isomorphism between the vector space $\mathscr{M}_{3}(\mathbb{R})$ and the vector space $\SS_{4}^{0}(\mathbb{R})$ of trace free symmetric matrices of dimension 4:

[TABLE]

such that for all $J\in\mathscr{M}_{3}(\mathbb{R})$, and $q\in\mathbb{H}/\pm 1$,

[TABLE]

The first dot product is defined by Equation (8) and the second one is the usual dot product in $\mathbb{R}^{4}$. 3. 3.

For all $q\in\mathbb{H}/\pm 1$, it holds that \phi\big{(}\Phi(q)\big{)}=q\otimes q-\frac{1}{4}I_{4}. 4. 4.

The isomorphism $\phi$ preserves the diagonal structure: if $D=\operatorname{diag}(d_{1},d_{2},d_{3})$ then,

[TABLE]

and if $Q=\operatorname{diag}(s_{1},s_{2},s_{3},s_{4})$ with $s_{1}+s_{2}+s_{3}+s_{4}=0$, then

[TABLE]

The proof of this proposition can be found in appendix A. We are now ready to prove Proposition 4.2.

Proof (of Proposition 4.2).

Using the first and second points of Proposition 4.3, it holds that

[TABLE]

The compatibility equation (22) then becomes:

[TABLE]

Applying the isomorphism $\phi$ defined in Proposition 4.3 to this last equation, we obtain thanks to the third point of Proposition 4.3 :

[TABLE]

Using the fourth point of Proposition 4.3, we then obtain the following equivalent problem: find all the trace-free diagonal matrices $Q=\operatorname{diag}(s_{1},s_{2},s_{3},s_{4})$ of dimension 4 such that

[TABLE]

where $Z$ is a normalisation constant:

[TABLE]

Equivalently, defining for $i\in\{1,2,3,4\}$ :

[TABLE]

we want to solve the system of compatibility equations:

[TABLE]

where $\mathbf{s}^{\prime}=(s_{1}^{\prime},s_{2}^{\prime},s_{3}^{\prime},s_{4}^{\prime})$ and $g_{\mathbf{s}^{\prime},2\rho}$ is given by (24). Thanks to Theorem 6, we conclude that if $\mathbf{s}^{\prime}$ is a solution of (25), then the coefficients $s_{1}^{\prime},s_{2}^{\prime},s_{3}^{\prime},s_{4}^{\prime}$ can take at most two distinct values. So, the same result holds for the coefficients $s_{1},s_{2},s_{3},s_{4}$. Now thanks to the fourth point of Proposition 4.3, we only have the following possibilities:

1. $\bullet$

if $s_{1}=s_{2}=s_{3}=s_{4}=0$, then

[TABLE] 2. $\bullet$

if $s_{1}=3\alpha/4$ and $s_{2}=s_{3}=s_{4}=-\alpha/4$ for $\alpha\in\mathbb{R}$, then

[TABLE] 3. $\bullet$

if $s_{2}=3\alpha/4$ and $s_{1}=s_{3}=s_{4}=-\alpha/4$ for $\alpha\in\mathbb{R}$, then

[TABLE]

and similarly by permuting the diagonal elements when $s_{3}=3\alpha/4$ and when the other elements are equal $s_{1}=s_{2}=s_{4}=-\alpha/4$ or when $s_{4}=3\alpha/4$ and $s_{1}=s_{2}=s_{3}=-\alpha/4$, 4. $\bullet$

if $s_{1}=s_{2}=\alpha/4$ and $s_{3}=s_{4}=-\alpha/4$ for $\alpha\in\mathbb{R}$, then

[TABLE]

and similarly by permuting the diagonal elements when $s_{1}=s_{3}=\alpha/4$ and when $s_{2}=s_{4}=-\alpha/4$ or $s_{1}=s_{4}=\alpha/4$ and $s_{2}=s_{3}=-\alpha/4$.

The computation of the SSVD for these matrices is an easy computation. This concludes the proof of Proposition 4.2.∎

4.3 Determination of the equilibria for each density $\rho$

In Theorem 5 we saw that the BGK operator can have three types of equilibria. The uniform equilibria $f=\rho$ (corresponding to $J=0$) is always an equilibrium. However, the existence of the other two types of equilibria depends on Equations (20) and (21) having a solution for a given $\rho$. Therefore the existence of these types of equilibria will depend on the value of $\rho$. In this section we will determine the existing equilibria for each value of $\rho$. In particular, we will draw the phase diagram for $\rho$ and $\alpha$, that is to say the parametrised curves defined by Equations (20) and (21) in the plane $(\rho,\alpha)$ (see Figure 2). We first prove the following proposition.

Proposition 4.4**.**

Let $\rho_{c}:=6$.

1. (i)

The function $\alpha\mapsto\alpha/c_{1}(\alpha)$ is well-defined on $\mathbb{R}$, its value at zero is $\rho_{c}$. Moreover, there exists $\alpha^{*}>0$ such that this function is decreasing on $(-\infty,\alpha^{*}]$ and increasing on $[\alpha^{*},+\infty)$. Defining $\rho^{*}:=\alpha^{*}/c_{1}(\alpha^{*})$, it holds that $\rho^{*}<\rho_{c}$. 2. (ii)

The function $\alpha\mapsto\alpha/c_{2}(\alpha)$ is even. It is decreasing on $(-\infty,0)$, increasing on $(0,\infty)$ and its value at zero is $\rho_{c}$. 3. (iii)

We have the following asymptotic behaviours:

[TABLE]

[TABLE]

Proof.

The idea of the proof is taken from [55].

1. (i)

Since

[TABLE]

an integration by parts shows that:

[TABLE]

It proves that $\alpha$ and $c_{1}(\alpha)$ have the same sign for all $\alpha\in\mathbb{R}$. Then we define function $m\leavevmode\nobreak\ :\alpha\mapsto\{\sin^{2}\theta\}_{\alpha}=3c_{1}(\alpha)/\alpha$ which satisfies the property:

[TABLE]

since

[TABLE]

and

[TABLE]

where $\operatorname{Var}_{\alpha}$ is the variance for the probability density (15). This property implies that $\alpha/c_{1}(\alpha)$ has only one critical point which is a global minimum. This minimum is attained at a point $\alpha^{*}>0$ as a simple computation shows that $m^{\prime}(0)>0$ and consequently $\rho^{*}<\rho_{c}$. A simple computation gives $m(0)=\frac{1}{2}$ so $\rho_{c}=6$. 2. (ii)

We have similarly:

[TABLE]

from which we can easily see that $\alpha\mapsto\alpha/c_{2}(\alpha)$ is even and has only one minimum attained at $\alpha=0$. A simple computation shows that its value at 0 is $\rho_{c}=6$. 3. (iii)

The behaviour at infinity is obtained by Laplace’s method: with the change of variable $s=1-\cos\theta$ on $[0,\pi]$, we get

[TABLE]

With the same method we have:

[TABLE]

∎

Thanks to Proposition 4.4 and Theorem 5 we can now fully describe the equilibria of the BGK operator. A graphical representation of this result is given by the phase diagram depicted in Figure 2 :

Corollary 4.2** (Equilibria of the BGK operator, depending on the density $\rho$).**

The set of equilibria of the BGK operator (13) depends on the value of $\rho$. In particular we need to distinguish three regions $\rho\in(0,\rho^{*})$, $\rho\in(\rho^{*},\rho_{c})$ and $\rho>\rho_{c}$ where $\rho^{*}$ and $\rho_{c}$ are defined in Proposition 4.4. We have the following equilibria in each region:

1. $\bullet$

For $0<\rho<\rho^{*}$, $\alpha=0$ is the unique solution of Equations (20) and (21) and therefore the only equilibrium is the uniform equilibrium $f^{\mathrm{eq}}=\rho$. 2. $\bullet$

For $\rho=\rho^{*}$, in addition to the uniform equilibrium, there is a family of anisotropic equilibria given by $f^{\mathrm{eq}}=\rho^{*}M_{\alpha^{*}\Lambda}$ where $\Lambda\in SO_{3}(\mathbb{R})$ and $\alpha^{*}=\rho^{*}c_{1}(\alpha^{*})$. 3. $\bullet$

For $\rho^{*}<\rho<\rho_{c}$, the compatibility equation (20) has two solutions $\alpha_{+}$ and $\alpha_{-}$ with $0<\alpha_{-}<\alpha_{+}$ which give, in addition to the uniform equilibrium, two families of anisotropic equilibria : $f^{\mathrm{eq}}=\rho M_{\alpha_{+}\Lambda}$ and $f^{\mathrm{eq}}=\rho M_{\alpha_{-}\Lambda}$ with $\Lambda\in SO_{3}(\mathbb{R})$. 4. $\bullet$

For $\rho=\rho_{c}$, we have $\alpha_{-}=0$. 5. $\bullet$

For $\rho>\rho_{c}$, Equation (20) has two solutions $\alpha_{3}<0<\alpha_{1}$ which give two families of anisotropic equilibria $f^{\mathrm{eq}}=\rho M_{\alpha_{3}\Lambda}$ and $f^{\mathrm{eq}}=\rho M_{\alpha_{1}\Lambda}$ with $\Lambda\in SO_{3}(\mathbb{R})$. Moreover, Equation (21) has two solutions $-\alpha_{2}<0<\alpha_{2}$ which give another family of equilibria: $f^{\mathrm{eq}}=\rho M_{\alpha_{2}B}$ where $B\in\mathscr{B}$. The uniform equilibrium is always an equilibrium.

When an equilibrium is of the form $f^{\mathrm{eq}}=\rho M_{\alpha\Lambda}$ with parameters $\alpha>0$ and $\Lambda\in SO_{3}(\mathbb{R})$ then these parameters can respectively be interpreted as a concentration parameter and a mean body-orientation. They are analogous to the equilibria found in [13] in the Vicsek case. However, in $SO_{3}(\mathbb{R})$, there exist other equilibria which are not of this form. We will see in Section 5 that these latter equilibria are always unstable.

Finally the following picture (Figure 3) is a representation in the space $\mathbb{R}^{3}$ of the diagonal parts of the SSVDs of the solution of the matrix compatibility equation (19) when $\rho>\rho_{c}$. They all belong to the domain $\mathscr{D}$ defined by (12) and depicted in orange in Figure 3.

5 Convergence to equilibria

Now that we know all the equilibria of the spatially homogeneous BGK equation (7) we proceed to investigate the asymptotic behaviour of $f(t,A)$ as $t\to+\infty$. This problem can be reduced to looking at the asymptotic behaviour of $J_{f}$ since, if $J_{f}\to J_{\infty}\in\mathscr{M}_{3}(\mathbb{R})$, then $f(t)$ will converge as $t\to+\infty$ towards $\rho M_{J_{\infty}}$ as it can be seen by writing Duhamel’s formula for equation (7) :

[TABLE]

The asymptotic behaviour of $J_{f}$ is much simpler than the one of $f$ since $J_{f}$ is the solution of the following ODE

[TABLE]

as it can be seen by multiplying (7) by $A\in SO_{3}(\mathbb{R})$ and integrating over $SO_{3}(\mathbb{R})$. Since $J\in\mathscr{M}_{3}(\mathbb{R})\mapsto M_{J}\in L^{\infty}(SO_{3}(\mathbb{R}))$ is locally Lipschitz, the flow of Equation (28) is defined globally in time as a bounded Lipschitz perturbation of the linear system $\frac{d}{dt}J=-J$.

Notice that the solutions of the compatibility equation (19) are exactly the equilibria of the dynamical system (28). We therefore obtain the following proposition:

Proposition 5.1** (Equilibria of the BGK operator, equilibria of the ODE (28)).**

A distribution $f^{\mathrm{eq}}=\rho M_{J}$ is an equilibrium of the BGK operator (13) if and only if $J\in\mathscr{M}_{3}(\mathbb{R})$ is an equilibrium of the dynamical system (28).

We will call stable/unstable an equilibrium of the BGK operator (13) such that the associated matrix $J\in\mathscr{M}_{3}(\mathbb{R})$ is a stable/unstable equilibrium of the ODE (28). This section is devoted to the proof of the following theorem:

Theorem 7** (Convergence towards equilibria).**

Let $\rho\in\mathbb{R}_{+}$ be such that $\rho\neq\rho^{*}$ and $\rho\leavevmode\nobreak\ \neq\leavevmode\nobreak\ \rho_{c}$ (as defined in Proposition 4.4). Let $f_{0}$ be an initial condition for (7) and let $J_{f_{0}}\leavevmode\nobreak\ =\leavevmode\nobreak\ PD_{0}Q$ be a SSVD. Let $f(t)$ be the solution at time $t\in\mathbb{R}_{+}$ of the spatially homogeneous BGK equation (7) with initial condition $f_{0}$. Let $D(t)$ be the solution at time $t\in\mathbb{R}_{+}$ of the ODE (28) with initial condition $D_{0}\in\mathscr{D}$. It holds that:

[TABLE]

is a SSVD and there exists a subset $\mathscr{N}_{\rho}\subset\mathbb{R}^{3}$ of zero Lebesgue measure such that:

if $D_{0}\notin\mathcal{N}_{\rho}$, then $f(t)$ converges as $t\to\infty$ towards an equilibrium $f^{\mathrm{eq}}$ of the BGK operator (13) of the form $f^{\mathrm{eq}}=\rho M_{J^{\mathrm{eq}}}$, where $J^{\mathrm{eq}}\in\mathscr{M}_{3}(\mathbb{R})$ is of one of the forms described in Theorem 5. The convergence is locally exponentially fast in the sense that there exists constants $\delta,K,\mu>0$ such that if $\|J_{f_{0}}-J^{\mathrm{eq}}\|\leq\delta$ then for all $t>0$,

[TABLE] 2. 2.

If $D_{0}\notin\mathcal{N}_{\rho}$, we have the following asymptotic behaviours depending on the density $\rho$ :

1. (i)

if $0<\rho<\rho^{*}$, then $D(t)\to 0$ as $t\to+\infty$ and, consequently, $f^{\mathrm{eq}}=\rho$, 2. (ii)

if $\rho^{*}<\rho<\rho_{c}$, then $D(t)\to 0$ or $D(t)\to\alpha_{+}I_{3}$ as $t\to+\infty$ and, consequently, $f^{\mathrm{eq}}=\rho$ or $f^{\mathrm{eq}}=\rho M_{\alpha_{+}\Lambda}$ respectively, where $\alpha_{+}>0$ is defined in corollary 4.2 and $\Lambda:=PQ\in SO_{3}(\mathbb{R})$, 3. (iii)

if $\rho>\rho_{c}$, then $D(t)\to\alpha_{1}I_{3}$ as $t\to+\infty$ and, consequently, $f^{\mathrm{eq}}=\rho M_{\alpha_{1}\Lambda}$ where $\alpha_{1}>0$ is defined in corollary 4.2 and $\Lambda:=PQ\in SO_{3}(\mathbb{R})$.

Remark 5.1**.**

The subset $\mathscr{N}_{\rho}\subset\mathbb{R}^{3}$ depends on the density $\rho$. This subset will be made explicit in the three cases $\rho<\rho^{*}$, $\rho^{*}<\rho<\rho_{c}$ and $\rho>\rho_{c}$ in Section 5.3. If $D_{0}\in\mathscr{N}_{\rho}$, then there is convergence towards an unstable equilibrium at a rate which may not be exponential.

Remark 5.2** (Phase transitions).**

Theorem 7 demonstrates a phase transition phenomenon triggered by the density of agents $\rho$ : when $\rho<\rho^{*}$, the system is disordered (asymptotically in time) in the sense that $J_{f}\to 0$ and we therefore cannot define a mean body-attitude. When the density increases and exceeds the critical value $\rho_{c}$, the system is self-organised (asymptotically in time and for almost every initial data), in the sense that $J_{f}\to\alpha\Lambda$ where $\alpha\in\mathbb{R}_{+}$ and $\Lambda\in SO_{3}(\mathbb{R})$ can be respectively interpreted as a concentration parameter and a mean body attitude. When $\rho^{*}<\rho<\rho_{c}$ the self-organised and disordered states are both asymptotically stable and the convergence towards one or the other state depends on the initial data. Such “transition region” also appears in the Vicsek model, as studied in [13], and gives rise to an hysteresis phenomenon.

The proof of this theorem can be found at the end of Section 5.2. It is based on a gradient flow structure for the flux $J_{f}$ studied in Section 5.1. This structure ensures the convergence of $J_{f}$ towards a matrix $J^{\mathrm{eq}}\in\mathscr{M}_{3}(\mathbb{R})$ as $t\to+\infty$ and consequently the convergence of $f(t)$ as $t\to+\infty$ towards an equilibrium. The stability of the equilibria determines which equilibrium can be attained. This question is addressed in Section 5.2. Additional details about the subset $\mathscr{N}_{\rho}$ as well as the study of the critical cases $\rho=\rho^{*}$ and $\rho=\rho_{c}$ are provided in Section 5.3.

5.1 A gradient-flow structure in $\mathbb{R}^{3}$

In this section we show that the ODE (28) can be reduced to a gradient-flow ODE in $\mathbb{R}^{3}$. We first show (Section 5.1.1) how (28) can be reduced to an ODE in $\mathbb{R}^{3}$, the equilibria of which are linked to the equilibria of (28) (and therefore of (7)). Then we show that this ODE in $\mathbb{R}^{3}$ has a gradient-flow structure which will allow us to conclude on the asymptotic behaviour of the solution of (7) (Section 5.1.2).

5.1.1 Reduction to a nonlinear ODE in $\mathbb{R}^{3}$ and equilibria

The ODE (28) is a matrix-valued nonlinear ODE (in dimension 9) but, as in the previous section (Proposition 4.1 and Corollary 4.1), we will use the left and right invariance of the Haar measure and the SSVD to reduce the problem to a vector-valued nonlinear ODE in dimension 3, as explained in Proposition 5.2 and Corollary 5.1.

Proposition 5.2** (Reduction to a nonlinear ODE in dimension 3).**

Let $J_{0}\in\mathscr{M}_{3}(\mathbb{R})$ be a given matrix and let $D_{0}\in\operatorname{Orb}(J_{0})$ be diagonal. Let $P,Q\in SO_{3}(\mathbb{R})$ such that $J_{0}=PD_{0}Q$. Let $J:[0,\infty)\to\mathscr{M}_{3}(\mathbb{R})$ be a $C^{1}$ curve in $\mathscr{M}_{3}(\mathbb{R})$ with $J(0)=J_{0}$. For all $t>0$, let $D(t)\in\mathscr{M}_{3}(\mathbb{R})$ be the matrix such that $J(t)=PD(t)Q$. It holds that:

1. (i)

$J=J(t)$* is the solution of the ODE (28) with initial condition $J(t=0)=J_{0}$ if and only if $D=D(t)$ is the solution of the same ODE (28) with initial condition $D(t=0)=D_{0}$.* 2. (ii)

Moreover, if (i) holds, then the matrix $D(t)\in\mathscr{M}_{3}(\mathbb{R})$ is diagonal for all time.

Proof.

1. (i)

Using the left and right invariance of the Haar measure, we see that if the matrix $J(t)\in\mathscr{M}_{3}(\mathbb{R})$ is a solution of (28), then for any $P,Q\in SO_{3}(\mathbb{R})$, $PJ(t)Q$ is also a solution (and conversely). 2. (ii)

Since $D_{0}$ is diagonal, the fact that $D(t)$ is also diagonal is a consequence of lemma 3.2 which states that $\langle A\rangle_{M_{D}}$ is diagonal when $D$ is diagonal.

∎

For any matrix $J_{0}$, such a diagonal matrix $D_{0}$ always exists: we can take the diagonal part of the SSVD of $J_{0}$. We therefore only have to study the following ODE for diagonal matrices $\mathscr{D}_{3}(\mathbb{R})$. Since this vector space is isomorphic to $\mathbb{R}^{3}$ through the isomorphism $\operatorname{diag}$ defined in the introduction, we obtain two equivalent ODEs:

[TABLE]

where $\widehat{D}_{0}=\operatorname{diag}^{-1}(D_{0})$ and the following equivalence holds : $\widehat{D}(t)=\operatorname{diag}^{-1}(D(t))$ is the solution of (29b) if and only if $D(t)=\operatorname{diag}(\widehat{D}(t))$ is the solution of (29a).

Note that it is not clear that if $J_{0}=PD_{0}Q$ is a SSVD, then $J(t)=PD(t)Q$ is a SSVD for all $t>0$. The following proposition and corollary ensure that the SSVD is preserved by the dynamical system which will allow us to restrict the domain on which the ODE (29a) is posed.

Proposition 5.3** (Invariant manifolds).**

The following subsets of $\mathbb{R}^{3}$ are invariant manifolds of the dynamical system (29b) :

1. $\bullet$

the planes \big{\{}(d_{1},d_{2},d_{3})\in\mathbb{R}^{3},\,\,d_{i}+d_{j}=0\big{\}} for $i\neq j\in\{1,2,3\}$, 2. $\bullet$

the planes \big{\{}(d_{1},d_{2},d_{3})\in\mathbb{R}^{3},\,\,d_{i}-d_{j}=0\big{\}} for $i\neq j\in\{1,2,3\}$, 3. $\bullet$

the intersections of two of these planes and in particular the lines

[TABLE]

Proof.

For $i=2$ and $j=3$, the result has already been proved in the second point of Lemma 4.1. The other cases are similar.∎

Corollary 5.1**.**

Let $J_{0}\in\mathscr{M}_{3}(\mathbb{R})$ and $J_{0}=PD_{0}Q$ be a SSVD with $P,Q\in SO_{3}(\mathbb{R})$ and $D_{0}\in\mathscr{D}_{3}(\mathbb{R})$ diagonal. Let $J(t)$ be the solution of the ODE (28) with initial condition $J(t=0)=J_{0}$. Let $D(t)$ the solution of the ODE (29a) with initial condition $D(t=0)=D_{0}$. Then the decomposition $J(t)=PD(t)Q$ is a SSVD for $J(t)$.

Proof.

The fact that $J(t)=PD(t)Q$ is a consequence of the first point of Proposition 5.2. The fact that it is a SSVD is a consequence of Proposition 5.3 which ensures that the conditions (12) remain true for all $t>0$ : $\mathscr{D}$ is stable in the sense that if $D_{0}\in\mathscr{D}$, then $D(t)\in\mathscr{D}$ for all time $t>0$. This follows from the fact that the image by the isomorphism $\operatorname{diag}$ of the invariant manifolds of (29b) described in Proposition 5.3 are invariant manifolds of the dynamical system (29a). These manifolds in $\mathscr{D}_{3}(\mathbb{R})$ form the boundary of the subset $\mathscr{D}$.∎

In conclusion, the study of the asymptotic behaviour of $f(t)$ as $t\to+\infty$ can be reduced to the study of the asymptotic behaviour of the solutions of the ODE (29a) posed on the domain $\mathscr{D}$ (see Figure 3).

The following Proposition describes the equilibria of the dynamical system (32a) and is a consequence of the results of Section 4.

Proposition 5.4** (Equilibria of the dynamical system (29a)).**

The equilibria of the dynamical system (29a) are, depending on the density $\rho$ :

1. $\bullet$

the matrix $D=0$ for any density $\rho$, 2. $\bullet$

the diagonal matrices of type (b) with parameters $\alpha_{+}$ and $\alpha_{-}$ when $\rho^{*}<\rho<\rho_{c}$, 3. $\bullet$

the diagonal matrices of type (b) with parameters $\alpha_{1}$ and $\alpha_{3}$ and the diagonal matrices of type (c) with parameters $\pm\alpha_{2}$ when $\rho>\rho_{c}$,

where the types (b) and (c) are defined in Proposition 4.2 and the elements $\alpha_{+}$, $\alpha_{-}$, $\alpha_{1}$, $\alpha_{2}$ and $\alpha_{3}$ are defined in Corollary 4.2.

Proof.

The equilibria of (29a) are the diagonal matrices $D$ such that

[TABLE]

which is exactly Equation (19) for the diagonal matrices. This equation has been solved in Theorem 5, Remark 4.7 and Corollary 4.2. ∎

Remark 5.3**.**

As in the previous section, we have found all the diagonal equilibria of (29a). However, thanks to Corollary 5.1, only the ones which belong to $\mathscr{D}$ are needed.

The following proposition is a straightforward consequence of the previous results.

Proposition 5.5** (Equilibria of the BGK operator, equilibria of the ODE (29a)).**

Let $J\in\leavevmode\nobreak\ \mathscr{M}_{3}(\mathbb{R})$ with a SSVD given by $J=PDQ$. The following assertions are equivalent.

1. (1)

The distribution $f^{\mathrm{eq}}=\rho M_{J}$ is an equilibrium of the BGK operator (13). 2. (2)

The matrix $D$ is an equilibrium of the dynamical system (29a) on the domain $\mathscr{D}$. 3. (3)

The matrix $D$ is of one of the forms:

1. $\bullet$

$D=0$, when $\rho<\rho^{*}$, 2. $\bullet$

$D=0$* or $D=\alpha_{-}I_{3}$ or $D=\alpha_{+}I_{3}$, when $\rho^{*}\leq\rho\leq\rho_{c}$,* 3. $\bullet$

$D=0$* or $D=\alpha_{1}I_{3}$ or $D=\alpha_{3}\operatorname{diag}(-1,-1,1)$ or $D=\alpha_{2}\operatorname{diag}(1,0,0)$, when $\rho>\rho_{c}$,*

where $\rho^{*}$, $\rho_{c}$, $\alpha_{-}$, $\alpha_{+}$, $\alpha_{1}$, $\alpha_{2}$ and $\alpha_{3}$ are defined in Proposition 4.4 and Corollary 4.2.

5.1.2 A gradient-flow structure

Lemma 5.1** (Gradient-flow structure).**

We define the partition function of a matrix $J\leavevmode\nobreak\ \in\leavevmode\nobreak\ \mathscr{M}_{3}(\mathbb{R})$ :

[TABLE]

and the potentials $V(J)$ and $\widehat{V}(\widehat{D})$ respectively on $\mathscr{M}_{3}(\mathbb{R})$ and on $\mathbb{R}^{3}$ :

[TABLE]

where $\|\cdot\|$ and $|\cdot|$ are the Euclidean norms respectively on $\mathscr{M}_{3}(\mathbb{R})$ and on $\mathbb{R}^{3}$. Then we can rewrite equations (29a) and (29b) into a gradient flow structure as follows:

[TABLE]

where $\nabla$ is the gradient operator in $\mathscr{M}_{3}(\mathbb{R})$ endowed with the Riemaniann structure (8) or the gradient operator in $\mathbb{R}^{3}$ endowed with the usual Euclidean structure.

Proof.

The partition function satisfies that for all $J\in\mathscr{M}_{3}(\mathbb{R})$,

[TABLE]

since $\nabla(e^{J\cdot A})=Ae^{J\cdot A}$. The result follows in $\mathscr{M}_{3}(\mathbb{R})$. The result in $\mathbb{R}^{3}$ follows from the fact that for any $w_{1},w_{2}\in\mathbb{R}^{3}$, it holds that:

[TABLE]

where $\cdot$ denotes the dot product on $\mathbb{R}^{3}$ and on $\mathscr{M}_{3}(\mathbb{R})$ as defined in (8). ∎

Remark 5.4**.**

This gradient-flow structure on $J_{f}$ is specific to the BGK equation and does not hold for the Fokker-Planck operator (as shown in [13], the differential equation satisfied by $J_{f}$ in the Vicsek case involves the spherical harmonics of degree 2 and higher of $f$; here the equation for $J_{f}$ is closed).

In particular the gradient-flow structure (32a) implies that the system (32a) will converge towards an equilibrium. When all the equilibria of the dynamical system (29a) are hyperbolic the convergence towards the stable equilibria is exponentially fast (see [38, Section 9.3]). The goal of the next section is to find which equilibria among the ones found in Section 4 are stable and to completely describe the asymptotic behaviour of the system depending on the initial condition and the local density $\rho$. We will see that phase transitions appear between ordered and disordered dynamics when the density $\rho$ increases.

5.2 Stability of the equilibria and conclusion

Since the flow of (29b) is the image by the isomorphism $\operatorname{diag}^{-1}$ of the flow of (29a), an equilibria $D^{\mathrm{eq}}$ of (29a) is stable (resp. unstable) if and only if $\operatorname{diag}^{-1}(D^{\mathrm{eq}})$ is a stable (resp. unstable) equilibria of (29b). The stability properties of the equilibria of (29b) are much simpler to study than the stability properties of the equilibria of (29a) since they are given by the signature of the Hessian matrix $\operatorname{Hess}\widehat{V}(\widehat{D}^{\mathrm{eq}})\in\mathscr{S}_{3}(\mathbb{R})$ of the potential $\widehat{V}$ given in equation (31b) (the linearisation of ODE (29b) around an equilibrium $\widehat{D}^{\mathrm{eq}}$ is indeed $\frac{d}{dt}\widehat{H}=-\operatorname{Hess}\widehat{V}(\widehat{D}^{\mathrm{eq}})\widehat{H}$). In particular, an equilibrium $\widehat{D}^{\mathrm{eq}}$ is stable if and only if the signature of $\operatorname{Hess}\widehat{V}(\widehat{D}^{\mathrm{eq}})$ is $(+++)$.

Note that in the matrix framework (32a), the Hessian of the potential (31a) in the Euclidean space $\mathscr{M}_{3}(\mathbb{R})$ would be a rank 4 tensor. Here we are reduced to the computation of the signature of $3\times 3$ symmetric matrices. For a diagonal matrix $D\in\mathscr{D}_{3}(\mathbb{R})$ and $\widehat{D}=\operatorname{diag}^{-1}(D)$ we will write with a slight abuse of notations:

[TABLE]

When $D\in\mathscr{D}_{3}(\mathbb{R})$ is an equilibrium of (29a) we call signature of the Hessian matrix $\operatorname{Hess}V(D)$ the signature of $\operatorname{Hess}\widehat{V}(\widehat{D})$ where $\widehat{D}=\operatorname{diag}^{-1}(D)$. A simple computation shows that the Hessian matrix $\operatorname{Hess}V(D)$ is given by :

[TABLE]

where $\Gamma=(\Gamma_{ij})_{i,j}$ with :

[TABLE]

The following theorem is an extension of Corollary 4.2 and gives a full description of the equilibria of the BGK operator with their stability.

Theorem 8** (Stability of the equilibria of the ODE (29a)).**

1. $\bullet$

For $0<\rho<\rho^{*}$, the only equilibrium is $D=0$. This equilibrium is stable. 2. $\bullet$

For $\rho^{*}<\rho<\rho_{c}$, the equilibrium $D=0$ and the equilibria of type (b) with parameter $\alpha_{+}$ are stable. The equilibria of type (b) with parameter $\alpha_{-}$ are unstable and the signature of the Hessian matrix is $(-++)$. 3. $\bullet$

For $\rho>\rho_{c}$, the stable equilibria are the equilibria of type (b) with parameter $\alpha_{1}$. The other equilibria ($D=0$, type (b) with parameter $\alpha_{3}$ and type (c) with parameter $\pm\alpha_{2}$) are unstable and the signatures of the Hessian matrix are respectively $(---)$, $(+--)$ and $(++-)$.

The proof of Theorem 8 will be based on the following lemma which states an important orbital invariance principle for the signature of the Hessian matrix.

Lemma 5.2** (Orbital invariance of the signature).**

Let $\widehat{D}\in\mathbb{R}^{3}$ be an equilibrium of (29b). The signature of $\operatorname{Hess}\widehat{V}(\widehat{D})$ depends only on the type (a), (b) or (c) of $\operatorname{diag}(\widehat{D})$ as defined in Proposition 4.2.

Proof.

Let $\widehat{D}_{1}$ and $\widehat{D}_{2}$ be two equilibria of (29b) such that $\operatorname{diag}(\widehat{D}_{1})$ and $\operatorname{diag}(\widehat{D}_{2})$ are of the same type. Then there exist $P,Q\in SO_{3}(\mathbb{R})$ and $R\in\mathcal{O}_{3}(\mathbb{R})$ such that

[TABLE]

Moreover for all $M\in\mathscr{M}_{3}(\mathbb{R})$, the matrices $P,Q,R$ satisfy (permutation of the diagonal coefficients):

[TABLE]

where $M_{ii}$ is the $(i,i)$-th coefficient of the matrix $M$. Therefore, one can check that:

[TABLE]

and

[TABLE]

The conclusion follows from Sylvester’s law of inertia ([45, Chapter 8 Theorem 1]).∎

Proof (of Theorem 8).

Thanks to lemma 5.2, we therefore do not have to compute the signatures of the Hessian matrix taken at all equilibria, it is enough to choose one matrix in each orbit: for the equilibria of type (b) (see Proposition 4.2) we will compute the signature of $\operatorname{Hess}(\alpha I_{3})$ where $\alpha=\rho c_{1}(\alpha)$ and for the equilibria of type (c) we will compute the signature of $\operatorname{Hess}(\alpha\operatorname{diag}(1,0,0))$ where $\alpha=\rho c_{2}(\alpha)$. We first start with the case of the uniform equilibrium.

**Case 1. Uniform equilibrium $D=0$

**

For the uniform equilibrium $D=0$, $\langle a_{ij}\rangle=0$ for all $(i,j)$. Moreover, by the change of variable $A^{\prime}=D^{ik}A$ where $k\neq i,j$, we have when $i\neq j$ :

[TABLE]

It proves that $\Gamma$ is diagonal. Then with the changes of variables $A^{\prime}=P^{ij}A$ or $A^{\prime}=AP^{ij}$ it can be seen that all the $3^{2}$ quantities $\langle a_{ij}^{2}\rangle$ are equal. Since their sum is equal to $n=3$ we get that

[TABLE]

where $\rho_{c}=6$. In conclusion, the signature of $\operatorname{Hess}V(0)$ is $(+++)$ if $\rho<\rho_{c}$ and $(---)$ if $\rho>\rho_{c}$.

**Case 2. Equilibria of type (b) : $D=\alpha I_{3}$

**

Let $D=\alpha I_{3}$ with $\alpha=\rho c_{1}(\alpha)$. We have $c_{1}(\alpha)=\langle a_{11}\rangle_{M_{D}}=\langle a_{22}\rangle_{M_{D}}=\langle a_{33}\rangle_{M_{D}}$ and a change of variable of the type $A^{\prime}=P^{ij}A(P^{ij})^{T}$ shows that all the $\langle a_{kk}a_{\ell\ell}\rangle_{M_{D}}$ are equal. The Hessian matrix is therefore equal to :

[TABLE]

with $\nu=\langle a_{11}^{2}\rangle_{M_{D}}-\langle a_{11}\rangle_{M_{D}}^{2}$ and $\gamma=\langle a_{11}a_{22}\rangle_{M_{D}}-\langle a_{11}\rangle_{M_{D}}\langle a_{22}\rangle_{M_{D}}$. The eigenvalues of $\operatorname{Hess}V(D)$ are :

1. $\bullet$

$1-\frac{1}{2}\rho\nu-\rho\gamma$ of order 1 with eigenvector $(1,1,1)^{T}$. But taking the derivative with respect to $\alpha$ of (14) with $\Lambda=I_{3}$ we obtain the relation :

[TABLE]

and using $\rho=\alpha/c_{1}(\alpha)$ we can rewrite :

[TABLE]

Its sign is then given by Proposition 4.4 and the fact that $c_{1}(\alpha)$ has the same sign as $\alpha$ : $c_{1}(\alpha)\left(\frac{id}{c_{1}}\right)^{\prime}(\alpha)>0$ when $\alpha<0$, $c_{1}(\alpha)\left(\frac{id}{c_{1}}\right)^{\prime}(\alpha)<0$ when $0\leq\alpha<\alpha^{*}$ and $c_{1}(\alpha)\left(\frac{id}{c_{1}}\right)^{\prime}(\alpha)>0$ when $\alpha>\alpha^{*}$. 2. $\bullet$

$1-\frac{1}{2}\rho\nu+\frac{1}{2}\rho\gamma$ of order 2 with eigenvectors $(1,-1,0)^{T}$ and $(0,1,-1)^{T}$. It can be rewritten as:

[TABLE]

To determine this sign, we use the explicit volume form of $SO_{3}(\mathbb{R})$ given by Rodrigues’ formula (9) to see that :

[TABLE]

Lemma 5.3**.**

The function

[TABLE]

satisfies $f(0)=0$, $f(x)\geq 0$ if $x\geq 0$ and $f(x)\leq 0$ if $x\leq 0$

Proof.

The value of $f(0)$ is given by the expansion of $\exp$. Note that :

[TABLE]

so that an integration by parts shows :

[TABLE]

We get that :

[TABLE]

We have :

[TABLE]

Linearizing $\sin^{2}(\theta/2)$ and expanding everything gives :

[TABLE]

so that :

[TABLE]

Now for $x\geq 0$, let $\theta_{0}=\arccos(-2/3)$. We cut the integral at $\theta_{0}$ and we get

[TABLE]

since the integrand is nonnegative and $x\geq 0$. Similarly when the integrand is nonpositive

[TABLE]

And finally :

[TABLE]

We find similarly that $f(x)\leq 0$ when $x\leq 0$. ∎

And therefore, we can deduce the sign of the eigenvalue: 1-\frac{\alpha}{4c_{1}(\alpha)}\big{\langle}(a_{11}-a_{22})^{2}\big{\rangle}_{M_{D}}>0 when $\alpha>0$ and 1-\frac{\alpha}{4c_{1}(\alpha)}\Big{\langle}(a_{11}-a_{22})^{2}\big{\rangle}_{M_{D}}<0 when $\alpha<0$.

**Case 3. Equilibria of type (c) : $D=\alpha\operatorname{diag}(1,0,0)$

**

For $D=\alpha\operatorname{diag}(1,0,0)$ with $\alpha=\rho c_{2}(\alpha)$, using the parametrisation (11), it holds that :

[TABLE]

where the first equality comes from the change of variable $A\mapsto P^{13}A(P^{13})^{T}$. Similarly,

[TABLE]

The Hessian matrix is therefore equal to :

[TABLE]

and since $\langle a_{22}^{2}\rangle_{M_{D}}=\langle a_{33}^{2}\rangle_{M_{D}}$ as it can be seen with the change of variable $A\mapsto P^{23}A(P^{23})^{T}$, its eigenvalues are :

[TABLE]

with eigenvector $(1,0,0)^{T}$,

[TABLE]

with eigenvector $(0,1,-1)^{T}$ and

[TABLE]

with eigenvector $(0,1,1)^{T}$. We have as before :

[TABLE]

so the first eigenvalue can be rewritten

[TABLE]

where we have used Proposition 4.4 to determine the sign. The two other eigenvalues are equal to:

[TABLE]

Using the change of variable $A\mapsto P^{13}A(P^{13})^{T}$ and the parametrisation (11), one can see that :

[TABLE]

so that :

[TABLE]

where $[\cdot]_{\alpha}$ is defined in Proposition 4.1. Using the relation $\rho=\frac{\alpha}{c_{2}(\alpha)}=\frac{4}{[\sin^{2}\phi]_{\alpha}}$ (see Formula (26)), the two eigenvalues are equal to :

[TABLE]

With the same technique used in the previous paragraph (Lemma 5.3) it is possible to show that the eigenvalue 1-\frac{1}{4}\rho\Big{\langle}(a_{22}-a_{33})^{2}\Big{\rangle}_{M_{D}} is nonpositive when $\alpha<0$ and nonnegative when $\alpha>0$. The contrary holds for the other eigenvalue (nonnegative when $\alpha<0$ and nonpositive when $\alpha>0$). Finally the signs of these two eigenvalues are always $(+-)$.

The conclusion of the proof follows from the study of the roots of Equations (20) and (21) provided by Corollary 4.2 and depicted in Figure 2. In particular, the cases $\alpha=\alpha^{*}$ and $\alpha=0$ correspond to the critical cases $\rho=\rho^{*}$ and $\rho=\rho_{c}$ which are studied in the next section. ∎

Remark 5.5**.**

The above calculations are similar to the computation of the equilibria and their stability of

[TABLE]

where $\mathcal{F}$ is the free-energy (6). In particular, it can be shown that the equilibria of $\Psi$ and their stability are the same as the ones described in Theorem 8. In particular, since (6) is also a free energy in the Fokker-Planck case, this analysis shows the instability of equilibria of the Fokker-Planck of the form $\rho M_{D}$ when $D$ is of one of the unstable equilibria of (29a) described in Theorem 8. This technique is similar to the one which was used in [57] in the case of 3D polymers. However, it does not provide global or local convergence of the solution of the Fokker-Planck equation towards an equilibrium. This requires further investigations which will be left for future work. In the Vicsek case [13], it was based on a LaSalle’s principle and on estimates for the dissipation term.

We can finally prove Theorem 7 :

Proof (of Theorem 7).

Thanks to the Duhamel’s formula (27), the asymptotic behaviour of $f(t)$ as $t\to+\infty$ is given by the asymptotic behaviour of $J_{f(t)}$. Thanks to Proposition 5.2 and Corollary 5.1, we only have to study the asymptotic behaviour of the solution of the ODE (29a) where the initial condition is the diagonal part of a SSVD of $J_{f_{0}}$. Since this equation has a gradient-flow structure (32a), we know that the solution $D(t)$ converges as $t\to+\infty$ towards an equilibrium $D^{\mathrm{eq}}$ and consequently $J_{f(t)}\to PD^{\mathrm{eq}}Q:=J^{\mathrm{eq}}$. Moreover the equilibrium $D^{\mathrm{eq}}$ is a stable equilibrium provided that $D_{0}$ does not belong to the stable manifold of an unstable equilibrium. Since these manifolds have dimension at most 2, the union of these manifolds, called $\mathcal{N}_{\rho}$ is of zero measure and there is convergence towards a stable equilibrium for all $D_{0}\notin\mathcal{N}_{\rho}$. In this case, the convergence is locally exponentially fast in the sense that there exist constants $\delta,\lambda,C>0$ such that if $\|J_{f_{0}}-J^{\mathrm{eq}}\|\leq\delta$ then for all $t>0$,

[TABLE]

Let $f^{\mathrm{eq}}=\rho M_{J^{\mathrm{eq}}}$. It follows from Duhamel’s formula (27) that for all $A\in SO_{3}(\mathbb{R})$ :

[TABLE]

Since $SO_{3}(\mathbb{R})$ is compact and $J_{f(t)}$ is bounded uniformly in $t$, there exists a constant $L>0$ such that for all $t>0$, the following Lipschitz bound holds:

[TABLE]

Reporting (35) into (34) and using (33), the first point of Theorem 7 follows with constants $K=CL/|1-\lambda|$ and $\mu=\min(1,\lambda)$ when $\lambda\neq 1$ and $K=CL$ and $\mu=1-\varepsilon$ for any $\varepsilon>0$ when $\lambda=1$.

The stability of the equilibria of the dynamical system (29a) is given in Theorem 8. Finally the conclusions of points 2.(ii) and 2.(iii) follow from the fact that the diagonal parts of the SSVD of the equilibria of type (b) with parameters $\alpha_{+}$ and $\alpha_{1}$ are respectively $\alpha_{+}I_{3}$ and $\alpha_{1}I_{3}$.∎

5.3 Final remark: characterisation of the stable manifold of the unstable equilibria and critical cases

This section is devoted to a more precise description of the subset $\mathscr{N}_{\rho}$ of zero measure of initial conditions which do not necessarily lead to one of the behaviours detailed in the previous section (see Theorem 7).

When $\rho\neq\rho^{*}$ and $\rho\neq\rho_{c}$, all the equilibria of (29b) are hyperbolic and we can apply the stable manifold theorem ([48] Section 2.7) which states that for a given hyperbolic equilibrium $\widehat{D}^{\mathrm{eq}}$ of the nonlinear equation (29b), the stable set of $\widehat{D}^{\mathrm{eq}}$ is a smooth manifold. Its dimension is given by the number of minuses in the signature of $-\operatorname{Hess}\widehat{V}(\widehat{D}^{\mathrm{eq}})$ and its tangent space at $\widehat{D}^{\mathrm{eq}}$ is the stable subspace of the linearized system around $\widehat{D}^{\mathrm{eq}}$. We are now ready to describe the behaviour of the system (29b) for a given density $\rho$ and an initial condition $\widehat{D}_{0}\in\mathbb{R}^{3}$. Note that we can restrict ourselves to the case $\widehat{D}_{0}\in\operatorname{diag}^{-1}(\mathscr{D})$ and that our goal is to describe $\operatorname{diag}^{-1}(\mathscr{N}_{\rho})$ (in the $\mathbb{R}^{3}$-framework). We write $D(t)=\operatorname{diag}\widehat{D}(t)$ the solution of (29a) with initial condition $D_{0}=\operatorname{diag}\widehat{D}_{0}$.

**Case 1. When $\rho<\rho^{*}$

**

When $\rho<\rho^{*}$, there is only the uniform equilibrium : $D(t)\underset{t\to+\infty}{\longrightarrow}0$ and $\mathscr{N}_{\rho}=\emptyset$.

Case 2. When $\rho^{*}<\rho<\rho_{c}$

There are three equilibria for (29b) : two are stable ([math] and $\alpha_{+}(1,1,1)^{T}$) and one is unstable ($\alpha_{-}(1,1,1)^{T}$) which satisfies that $-\operatorname{Hess}(\widehat{V})(\alpha_{-}(1,1,1)^{T})$ has signature $(--+)$. The stable manifold theorem ensures that $\operatorname{diag}^{-1}(\mathscr{N}_{\rho})$ is the stable set of the unstable equilibrium and is a smooth manifold of dimension 2 and its tangent plane at $\alpha_{-}(1,1,1)^{T}$ admits the eigenvectors $(1,-1,0)^{T}$ and $(0,1,-1)^{T}$ for basis. The unstable manifold of the unstable equilibrium has dimension : it is the line of direction $(1,1,1)^{T}$.

Finally, the space $\mathbb{R}^{3}$ is partitioned in two domains by $\operatorname{diag}^{-1}(\mathscr{N}_{\rho})$ : depending on which domain $\widehat{D}_{0}$ belongs to in $\mathbb{R}^{3}\setminus\operatorname{diag}^{-1}(\mathscr{N}_{\rho})$, $\widehat{D}(t)$ will converge towards either the uniform equilibrium or to $\alpha_{+}(1,1,1)^{T}$ (or towards $\alpha_{-}(1,1,1)^{T}$ if $\widehat{D}_{0}$ belongs to its stable manifold).

**Case 3. When $\rho>\rho_{c}$

**

When $\rho>\rho_{c}$, $\widehat{D}(t)$ can converge towards one of the four equilibria : [math], $\alpha_{3}(-1,-1,1)^{T}$, $\alpha_{2}(1,0,0)^{T}$ or $\alpha_{1}(1,1,1)^{T}$ (see Figure 2).

1. $\bullet$

The uniform equilibrium is unstable and $-\operatorname{Hess}(\widehat{V})(0)$ has signature $(+++)$ therefore, there is no stable direction and $D(t)$ cannot converge towards 0 unless $D_{0}=0$ (and then $D(t)=0$ for all $t\in\mathbb{R}_{+}$). 2. $\bullet$

The equilibrium $\alpha_{3}\operatorname{diag}(-1,-1,1)$ is unstable and its stable manifold $\mathcal{M}_{1}$ has dimension 1: it is the half line \mathbb{R}_{+}^{*}\left(\begin{array}[]{c}1\\ 1\\ -1\end{array}\right) (it is indeed an invariant manifold of dimension 1 and a solution on this half-line cannot converge towards an other equilibrium). Therefore $D(t)$ converges towards $\alpha_{3}\operatorname{diag}(-1,-1,1)$ if and only if $\widehat{D}_{0}$ lies on this half line. 3. $\bullet$

The equilibrium $\alpha_{2}(1,0,0)^{T}$ is unstable and its stable manifold has dimension 2 with a tangent plane generated by the two vectors $(1,0,0)^{T}$ and $(0,1,-1)^{T}$ : it is the plane $\{d_{2}+d_{3}=0\}$. We note that this plane is an invariant manifold, so if $\widehat{D}_{0}$ belongs to this plane, the limit when $t\to+\infty$ also belongs to this plane. But since the half-lines \mathbb{R}_{+}^{*}\left(\begin{array}[]{c}1\\ 1\\ -1\end{array}\right) and \mathbb{R}_{+}^{*}\left(\begin{array}[]{c}1\\ -1\\ 1\end{array}\right) are the stable manifolds of the equilibria $\alpha_{3}(-1,-1,1)^{T}$ and $\alpha_{3}(-1,1,-1)^{T}$, the limit when $t\to+\infty$ is $\alpha_{2}(1,0,0)^{T}$ if and only if $D_{0}$ belongs to the (open) quarter of the plane $\{d_{2}+d_{3}=0\}$ delimited by these two half-lines: this is the stable manifold of $\alpha_{2}\operatorname{diag}(1,0,0)$. In conclusion, $D(t)$ converges towards $\alpha_{2}\operatorname{diag}(1,0,0)$ if and only if $\widehat{D}_{0}$ is of the form $(d_{1},d_{2},-d_{2})^{T}$ with $d_{1}>d_{2}\geq 0$. Since $\widehat{D}_{0}\in\operatorname{diag}^{-1}(\mathscr{D})$, we can restrict ourselves to the eighth of plan $\mathcal{M}_{2}$ delimited by the half lines $\mathbb{R}_{+}^{*}(1,1,-1)^{T}$ excluded and $\mathbb{R}_{+}^{*}(1,0,0)^{T}$ included. 4. $\bullet$

In every other cases, $D(t)$ converges towards $\alpha_{1}I_{3}$ which is the only stable equilibrium.

Finally, in the case $\rho>\rho_{c}$, the subset $\mathscr{N}_{\rho}\cap\mathscr{D}$ is equal to \operatorname{diag}\big{(}\{0\}\cup\mathcal{M}_{1}\cup\mathcal{M}_{2}\big{)}.

**The critical cases

**

We end this section by an informal description of the expected behaviour in the critical cases.

1. $\bullet$

When $\rho=\rho^{*}$ the uniform equilibrium is stable and there is an other equilibrium of the form $\alpha^{*}I_{3}$ with $\alpha^{*}>0$. This equilibrium is non hyperbolic: the kernel of $-\operatorname{Hess}\widehat{V}(\alpha^{*}(1,1,1)^{T})$ is one dimensional, spanned by the vector $(1,1,1)^{T}$. The two other eigenvalues are negative. Therefore, for the system (29b), there exists a center manifold of dimension 1 and a stable manifold of dimension 2, which tangent plane being the orthogonal of $(1,1,1)^{T}$.

If $\widehat{D}_{0}$ belongs to the stable manifold, $D(t)$ converges exponentially fast towards $\alpha^{*}I_{3}$ ([33] Theorem 3.22). The stable manifold of dimension 2 delimits two domains and one of them is included in the subset $\{x\in\mathbb{R}^{3},x\cdot(1,1,1)^{T}\geq\alpha^{*}\}$. If $\widehat{D}_{0}$ belongs to this domain, then it belongs to a center manifold and $\widehat{D}(t)$ is attracted exponentially fast towards the the line $\mathbb{R}(1,1,1)^{T}$ and converges at rate $1/t$ towards $\alpha^{*}(1,1,1)$. In every other cases, $\widehat{D}(t)$ converges exponentially fast towards 0. 2. $\bullet$

When $\rho=\rho_{c}$ there is one stable equilibrium of the form $\alpha I_{3}$ with $\alpha>0$. The uniform equilibrium is non hyperbolic and $-\operatorname{Hess}\widehat{V}(0)=0$.

In the case when $\widehat{D}(t)$ converges towards 0, with a formal computation, the rate of convergence is expected to be either $1/t$ or $1/\sqrt{t}$ depending on $\widehat{D}_{0}$. The convergence towards the stable anisotropic equilibrium is exponentially fast.

6 Macroscopic limit for the stable equilibria

We now go back to the spatially inhomogeneous model. We want to investigate (at least formally) the hydrodynamic models derived from the BGK equation (4). To do so, we introduce the scaling $t^{\prime}=\varepsilon t$ and $x^{\prime}=\varepsilon x$ for $\varepsilon>0$ and we define $f^{\varepsilon}(t^{\prime},x^{\prime},A):=f(t,x,A)$. After this change of variables in (4) and dropping the primes, we see that $f^{\varepsilon}$ satisfies the following equation:

[TABLE]

We want to investigate the macroscopic limit $\varepsilon\to 0$ with the assumption that $f^{\varepsilon}$ converges towards a stable equilibrium. Thanks to the results of the last section, we will assume that $f^{\varepsilon}\to\rho M_{\alpha\Lambda}$ where $\rho=\rho(t,x)$, $\alpha=\rho c_{1}(\alpha)$, $\alpha\in\mathbb{R}_{+}$ and $\Lambda=\Lambda(t,x)\in\mathscr{M}_{3}(\mathbb{R})$ (with a notion of convergence as strong as needed). Since the equilibrium is assumed to be stable there are two cases: either $\Lambda\in SO_{3}(\mathbb{R})$ (and therefore $\rho>\rho^{*}$) or $\Lambda=0$ that is to say $f^{\varepsilon}$ is uniform in the body-orientation variable and converges towards $\rho=\rho(t,x)$ (and therefore $\rho<\rho_{c}$). For a given time $t\in\mathbb{R}_{+}$, we will say that $x\in\mathbb{R}^{3}$ belongs to a disordered region when $\Lambda(t,x)=0$. Otherwise, when $\Lambda(t,x)\in SO_{3}(\mathbb{R})$, we will say that $x\in\mathbb{R}^{3}$ belongs to an ordered region.

The purpose of the two next sections is to write at least formally the hydrodynamic equations satisfied by $\rho=\rho(t,x)$ and $\Lambda=\Lambda(t,x)$. First notice that integrating (36) over $SO_{3}(\mathbb{R})$ leads to the conservation law:

[TABLE]

where $\rho^{\varepsilon}\equiv\rho_{f^{\varepsilon}}$ and

[TABLE]

The macroscopic model then depends on the region considered.

In a disordered region, $j[f^{\varepsilon}]\to 0$ and assuming that the convergence is sufficiently strong, we get that $\partial_{t}\rho=0$. To obtain more information we will look at the next order in the Chapman-Enskog expansion (Section 6.1). 2. 2.

In an ordered region, $j[f^{\varepsilon}]\to\alpha\Lambda e_{1}$ where $\alpha=\rho c_{1}(\alpha)$ and therefore, assuming that the convergence is strong enough:

[TABLE]

However due to the lack of conserved quantities, we will need specific tools to write an equation satisfied by $\Lambda=\Lambda(t,x)$ in order to obtain a closed system of equations on $(\rho,\Lambda)$. This is the purpose of Section 6.2.

6.1 Diffusion model in a disordered region

We consider a region where $f^{\varepsilon}$ converges as $\varepsilon\to 0$ to a density $\rho(t,x)$ uniform in the body-attitude variable. The following proposition gives the diffusion model obtained by looking at the next order in the Chapman-Enskog expansion.

Proposition 6.1** (Formal).**

In a disordered region, the density $\rho^{\varepsilon}$ satisfies formally at first order the following diffusion equation:

[TABLE]

Proof.

We follow the same calculations as in [13] : we write $f^{\varepsilon}=\rho^{\varepsilon}+\varepsilon f_{1}^{\varepsilon}$ (where $f^{\varepsilon}_{1}$ is defined by this relation) and notice that:

[TABLE]

Inserting this in (36), multiplying by $A$ and integrating over $SO_{3}(\mathbb{R})$ leads to:

[TABLE]

Using Lemma 3.3, it holds that

[TABLE]

To compute the second term, we note that $Ae_{1}\cdot\nabla_{x}\rho^{\varepsilon}=2A\cdot R_{\rho^{\varepsilon}}$ where $R_{\rho^{\varepsilon}}$ is the matrix, the first column of which is equal to $\nabla_{x}\rho^{\varepsilon}$ and the others are equal to zero. Using Lemma 3.3 we obtain

[TABLE]

By multiplying (39) by $e_{1}$, it follows from (40) and (41) that :

[TABLE]

which gives the result by inserting this in (37). ∎

Remark 6.1**.**

This analysis does not depend on the dimension. In $SO_{n}(\mathbb{R})$ the same formal result holds:

[TABLE]

6.2 Self-organised hydrodynamics in an ordered region

In the following, for a given density $\rho\in\mathbb{R}_{+}$, $\alpha(\rho)$ denotes the maximal nonnegative root of $\alpha=\rho c_{1}(\alpha)$. We are going to prove the following theorem.

Theorem 9** (Formal).**

We suppose that $f^{\varepsilon}\to\rho(x,t)M_{J(x,t)}$ (as strongly as necessary) as $\varepsilon\to 0$ where $J(x,t)=\alpha(\rho(x,t))\Lambda(x,t)$ and $\Lambda(x,t)\in SO_{3}(\mathbb{R})$. Then $\rho$ and $\Lambda$ satisfy the following system of partial differential equations:

[TABLE]

where $\tilde{c}_{2}$, $\tilde{c_{3}}$, $c_{4}$ are functions of $\rho$ to be defined later and $\delta$ and $\mathbf{r}$ are the “divergence” and “rotational” operators defined in [14] : if $\Lambda(x)=\exp([\mathbf{b}(x)]_{\times})\Lambda(x_{0})$ with $\mathbf{b}$ smooth around $x_{0}$ and $\mathbf{b}(x_{0})=0$, then

[TABLE]

where $\nabla_{x}\times$ is the curl operator.

The first equation (42a) is the conservation law (37) and the goal is to obtain the equation (42b) for $\Lambda=\Lambda(t,x)$. However, here and contrary to the classical gas dynamics, the total momentum is not conserved:

[TABLE]

and we therefore cannot deduce easily a closed system of equations. This lack of conserved quantities is specific to self-propelled particle models such as the Vicsek model. The main tool to tackle the problem will be the Generalised Collision Invariants method introduced in [18] for the study of the hydrodynamic limit of the continuum Vicsek model. Its precise setting in the context of our body-attitude model is detailed Section 6.2.1. The formal proof of Theorem 9 can be found in Section 6.2.2.

6.2.1 Generalised collision invariants

To obtain an equation on $\Lambda$, the main tool are the Generalised Collisional Invariants (GCI) first introduced in [18]. For a given $J\in\mathscr{M}_{3}(\mathbb{R})$, we define first the linear collision operator:

[TABLE]

so that $Q_{BGK}(f)=\mathcal{L}_{J_{f}}(f)$. Let $J\in\mathscr{M}_{3}(\mathbb{R})$ with $\det J>0$ and let $\Lambda\leavevmode\nobreak\ :=PD(J)\in SO_{3}(\mathbb{R})$ be the orthogonal part of its polar decomposition. The set of GCI associated to $J$ is defined as:

[TABLE]

The condition $\psi\in\mathscr{C}_{J}$ is equivalent to:

[TABLE]

Therefore, following the ideas of the proof of [14, Proposition 4.3], we have:

[TABLE]

that is to say:

[TABLE]

or equivalently since $B\in T_{\Lambda}$ means that there exists $P\in\mathscr{A}_{3}(\mathbb{R})$ such that $B=\Lambda P$ :

[TABLE]

where

[TABLE]

Now for any $P\in\mathscr{A}_{3}(\mathbb{R})$, denoting $\Lambda_{f^{\varepsilon}}\equiv\Lambda^{\varepsilon}=PD(J_{f^{\varepsilon}})$, we get by multiplying the equation (36) by $\psi_{P}^{\Lambda^{\varepsilon}}$ :

[TABLE]

The right-hand side vanishes by Definition (43) of the GCI. Note that if $f^{\varepsilon}\to\rho M_{J}$ with $\det(J)>0$, we have also $\det(J_{f^{\varepsilon}})>0$ for $\varepsilon$ sufficiently small and $\Lambda^{\varepsilon}\in SO_{3}(\mathbb{R})$.

6.2.2 Hydrodynamic limit (formal proof of Theorem 9)

Taking formally the limit $\varepsilon\to 0$ in (44) and since it is true for all $P\in\mathscr{A}_{3}(\mathbb{R})$, we obtain:

[TABLE]

We have:

[TABLE]

where $\alpha^{\prime}$ denotes the derivative of $\alpha(\rho)$ with respect to $\rho$ and similarly for $\partial_{i}(M_{\alpha\Lambda})$. With this we compute the term:

[TABLE]

where we have used that:

[TABLE]

and

[TABLE]

Most of the terms that appear in $X$ are computed in [14]. Precisely,

[TABLE]

where $X_{1}$, $X_{2}$, $X_{3}$ and $X_{4}$ are computed in [14] :

[TABLE]

where the coefficients

[TABLE]

and the matrix

[TABLE]

are the same as in [14]. The matrix $\mathcal{D}_{x}(\Lambda)\in\mathscr{M}_{3}(\mathbb{R})$ is defined as the unique matrix such that for all $\mathbf{w}\in\mathbb{R}^{3}$, and smooth functions $\Lambda:\mathbb{R}^{3}\to SO_{3}(\mathbb{R})$,

[TABLE]

(see [14, Section 4.5]). Note that $C_{3}=C_{2}$ since the noise and alignement parameters which were denoted by $\nu$ and $d$ in [14] have been taken equal to 1 here. Note also that these coefficients are functions of $\rho$ (through $\alpha$ only). The term $Y$ is an additional term which appears here due to the presence of the parameter $\alpha=\alpha(\rho)$ which is a function of $\rho$. It depends also on the derivative $\alpha^{\prime}$ of $\alpha$ :

[TABLE]

All the terms that involve the time derivative of $\rho$ are equal to zero since $\partial_{t}\rho$ does not depend on $A$ and with the change of variable $A^{\prime}=\Lambda A^{T}\Lambda$ which has unit jacobian, it holds that $A^{\prime}\cdot\Lambda=A\cdot\Lambda$ and therefore

[TABLE]

We thus have:

[TABLE]

where

[TABLE]

and

[TABLE]

With the change of variable $A\mapsto\Lambda^{T}A$ these terms become

[TABLE]

[TABLE]

and they can be computed using the same techniques as in [14] or lemma 3.4 in the appendix. More precisely, we can write

[TABLE]

where $R_{1,\rho}$ is the matrix, the first column of which is equal to $2\Lambda^{T}\nabla_{x}\rho$ and the others are all equal to zero. It satisfies:

[TABLE]

By the change of variable $B\mapsto B^{T}$, we have

[TABLE]

This integral is of the form (46) where

[TABLE]

is invariant by transposition and conjugation and

[TABLE]

is a skew-symmetric matrix. From (48) and (49) we get

[TABLE]

with

[TABLE]

The first term $Y_{1}$ can therefore be written:

[TABLE]

and using Rodrigues’ formula we obtain:

[TABLE]

where $\{\cdot\}_{\alpha}$ has been defined in Proposition 4.1. Similarly the second term $Y_{2}$ can be written:

[TABLE]

where the coefficient $C_{3}$ is the same as in [14] :

[TABLE]

Finally, we obtain:

[TABLE]

Putting all the terms together, we can conclude as in [14]. First we notice that:

[TABLE]

Therefore we obtain by multiplying (45) by $\Lambda$ and dividing by $\alpha C_{2}$ :

[TABLE]

where the coefficients

[TABLE]

are respectively equal to the coefficients $c_{2}$ and $c_{4}$ in [14] and the coefficient $c_{3}$ in [14] (which is equal to $1$) becomes:

[TABLE]

7 Conclusion

In this work, we have presented a new BGK model of body-attitude coordination where agents are described by a rotation matrix. Starting from the kinetic level (a space homogeneous BGK equation) we have drawn a parallel between our Vicsek-type model and the models of nematic alignment of polymers. We then have deduced the equilibria of the system and have shown a phase transition phenomenon triggered by the density of agents. Thanks to a gradient-flow structure specific to the BGK equation we have been able to describe the asymptotic behaviour of the system. Finally, we have derived the macroscopic models (SOHB) in the spatially inhomogeneous case.

On the modelling side, a rigorous mean-field limit which leads to the BGK equation is currently under study and will be the object of future work. However, many other questions remain open. At the kinetic level, our study relies on the dimension 3 and it would be interesting to extend the ideas developed here in $SO_{n}(\mathbb{R})$, $n\geq 3$, for example by drawing new parallels with higher dimensional polymers models or other similar models [28, 27]. In addition, the mathematical and numerical analyses of the macroscopic SOHB model are still in progress.

The BGK model studied here is a step towards the full description of the models of collective behaviour depicted in Figure 1. The tools and ideas that we have presented here may help to analyse other models of body-attitude coordination such as the non-normalised Fokker-Planck model in $SO_{3}(\mathbb{R})$. Other models which take into account curvature control in addition to body-orientation, in the spirit of [19], could also be considered.

Appendix A Quaternions and rotations

These appendix is devoted to the proof of Proposition 4.3. We also give additional results about quaternions. The following lemma gives a link between quaternions and the theory of $Q$-tensors.

Lemma A.1**.**

Let $\SS_{4}^{0}(\mathbb{R})$ be the space of symmetric $4\times 4$ trace free matrices. If $Q\in\SS_{4}^{0}(\mathbb{R})$ has two eigenvalues with eigenspaces of dimensions 1 and 3, then $Q$ can be written

[TABLE]

for a given unit quaternion $q$ seen as a vector of $\mathbb{R}^{4}$. A matrix of this form is called a uniaxial $Q$-tensor. When $\alpha=1$ we will say that $Q$ is a normalised uniaxial $Q$-tensor.

Proof.

Let $Q\in\SS_{4}^{0}(\mathbb{R})$ such that $Q$ has two eigenvalues with eigenspaces of dimensions 1 and 3. By the spectral theorem, there exists $P\in\mathcal{O}_{3}(\mathbb{R})$ such that for a given $\alpha>0$ :

[TABLE]

and the result follows by taking $q$ equals to the first column of $P$. ∎

Proof of Proposition 4.3.

The group isomorphism $\Phi$ is explicitly computed in [52]. In particular, let $q\in\mathbb{S}^{3}/\pm 1$. The matrix $A=\Phi(q)$ is defined for all purely imaginary quaternion $u\in\mathbb{H}$ by $A[u]=[quq^{*}]$ where $(u_{1},u_{2},u_{3})^{T}=:[u]\in\mathbb{R}^{3}$ is the vector associated to $u=u_{1}i+u_{2}j+u_{3}k$. More explicitly, if $q=x+iy+zj+tk$, then

[TABLE]

Note that we have identified $q\in\mathbb{H}$ and its equivalence class in $\mathbb{S}^{3}/\pm 1$ and that $\Phi$ is well defined since only quadratic expressions are involved. The fact that this group isomorphism is an isometry follows from [17, Proposition A.3] and [14, Lemma 4.2]. 2. 2.

The expression

[TABLE]

is a quadratic form for $q$. We take $Q$ the matrix associated to this quadratic form. For $J=(J_{ij})_{i,j}$, using the explicit form of $A=\Phi(q)$ with $q=x+yi+zj+tk$ we obtain:

[TABLE]

This is an isomorphism since $\dim\mathscr{M}_{3}(\mathbb{R})=\dim\SS_{4}^{0}(\mathbb{R})$ and $J$ can be obtained from $Q$ similarly. Moreover, since the bilinear matrix associated to a quadratic form is uniquely defined, if $Q\in\SS_{4}^{0}(\mathbb{R})$ is such that $\frac{1}{2}J\cdot\Phi(q)=q\cdot Qq$ for all $q\in\mathbb{S}^{3}/\pm 1$, then $Q=\phi(J)$. 3. 3.

To prove the third point, we note that a unit quaternion can be seen as a rotation in $\mathbb{R}^{3}$ in a more geometrical way ([17, Section 5.1]): for $\theta\in[0,\pi]$ and $\mathbf{n}=(n_{1},n_{2},n_{3})^{T}\in\mathbb{S}^{2}$, let us define the unit quaternion $q$ by

[TABLE]

The unit quaternion $q$ represents the rotation of angle $\theta\in[0,\pi]$ and axis $\mathbf{n}\in\mathbb{S}^{2}$ in the sense that if $R(\theta,\mathbf{n})\in SO_{3}(\mathbb{R})$ denotes the matrix associated to the rotation of angle $\theta\in[0,\pi]$ around the axis $\mathbf{n}\in\mathbb{S}^{2}$ then $\Phi(q)=R(\theta,\mathbf{n})$. Note that $q$ and $-q$ represent the same rotation so $\Phi(q)$ is well defined by identifying $q$ with its equivalence class in $\mathbb{S}^{3}/\pm 1$.

In $\mathbb{R}^{3}$ the composition of two rotations of respective angles and axis $(\theta,\mathbf{n})\in[0,\pi]\times\mathbb{S}^{2}$ and $(\theta^{\prime},\mathbf{n}^{\prime})\in[0,\pi]\times\mathbb{S}^{2}$ is itself a rotation: we have $R(\theta,\mathbf{n})R(\theta^{\prime},\mathbf{n}^{\prime})=R(\hat{\theta},\hat{\mathbf{n}})$ where the angle $\hat{\theta}\in[0,\pi]$ is defined by

[TABLE]

Note that $\cos(\hat{\theta}/{2})=q\cdot\bar{q}^{\prime}$ where $q$ and $q^{\prime}$ are the associated unit quaternions seen as vectors of dimension 4. In particular the dot product of two rotations matrices is

[TABLE]

where

[TABLE]

Besides, for the quaternions $q$ and $q^{\prime}$ respectively associated to the rotations $R(\theta,\mathbf{n})$ and $R(\theta^{\prime},\mathbf{n}^{\prime})$, we have:

[TABLE]

where $Q$ is the normalised uniaxial $Q$-tensor:

[TABLE]

Finally $\frac{1}{2}R(\theta,\mathbf{n})\cdot R(\theta^{\prime},\mathbf{n}^{\prime})=q^{\prime}\cdot Qq^{\prime}$ and we obtain thanks to the previous point:

[TABLE]

that is to say: if $J\in SO_{3}(\mathbb{R})$ then $\phi(J)$ is a normalised uniaxial $Q$-tensor. 4. 4.

If $D=\operatorname{diag}(d_{1},d_{2},d_{3})$ then using the explicit form of $\phi$ given in the second point:

[TABLE]

and if $Q=\operatorname{diag}(s_{1},s_{2},s_{3},s_{4})$ with $s_{1}+s_{2}+s_{3}+s_{4}=0$ then

[TABLE]

∎

Appendix B More about $SO_{3}(\mathbb{R})$ and $SO_{n}(\mathbb{R})$

Lemma B.1**.**

For all $n\geq 3$ :

[TABLE]

Proof.

First we prove that the diagonal matrices form a subset of $\operatorname{Span}(SO_{n}(\mathbb{R}))$: it is enough to show that

[TABLE]

(the other diagonal matrices with only one nonzero coefficient can be obtained in a similar way). When $n$ is odd:

[TABLE]

and both matrices in the sum are in $SO_{n}(\mathbb{R})$. When $n\geq 4$ is even,

[TABLE]

thus:

[TABLE]

and similarly:

[TABLE]

The SSVD (Definition 3.2) gives the result for any matrix. ∎

Corollary B.1**.**

For $n\geq 3$, a matrix that commutes with any matrix of $SO_{n}(\mathbb{R})$ is of the form $\lambda I_{n}$.

We can now prove Lemma 3.3 and its generalisation Lemma 3.4.

Proof (of Lemma 3.3).

The linear map $\Phi\leavevmode\nobreak\ :\mathscr{M}_{n}(\mathbb{R})\to\mathscr{M}_{n}(\mathbb{R})$ defined by

[TABLE]

satisfies for all $P\in SO_{n}(\mathbb{R})$,

[TABLE]

which by corollary B.1 means that

[TABLE]

Therefore, $\Phi(P)=\lambda P$ for any $P\in SO_{n}(\mathbb{R})$ and the same is true for any matrix $J\in\mathscr{M}_{n}(\mathbb{R})$ by lemma B.1. To compute $\lambda$, notice that for the matrix $e_{i}\otimes e_{j}$ :

[TABLE]

Summing this equality for all $i,j$ gives:

[TABLE]

and $\lambda=1/2n$. ∎

Proof (of Lemma 3.4).

Let us define the linear map $\psi\leavevmode\nobreak\ :\mathscr{M}_{n}(\mathbb{R})\longrightarrow\mathscr{M}_{n}(\mathbb{R})$ by

[TABLE]

We first note that $\psi$ is self-adjoint for the dot product $A\cdot B=\operatorname{Tr}(A^{T}B)$: for any $K\in\mathscr{M}_{n}(R)$,

[TABLE] 2. 2.

We prove that $\operatorname{Span}(I_{n})$ is a stable supspace for $\psi$ : for any $P\in SO_{n}(\mathbb{R})$ we have:

[TABLE]

and we conclude with corollary B.1 that $\psi(I_{n})=\alpha I_{n}$ with:

[TABLE] 3. 3.

Since $\psi$ is a self adjoint operator, the orthogonal subspace $\operatorname{Span}(I_{n})^{\perp}$ is also a stable subspace. Moreover, using the change of variable $A\mapsto A^{T}$, we see that $\psi(J^{T})=\psi(J)^{T}$ and we have the decomposition:

[TABLE]

where $\SS_{n}^{0}(\mathbb{R})$ and $\mathscr{A}_{n}(\mathbb{R})$ are respectively the subspace of trace free symmetric matrices and the subspace of skew-symmetric matrices. They are both stable subspaces. 4. 4.

We prove now that $\psi\leavevmode\nobreak\ :\SS_{n}^{0}(\mathbb{R})\to\SS_{n}^{0}(\mathbb{R})$ is a uniform scaling. By the spectral theorem, every matrix $J\in\SS_{n}^{0}(\mathbb{R})$ can be written

[TABLE]

where $P\in SO_{n}(\mathbb{R})$ and $D$ is diagonal. Since

[TABLE]

it is enough to prove that there exists $\lambda\in\mathbb{R}$ such that for all diagonal matrices $D$, $\psi(D)=\lambda D$. For $D=\operatorname{diag}(d_{1},\ldots,d_{n})$, we have:

[TABLE]

Now, if $i\neq j$, let $D^{ik}\in SO_{n}(\mathbb{R})$ be the diagonal matrix such that all the coefficients are equal to 1 except $D^{ik}_{ii}$ and $D^{ik}_{kk}$ which are equal to $-1$. Then the change of variable $A\mapsto D^{ik}A(D^{ik})^{T}$ gives

[TABLE]

which proves that $\psi(D)$ is diagonal and the $i$-th coefficient of $\psi(D)$ is:

[TABLE]

Using the fact that $\operatorname{Tr}(D)=0$ and that all the $\int_{SO_{n}(\mathbb{R})}a_{ii}\,a_{kk}\,g(A)\,dA$ are equal for $i\neq k$ (by using conjugation by the matrices $P^{ij}$, $i,j\neq k$, see Definition 3.1), we obtain:

[TABLE]

and we conclude that for all $J\in\SS_{n}^{0}(\mathbb{R})$ :

[TABLE]

with

[TABLE] 5. 5.

We prove similarly that $\psi\leavevmode\nobreak\ :\mathscr{A}_{n}(\mathbb{R})\to\mathscr{A}_{n}(\mathbb{R})$ is a uniform scaling. Every $J\in\mathscr{A}_{n}(\mathbb{R})$ can be written

[TABLE]

where

[TABLE]

is a block diagonal matrix with blocks

[TABLE]

so that,

[TABLE]

and

[TABLE]

When $n\geq 3$, by using conjugation by the matrices $D^{2j-1,2j}$ for $j\in\{1,\ldots,\left\lfloor n/2\right\rfloor\}$ (Definition 3.1) we can see that for each $k\in\{1,\ldots,p\}$, the matrix

[TABLE]

is of the form (47). Moreover, when $n\geq 5$, by using conjugation by matrices $D^{2j-1,2\ell-1}$ for $j\neq\ell$, $j,\ell\in\{1,\ldots,\left\lfloor n/2\right\rfloor\}$ and $j,\ell\neq k$, we can see that all the diagonal blocks of $M_{k}$ are equal to zero except the one in position $2k-1$. When $n=3$ there is only one block in position 1 so the result holds but when $n=4$ such $j$ and $\ell$ do not exist. In conclusion, when $n\neq 4$, $\psi(C)$ is of the form (47) and each diagonal block $C_{2k-1}^{\prime}$ of $\psi(C)$ is written

[TABLE]

with

[TABLE]

where this equality follows by using conjugation by “block-permutation matrices”:

[TABLE]

where the first zero on the diagonal is in position $2k-1$. Therefore,

[TABLE]

with

[TABLE] 6. 6.

Finally, for $J\in\operatorname{Span}(I_{n})^{\perp}$, writing

[TABLE]

we have:

[TABLE]

with

[TABLE]

and

[TABLE] 7. 7.

In conclusion, writing the decomposition

[TABLE]

where $K\in\operatorname{Span}(I_{n})^{\perp}$, we obtain

[TABLE]

with

[TABLE]

and

[TABLE]

and

[TABLE]

And there are of course many other ways to write the coefficients $a$, $b$ and $c$.

∎

Remark B.1**.**

In dimension 4, the result still holds for symmetric matrices. For general matrices, the result can be proved in particular cases, for instance when $g$ is a function of the trace, by using an explicit parametrisation of $SO_{4}(\mathbb{R})$ such as the 4-dimensional version of (10).

Acknowledgments

PD acknowledges support by the Engineering and Physical Sciences Research Council (EPSRC) under grants no. EP/M006883/1 and EP/P013651/1, by the Royal Society and the Wolfson Foundation through a Royal Society Wolfson Research Merit Award no. WM130048 and by the National Science Foundation (NSF) under grant no. RNMS11-07444 (KI-Net). PD is on leave from CNRS, Institut de Mathématiques de Toulouse, France.

A. D. acknowledges the hospitality of CEREMADE, Université Paris-Dauphine and the University of Sussex where part of this research was carried out. A.D. was supported for this work by the École Normale Supérieure de Rennes through an Erasmus+ grant. This work was initiated at the École Normale Supérieure de Lyon as part of a Master degree delivered to A.D.

A.F. acknowledges support from the EFI project ANR-17-CE40-0030 and the Kibord project ANR-13-BS01-0004 of the French National Research Agency (ANR), from the project Défi S2C3 POSBIO of the interdisciplinary mission of CNRS, and the project SMS co-funded by CNRS and the Royal Society.

SMA is supported by the Vienna Science and Technology Fund (WWTF) with a Vienna Research Groups for Young Investigators, grant VRG17-014.

Data statement: no new data were collected in the course of this research.

Conflict of interest: the authors have no conflict of interest to declare.