Rinaldo M Colombo

Mauro Garavello

Magali Lécureux-Mercier

A Class of Non-Local Models for Pedestrian Traffic

Keywords: 2000 Mathematics Subject Classification: 35L65, 90B20 Crowd Dynamics, Macroscopic Pedestrian Model, Non-Local Conservation Laws

We present a new class of macroscopic models for pedestrian flows. Each individual is assumed to move towards a fixed target, deviating from the best path according to the instantaneous crowd distribution. The resulting equation is a conservation law with a nonlocal flux. Each equation in this class generates a Lipschitz semigroup of solutions and is stable with respect to the functions and parameters defining it. Moreover, key qualitative properties such as the boundedness of the crowd density are proved. Specific models are presented and their qualitative properties are shown through numerical integrations.

Introduction

From a macroscopic point of view, a moving crowd is described by its density ρ = ρ(t, x), so that for any subset A of the plane, the quantity A ρ(t, x) dx is the total number of individuals in A at time t. In standard situations, the number of individuals is constant, so that conservation laws of the type ∂ t ρ + div x (ρ v) = 0 are the natural tool for the description of crowd dynamics. A key issue is the choice of the speed v, which should describe not only the target of the pedestrians and the modulus of their speed, but also their attitude to adapt their path choice to the crowd density they estimate to find along this path.

Our starting point is the following Cauchy problem for the conservation law

   ∂ t ρ + div ρ v(ρ) ν(x) + I(ρ) = 0 , ρ(0, x) = ρ 0 (x) . (1.1)
The scalar function ρ → v(ρ) describes the modulus of the pedestrians' speed, independently from geometrical considerations. In other words, an individual at time t and position x ∈ R N moves at the speed v ρ(t, x) that depends on the density ρ(t, x) evaluated at the same time t and position x. Given that the density is ρ, the vector ν(x) + I(ρ) describes the direction that the individual located at x follows and has norm (approximately) 1. More precisely, the individual at position x and time t is assumed to move in the direction ν(x) + I ρ(t) (x).

In situations like the evacuation of a closed space Ω ⊂ R N , it is natural to assume that the first choice of each pedestrian is to follow a path optimal with respect to the visible geometry, for instance the geodesic. As soon as walls or obstacles are relevant, it is necessary to take into consideration the discomfort felt by pedestrians walking along walls or too near to corners, see for instance [START_REF] Helbing | Pedestrian, crowd and evacuation dynamics[END_REF][START_REF] Hoogendoorn | Pedestrian route-choice and activity scheduling theory and model[END_REF] and the references therein.

The vector I ρ(t) (x) describes the deviation from the direction ν(x) due to the density distribution ρ(t) at time t. Hence, the operator I is in general nonlocal, so that I ρ(t) (x) depends on all the values of the density ρ(t) at time t in a neighborhood of x. More formally, it depends on all the function ρ(t) ∈ L 1 (R N ; [0, R]) and not only on the value ρ(t, x) ∈ [0, R].

The case in which I = 0 is equivalent to assume that the paths followed by the individuals are chosen a priori, independently from the dynamics of the crowd.

Here we present two specific choices that fit in (1.1). A first criterion assumes that each individual aims at avoiding high crowd densities. Fix a mollifier η. Then, the convolution (ρ * η) is an average of the crowd density around x. This leads to the natural choice

I(ρ) = -ε ∇(ρ * η) 1 + ∇(ρ * η) 2 , (1.2)
related to [START_REF] Bressan | models of pedestrian flow[END_REF], which states that individuals deviate from the optimal path trying to avoid entering regions with higher densities. Through numerical integrations, below we provide examples of solutions to (1.1)-(1.2). They show the interesting phenomenon of pattern formation. In the case of a crowd walking along a corridor, coherently with the experimental observation described in the literature, see for instance [START_REF] Helbing | Simulation of pedestrian crowds in normal and evacuation situations, volume Pedestrian and Evacuation Dynamics[END_REF][START_REF] Helbing | Quantitative agent-based modeling of human interactions in space and time[END_REF][START_REF] Hoogendoorn | Simulation of pedestrian flows by optimal control and differential games[END_REF][START_REF] Piccoli | Pedestrian flows in bounded domains with obstacles[END_REF], the solution to (1.1)-(1.2) self-organizes into lanes. The width of these lanes depends on the size of the support of the averaging kernel η. This feature is stable with respect to strong variations in the initial datum and also in the geometry. Indeed, we have lane formation also in the case of the evacuation of a room, when the crowd density sharply increases in front of the door. Section 4.1 is devoted to this property.

From the analytical point of view, we note that the convolution term in (1.2) seems not sufficient to regularize solutions. Indeed, the present analytical framework is devised to consider solutions in L 1 ∩ BV. Both in the case of a crush in front of an exit (Section 4.2) and in the specific example in Section 4.3, numerical simulations highlight that the space gradient of ρ increases dramatically.

According to (1.2), pedestrians evaluate the crowd density all around their position. When restrictions on the angle of vision are relevant, the following choice is reasonable:

I(ρ) = -ε ∇ R N ρ(y) η(x -y) ϕ (y -x) • g(x) dy 1 + ∇ R N ρ(y) η(x -y) ϕ (y -x) • g(x) dy 2 . (1.3)
Here, η is a fixed mollifier as above and the smooth function ϕ weights the deviation from the preferred direction g(x).

We note that the constructions in [START_REF] Hughes | A continuum theory for the flow of pedestrians[END_REF][START_REF] Hughes | The ow of human crowds[END_REF] and [START_REF] Di Francesco | On the hughes' model for pedestrian flow: The one-dimensional case[END_REF] fit in the present setting. Indeed, there the following choices were considered:

∂ t ρ -div ρ f 2 (ρ) ∇ϕ = 0 |∇ϕ| = 1/f (ρ) and ∂ t ρ -div ρ f 2 (ρ) ∇ϕ = 0 -ε∆ϕ + |∇ϕ| 2 = 1/ f (ρ) + ε 2 .
The former admits an immediate interpretation: the direction ∇ϕ of the speed v is chosen along the solutions to the eikonal equations, i.e. all pedestrians follow the shortest path, weighing at every instant the length of paths with the amount of people that are moving along it. In the former case, ν(x) = 0 and I(ρ) is the gradient of the solution to the eikonal equation with 1/f (ρ) in the right hand side, while in the latter case ν(x) = 0 and I(ρ) is the gradient of the solution to an elliptic partial differential equation. The model introduced in [START_REF] Cristiani | Multiscale modeling of granular flows with application to crowd dynamics[END_REF][START_REF] Piccoli | Time-evolving measures and macroscopic modeling of pedestrian flow[END_REF] relies on this measure valued conservation law:

∂ t µ + div (µ v) = 0 where v = ν(x) + R 2 f (|x -y|) ϕ (y -x) • ν(x) x -y x -y dµ t (y) .
Here the unknown is a map µ : [0, T ] → M(R N ; R) assigning at every time t a positive measure µ(t) which substitutes the crowd density, in the sense that the amount of people that at time t are in A is µ(t) (A). Contrary to this model, the present framework is the space of L 1 densities and physical a priori L ∞ bounds on the solutions to (1.1) are rigorously proved, preventing any focusing effect as well as the rise of any Dirac delta. Moreover, below we prove global in time existence of solutions, their continuous dependence from the initial data and their stability with respect to variations in the speed law v.

In [START_REF] Coscia | First-order macroscopic modelling of human crowd dynamics[END_REF], the geometric part ν of the speed is chosen a priori, while its modulus depend on the density as well as on the gradient of the density:

∂ t ρ + div ρ ϕ(ρ, ∇ρ) ν(x) = 0 .
On the contrary, the model introduced in [START_REF] Colombo | Control of the continuity equation with a non local flow[END_REF] postulates a nonlocal dependence of the speed from the density:

∂ t ρ + div ρ v(ρ * η) ν(x) = 0 ,
which amounts to assume that pedestrians choose their behavior according to evaluations of an average of the density around their position, rather than according to the density at their place. The following second order model was presented in [START_REF] Bellomo | On the modelling crowd dynamics from scaling to hyperbolic macroscopic models[END_REF][START_REF] Dogbé | On the Cauchy problem for macroscopic model of pedestrian flows[END_REF] and does not fit in (1.1):

∂ t ρ + div (ρv) = 0 ∂ t v + (v • ∇)v = F (ρ, ∇ρ, v) ,
we refer to the review [START_REF] Bellomo | On the modelling of traffic and crowds -a survey of models, speculations, and perspectives[END_REF] for further details. To underline the variety of analytical techniques with which crowd dynamics has been tackled, we recall the further approaches: optimal transport in [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF], the mean field limit in [START_REF] Dogbé | Modeling crowd dynamics by the mean-field limit approach[END_REF], the functional analytic one in [START_REF] Maury | A macroscopic crowd motion model of the gradient-flow type[END_REF][START_REF] Maury | Handling congestion in crowd motion modeling[END_REF] and the nonclassical shocks used in [START_REF] Chalons | Numerical approximation of a macroscopic model of pedestrian flows[END_REF][START_REF] Colombo | Pedestrian flows and non-classical shocks[END_REF][START_REF] Colombo | Existence of nonclassical solutions in a pedestrian flow model[END_REF]. Throughout this work, we set all statements in all R N . In Section 3 we show how the present framework is able to take into consideration the presence of various constraints. In the case of pedestrian dynamics, for instance, this amounts to prove that no individual passes through the walls of a given room, provided the initial datum is assigned inside it. Concerning the dimension, our main applications are referred to the case of crowd dynamics, i.e. N = 2. Nevertheless, from the analytical point of view, considering the case of a general N does not add any difficulty. Furthermore, we believe that the present setting can be reasonably applied also to, say, fishes and birds moving in R 3 , for example in predator-prey like situations as described in [START_REF] Colombo | An analytical framework to describe the interactions between individuals and a continuum[END_REF].

The next section is devoted to the analytical properties of (1.1): well posedness and stability. The general theory is particularized to specific examples in Section 3, where the presence of obstacle (walls) is considered. Sample numerical integrations are provided in Section 4. The final Section 5 collects the analytical proofs. We defer to Appendix A further remarks related to the geometry of the physical domain.

Analytical Results

In the following, N ∈ N \ {0} is the (fixed) space dimension. We denote R + = [0, +∞[ ; the open ball in R N centered at x and with radius r > 0 is B(x, r) and we let W N = π/2 0 (cos ϑ) N dϑ. As usual, we denote

S N -1 = x ∈ R N : x = 1 .
The density ρ can be defined as the fraction of space occupied by pedestrians, so that ρ turns out to be a nondimensional scalar in [0, 1]. Otherwise, it can be useful to think at ρ as measured in, say, individuals/m 2 and varying in [0, R], with R > 0 being a given maximal density, for example R = 8 individuals/m 2 .

Our first step in the study of (1.1) is the formal definition of solution.

Definition 2.1. Fix a positive T and an initial datum

ρ 0 ∈ L 1 (R N ; [0, R]). A function ρ ∈ C 0 [0, T ]; L 1 (R N ; R)
∂ t ρ + div ρ v(ρ) w(t, x) = 0 ρ(0, x) = ρ 0 (x) (2.1)
where w(t, x) = ν(x) + I ρ(t) (x).

In other words, recalling [24, Definition 1], for all k ∈ R,

for all ϕ ∈ C ∞ c (]-∞, T ] × R N ; R + ), T 0 R N |ρ -k| ∂ t ϕ + ρ v(ρ) -k v(k) w(t, x) • ∇ϕ -k v(k) div w(t, x) ϕ sgn(ρ -k) dx dt + R N ρ 0 (x) -k ϕ(0, x) dx 0.
On the functions defining the general model (1.1), we introduce the following hypotheses:

(v) v ∈ C 2 (R; R) is non increasing, v(0) = V and v(R) = 0 for fixed V, R > 0. (ν) ν ∈ (C 2 ∩ W 1,∞ )(R N ; R N ) is such that div ν ∈ (W 1,1 ∩ W 1,∞ )(R N ; R). (I) I ∈ C 0 L 1 (R N ; [0, R]); C 2 (R N ; R N
) satisfies the following estimates:

(I.1) There exists an increasing

C I ∈ L ∞ loc (R + ; R + ) such that, for all r ∈ L 1 (R N ; [0, R]), I(r) W 1,∞ C I ( r L 1 ) and div I(r) L 1 C I ( r L 1 )
.

(I.2)

There exists an increasing

C I ∈ L ∞ loc (R + ; R + ) such that, for all r ∈ L 1 (R N ; [0, R]), ∇ div I(r) L 1 C I ( r L 1 ) . (I.
3) There exists a constant K I such that for all r 1 , r 2 ∈ L 1 (R N ; [0, R]),

I(r 1 ) -I(r 2 ) L ∞ K I • r 1 -r 2 L 1 , I(r 1 ) -I(r 2 ) L 1 + div (I(r 1 ) -I(r 2 )) L 1 K I • r 1 -r 2 L 1 .
Furthermore, throughout we denote by q the map q : [0, R] → R defined by q(ρ) = ρ v(ρ).

As a first justification of these conditions, note that they make Definition 2.1 acceptable.

Lemma 2.2. Fix a positive T . Let (v), (ν) and (I.1) hold. Choose an arbitrary density

r ∈ C 0 [0, T ]; L 1 (R N ; [0, R])
. Then, the Cauchy problem

∂ t ρ + div ρ v(ρ) w(t, x) = 0 ρ(0, x) = ρ 0 (x) with w(t, x) = ν(x) + I r(t) (x) (2.2)
satisfies the assumptions of Kružkov Theorem [START_REF] Kružkov | First order quasilinear equations with several independent variables[END_REF]Theorem 5].

This lemma is proved in Section 5 below. In Section 3 we show that the above assumption (I) allows to comprehend physically reasonable cases.

The next result is devoted to the proof of existence and uniqueness for (1.1). It is obtained through Banach Fixed Point Theorem.

Theorem 2.3. Let (v), (ν) and (I) hold. Choose any ρ 0 ∈ (L 1 ∩ BV)(R N ; [0, R]). Then, there exists a unique weak entropy solution ρ ∈ C 0 R + ; L 1 (R N ; [0, R]) to (1.1). Moreover, ρ satisfies the bounds

ρ(t) L 1 = ρ 0 L 1 , for a.e. t ∈ R + , TV (ρ(t)) TV (ρ 0 ) e kt + te kt N W N q L ∞ ([0,R]) ∇ div ν L 1 + C I ( ρ 0 L 1 )
,

where k = (2N + 1) q ′ L ∞ ([0,R]) ∇ν L ∞ + C I ( ρ 0 L 1 ) .
Using the techniques in [START_REF] Colombo | Control of the continuity equation with a non local flow[END_REF][START_REF] Colombo | Stability and total variation estimates on general scalar balance laws[END_REF], we now obtain the continuous dependence of the solution to (1.1) from the initial datum and its stability with respect to v, ν and I in the natural norms.

Theorem 2.4. Let (v), (ν) and (I) be satisfied by both systems

   ∂ t ρ + div ρ v 1 (ρ) ν 1 (x) + I 1 (ρ) = 0 ρ(0, x) = ρ 0,1 (x)
and

   ∂ t ρ + div ρ v 2 (ρ) ν 2 (x) + I 2 (ρ) = 0 ρ(0, x) = ρ 0,2 (x) with ρ 0,1 , ρ 0,2 ∈ (L 1 ∩ BV)(R N ; [0, R]).
Then, the two solutions ρ 1 and ρ 2 satisfy

ρ 1 (t) -ρ 2 (t) L 1 C(t) ρ 0,1 -ρ 0,2 L 1 + q 1 -q 2 W 1,∞ + ν 1 -ν 2 L ∞ + div (ν 1 -ν 2 ) L 1 + d (I 1 , I 2 )
where

d(I 1 , I 2 ) = sup I 1 (ρ) -I 2 (ρ) L ∞ + div I 1 (ρ) -div I 2 (ρ) L 1 : ρ ∈ L 1 (R N ; [0, R]) , the map C ∈ C 0 (R + ; R + ) vanishes at t = 0 and depends on TV (ρ 0,1 ), ρ 0,1 L 1 , ν 1 L ∞ , div ν 1 W 1,1 , q 1 W 1,∞ , q 2 W 1,∞ .
Thanks to these stability results, several control problems can be proved to admit a solution through a direct application of Weierstraß Theorem. A possible standard application could be the minimization of the evacuation time from a given room. Describing the actions of a controller able to determine the initial pedestrians' distribution and/or their preferred paths, one is lead to an optimal control problem with the initial datum ρ 0 and the vector field ν as control parameters, for instance. Without any loss of generality, the compact sets on which the optimization is made can be ρ

∈ (L 1 ∩ BV)(R N ; [0, 1]) : TV (ρ) < M for ρ 0 and ν ∈ C 2 R N ; B(0, 1) : ∇ 3 ν
M for ν, where M > 0 is arbitrary.

A different problem solved by the same analytical techniques is that of the dynamic parameter estimation. Once real data are available, one is left with the problem of determining the various parameters entering ν, v or I. Theorem 2.4 ensures the existence of the parameters that allows a best agreement between the solutions to (1.1) and the data.

The Models

This section is devoted to the study of specific cases of (1.1). Aiming at real applications, it is necessary to take into consideration the various constraints that are present on the movement of pedestrians. Therefore, we introduce the subset Ω of R N which characterizes the region reachable to any individual. The boundary ∂Ω consists of walls that can not be crossed by any individual. On the set Ω we require that (Ω) Ω ⊆ R N is the closure of a non empty connected open set. If ∂Ω is not empty, there exists a positive r Ω such that the function

d ∂Ω : B(∂Ω, r Ω ) ∩ Ω → R + x → inf d(x, w) : w ∈ ∂Ω is of class C 2 B(∂Ω, r Ω ) ∩ Ω; R + .
Note that we do not require Ω to be bounded. This assumptions allows us to introduce the inward normal

n(x) = ∇d ∂Ω (x) (3.1)
on all the strip B(∂Ω, r Ω ) ∩ Ω. Moreover, n(x) = 1.

For the present models to be acceptable, it is mandatory that no individual enters any wall, provided the initial datum is supported inside the physically admissible space. Analytically, this is described by the following invariance property:

(P) The model (1.1) is invariant with respect to Ω if spt ρ 0 ⊂ Ω =⇒ spt ρ(t) ⊂ Ω for all t 0 .
Below, we show that theorems 2.3 and 2.4 can be applied to reasonable crowd dynamics models, in the sense that the presence of a physically admissible set Ω is also considered and its invariance in the sense of (P) is proved. More precisely, we show below the following sufficient condition for invariance.

L 1 (R N ; [0, R]) with spt ρ ⊆ Ω ν(x) + I(ρ) (x) • n(x) 0 (3.2)
then, property (P) holds.

The Model

(1.1)-(1.
2)

The starting point for (1.1)-(1.2) is provided by the degenerate parabolic model

∂ t ρ + div   ρ v(ρ) ν(x) -ε ∇ψ(ρ) 1 + ∇ψ(ρ) 2    = 0 introduced in [4], which fits in (1.1) with I(ρ) = -ε ∇ψ(ρ) 1 + ∇ψ(ρ) 2 , motivated
by the desire of each individual to avoid entering regions occupied by a high crowd density.

Here, ψ is a suitable weight function. Assuming that each pedestrian reacts to evaluations of averages of the density, we obtain

∂ t ρ + div    ρ v(ρ) ν(x) -ε ∇(ψ(ρ) * η) 1 + ∇(ψ(ρ) * η) 2    = 0 .
Here, we avoid the introduction of ψ to limit the analytical technicalities, obtaining (1.1)-(1.2). We assume throughout that ε > 0 is fixed and that the mollifier η satisfies

(η) η ∈ C 3 c (R N ; R + ) with R N η(x) dx = 1.
This convolution kernel has a key role: the individual at x deviates from the optimal path considering the crowd present in the region x -spt η, when no walls are present. The value η(ξ) is the relevance that the individual at x gives to the density ρ(x -ξ) located at x -ξ.

The next result shows that the present model fits in the framework described in Section 2 when applied on all of R N . Lemma 3.2. Fix ε > 0 and let η satisfy (η). Then, the operator I in (1.2) satisfies (I).

When the region Ω reachable by the crowd is restricted by the walls ∂Ω, we intend the convolution restricted to Ω

(ρ * η)(x) = Ω ρ(y) η(y -x) dy (3.3)
which coincide with the previous definition in the case spt ρ ⊂ Ω. A better choice is described in (A.1). The vector ν is here chosen as a sum ν = g + δ .

(3.4)

The former vector g is tangent to the "optimal" path, depending only on the geometry of the environment and coherent with it. Hence we assume that

(g) g ∈ C 2 (R N ; S N -1 ) satisfies ∇g ∈ L ∞ (R N ; R N ×N ), div g ∈ (W 1,1 ∩ W 1,∞ )(R N ; R)
and the invariance condition g(x)

• n(x) 0 holds for all x ∈ ∂Ω.

The latter vector δ describes the discomfort felt when passing too near to walls or obstacles. Below, we choose δ(x) = λ α(x) n(x) .

(3.5)

Here, λ ∈ R + is a suitable constant and n(x) is the inward normal (3.1). The role of the function α is to confine this discomfort to the region near the walls, i.e. we require

(α) α ∈ C 2 (R N ; [0, 1]) is such that α(x) = 0 whenever B(x, r Ω ) ⊆ Ω and α(x) = 1 for x ∈ ∂Ω.
A better choice for the discomfort is discussed in Appendix A. The present setting (1.1)-(1.2) can be effectively applied also in presence of walls, as shown by the next result.

2). Assume moreover that either

∂Ω is compact, or α ∈ (W 2,1 ∩ W 2,∞ ) and d ∂Ω ∈ (W 3,1 ∩ W 3,∞ ). Then, (1.1)-(1.
2) satisfies the assumptions of Theorem 2.3 and Theorem 2.4. Furthermore, if λ ε R ∇η L 1 , then property (P) holds.

The Model (1.1)-(1.3)

In the framework of (1.1), we now extend the model [13, (2.1)-(2.4)] to take into account the effects of crowd density on the speed modulus, i.e. of v(ρ). Moreover, we interpret the nonlocal term as a weighted average of the density, where the amount of crowd in the direction of the optimal path g is given more importance.

Using the same notation as above, we thus consider the model (1.1) with

I(ρ) = -ε ∇ R N ρ(y) η(x -y) ϕ (y -x) • g(x) dy 1 + ∇ R N ρ(y) η(x -y) ϕ (y -x) • g(x) dy 2 . (3.6)
Here, ε, η and v are as in § 3.1 and, in particular, ν is as in (3.4). The (almost) isotropic convolution (3.3) is here weighted by ϕ, whose argument is essentially the angle between the preferred path g and y -x, the point x being the position of the individual an the point y being the location where the density is evaluated. This term takes into account the preference of each individual to deviate little from the preferred path defined by g. For example, let ϕ ∈ C ∞ (R, [0, 1]) be such that ϕ ≡ 0 on ]-∞, 0] and ϕ ≡ 1 on [ϑ, +∞[, where ϑ > 0 is a given parameter. Then, adding the function ϕ (y -x) • g(x) into the nonlocal term means that the individual at x reacts to the average density evaluated in the prefer ed direction g(x).

The denominator is a regularized normalization. Its presence is necessary from the modeling point of view, for coherence with the presence of v(ρ). From the analytical point of view, this normalization makes various expressions slightly more complicate, but all estimates remain doable.

We are thus lead to the equation

∂ t ρ + div       ρ v(ρ)      ν(x) -ε ∇ R N ρ(y) η(x -y) ϕ (y -x) • g(x) dy 1 + ∇ R N ρ(y) η(x -y) ϕ (y -x) • g(x) dy 2            = 0 (3.7)
and we verify that it fits in the analytical framework provided in Section 2.

Lemma 3.4. Fix ε > 0, η ∈ C 3 c (R N ; R + ) with R N η(x) dz = 1, g ∈ W 3,∞ (R; [0, 1]) and ν ∈ W 3,∞ (R N , R).
Then, the operator I defined by (1.3) satisfies (I).

Passing now to the case in which the crowd's movement is constrained by the walls ∂Ω,as in the preceding section, we intend all integrals in (3.6)-(3.7) restricted to Ω. In particular, we consider now

I(ρ) = -ε ∇ Ω ρ(y) η(x -y) ϕ (y -x) • g(x) dy 1 + ∇ Ω ρ(y) η(x -y) ϕ (y -x) • g(x) dy 2 , (3.8)
with v as in (v) and ν as in (3.4). The applicability of theorems 2. Moreover, assume that ϕ ′ 0 and call ℓ = diam(spt η). Then,

(C 3 ∩ W 3,∞ )(R N ; S N -1 ). Let ϕ ∈ (C 3 ∩ W 3,∞ )(R; R). Define ν by (3.
ϕ ′ 0 and λ R η W 1,1 ϕ W 1,∞ 1 + ℓ ∇g L ∞
imply that property (P) holds.

The Model [13, (2.1)-(2.4)]

The model in [START_REF] Cristiani | Multiscale modeling of granular flows with application to crowd dynamics[END_REF][START_REF] Piccoli | Time-evolving measures and macroscopic modeling of pedestrian flow[END_REF], although set therein in the space M(R N ; R + ) of positive Radon measures, can be seen as a particular case of (1.1) setting

v(ρ) = 1 ν(x) = v des (x) and I(ρ) = ε R N ρ(y) ∇η(x -y) ϕ (y -x) • ν(x) dy , (3.9)
v des being the desired speed, see [13, formula (2.4)]. When η is radially symmetric, we recover [13, formula (2.6)] with η(x) = η( x ) and η′ = f . We leave to [START_REF] Cristiani | Multiscale modeling of granular flows with application to crowd dynamics[END_REF] the motivations of this model.

Proposition 3.6. Let ε > 0, (v), (ν) and (η) hold. Assume that ϕ ∈ (C 2 ∩ W 2,∞ )(R; R).
Then, (1.1)-(3.9) satisfies the assumptions of Theorem 2.3 and Theorem 2.4.

In particular, above we prove that if the initial datum ρ 0 is in L 1 (R N ; [0, R]), then the corresponding solution satisfies the same bounds. This ensures that neither focusing to any Dirac delta takes place, nor values of the density above R can be expected. (For the sake of completeness, we note that the conditions η 0 and R N η dx = 1 is in the case of (1.1)-(3.9) neither necessary, nor meaningful and can be replaced by η ∈ (C 3 c ∩ W 2,∞ ∩ W 3,1 )(R N ; R), see the proof in Section 5 for more details.).

Qualitative Properties

This section is devoted to sample numerical integrations of (1.1). In all the examples below, we choose as vector field ν = ν(x) the geodesic one, computed solving numerically the eikonal equation. This leads to the formation of congested queues near to the door jambs. From the modeling point of view, a more refined choice would consist in choosing ν so that most pedestrians are directed towards the central part of the exit. This choice increases the difficulties neither of the analytical treatment nor of the numerical integration but imposes the introduction of several further parameters.

The algorithm used is the classical Lax-Friedrichs method with dimensional splitting. As usual, a fixed grid (x i , y j ) for i 0 , . . . , n x and j = 1, . . . , n y is introduced and the density ρ is approximated through the values ρ ij on this grid. At every time step, the convolution in vector I(ρ) is then computed through products of the type A ih ρ hk B kj , for two fixed matrices A and B depending only on the kernel η.

All the examples below are set in R 2 , due to obvious visualization problems in higher dimensions. As is well known, the analytical techniques are essentially independent from the dimension as also the numerical algorithm. The time of integration obviously increases with the dimension.

Lane Formation

A widely detected pattern formed in the context of crowd dynamics is that of lane formation, see for instance [START_REF] Helbing | Simulation of pedestrian crowds in normal and evacuation situations, volume Pedestrian and Evacuation Dynamics[END_REF][START_REF] Helbing | Quantitative agent-based modeling of human interactions in space and time[END_REF][START_REF] Hoogendoorn | Simulation of pedestrian flows by optimal control and differential games[END_REF][START_REF] Piccoli | Pedestrian flows in bounded domains with obstacles[END_REF]. This feature has been often related to the specific qualities of each individual, i.e. it has usually been explained from a microscopic point of view. Here, in a purely macroscopic setting, we show that the solutions to (1.1)-(1.2) also display this pattern formation phenomenon, with pedestrian self organizing along lanes.

ν(x) = 1 0 + δ(x) , η(x) = 1 -x 1 r 2 3 1 -x 2 r 2 3 χ [-r,r] 2 (x) , r = 4 5 , v(ρ) = 1 2 (1 -ρ) , ρ 0 (x) = χ [3/5,4]×[-3/5,3/5] (x) , ε = 2 5 , (4.1)
where δ = δ(x) describes the discomfort due to walls: it is a vector normal to the walls, pointing inward, with intensity 3/2 along the walls, decreasing linearly to 0 at a distance 3/10 from the walls. As Figure 1 shows, the initially uniform crowd distribution evolves into a patterned configuration, first with 4 lanes and then with 5 lanes. The number of lanes depends on the size of the support of the convolution kernel η. Indeed, keeping all functions and parameters fixed, but not the parameter r, which determines the size of spt η, we obtain patterns differing in the number of lanes, see Figure 2. 1) computed at time t = 9.014 and with spt η with radius r = 0.5, 0.6, 0.8, 0.9, 1.0, 1.4. Note that as r increases, the number of lanes diminishes.

The formation of lanes is a rather stable phenomenon. Indeed, Figure 3 shows the result of the integration of (1.1)-(1.2)-(3.4)-(4.1) computed at time t = 0, 5, 10 with r = 3/5 (above) and r = 9/10 (below) with initial data different from that in (4.1). In both cases, lanes are formed similar to the corresponding situations in Figure 2.

We also note that in the present framework, using the terms in [28, Section 5.4], lanes form also in an isotropic setting. Indeed, the integrations in figures 1-2 were obtained with individuals able to see both forward and behind.

Evacuation of a Room

A standard application of macroscopic models for crowd dynamics is the minimization of evacuation times. The present setting applies to general geometries, see the assumption (Ω).

Here we show that (1.1) captures reasonable features of the escape dynamics. We consider a room with an exit, as in Figure 4. The vector ν = ν(x) is chosen as the unit vector tangent at x to the geodesic connecting x to the exit. The discomfort d = d(x) is a vector normal to the walls, pointing inward, with intensity 1 along the walls, decreasing linearly to 0 at a distance 1/2 from the walls. The other quantities are in (4.2). Keeping the

v(ρ) = 6(1 -ρ) η(x) = 1 -x 1 r 2 3 1 -x 2 r 2 3 χ [-r,r] 2 (x) , r = 0.6 ρ 0 (x) = 0.75 χ [2,7]×[-2,2] (x) ε = 0.4 (4.2)
above parameters fixed, as well as the outer walls of the room, we insert various obstacles (columns) to direct the movement of the crowd. First, Figure 5, first line, shows an integration of the case with two columns that direct people towards the exit. The number of lanes self adapts to the available space, with three lanes merging into one before the bottleneck. On the second line of Figure 5, the insertion of three columns in these positions delays but does not avoid the congestion at the exit. Those individuals that pass through the bottleneck are favored in exiting the room.

In Figure 6, first line, the insertion of four columns is more successful. Note, in the first two diagrams, that the number of lanes changes from 7 before the bottleneck, to 4 in between it. The individuals that do not pass through the bottleneck are penalized by the high density they meet near the door jambs. Finally, on the second line in Figure 6, an asymmetric layout hinders the lanes' pattern.

On the Rise of Singularities

A typical feature of conservation laws is the possible rise of singularities, see for instance [START_REF] Bressan | Hyperbolic systems of conservation laws[END_REF]Example 1.4]. The nonlocal equation (1.1) shares this characteristic. Indeed, assume ρ = ρ(t, x) is a given solution to (1.1), smooth up to time T > 0. Then, setting w(t, x) = ν(x) + I(t) (x) , simple computations lead to the following equation for the space derivative ρ j of ρ in the direction x j , for j = 1, . . . , N :

∂ t ρ j + (div ρ j ) q ′ w = (ρ j ) 2 q ′′ w + ρ j q ′′ i =j (∂ i ρ) w + ρ j q ′ div w + q ′ (∂ i ρ) ∂ j w + q ∂ j div w .
The first term in the right hand side is quadratic in ρ j , showing that a blow up of ρ j may take place in finite time. Moreover, we consider (1.1) with I as in (1.2) and

v(ρ) = 1 -ρ , η(x, y) = (1 -16x 2 ) 3 (1 -16y 2 ) 3 χ [-1/4,1/4] 2 (x, y) , ν(x) = 0 , ρ 0 (x, y) = (1 -4x 2 /9) 2 (1 -4y 2 /9) 2 χ [-3/2,3/2] 2 (x, y) , ε = -1 (4.
3) and we obtain the solution in Figure 7. 5 Technical Details

Existence and uniqueness

Proof of Lemma 2.2. Let q(ρ) = ρ v(ρ). Thanks to (ν), the assumptions on r and (I.1), we have for all M > 0

∂ ρ f (t, x, ρ) = q ′ (ρ) ν(x) + I r(t) (x) ; ∂ ρ f ∈ L ∞ ([0, T ] × R N × [-M, M ]; R N ) ; div f (t, x, ρ) = q(ρ) div ν(x) + I r(t) (x) ; div f ∈ L ∞ ([0, T ] × R N × [-M, M ]; R) .
Thus, we can apply Kružkov Theorem [24, Theorem 5 & 5.4] and ensure that (2.2) admits a unique Kružkov solution ρ ∈ L ∞ R + ; L 1 (R N ; R) , which is continuous from the right in time.

Lemma 5.1.

Let ρ 0 ∈ L 1 (R N ; [0, R]), r ∈ C 0 R + ; L 1 (R N ; [0, R]) . Under hypotheses (v)- (ν)-(I.1), the Cauchy problem (2.
2) admits a unique weak entropy solution ρ, with ρ ∈

C 0 R + ; L 1 (R N ; [0, R]) satisfying, for all t ∈ R + , ρ(t) L 1 = ρ 0 L 1 . (5.1)
Assume, in addition, that (I.2) is satisfied. Then,

ρ 0 ∈ BV(R N ; [0, R]) implies ρ(t) ∈ BV(R N ; [0, R]
) for all time t 0. Moreover, the following bound is satisfied

TV ρ(t) TV (ρ 0 ) + t N W N q L ∞ ([0,R]) ∇ div ν L 1 + C I ( r L ∞ ([0,T ];L 1 ) ) e κ * 0 t , (5.2)
where W N = π/2 0 (cos ϑ) N dϑ and the constant κ * 0 is bounded above as follows:

κ * 0 (2N + 1) q ′ L ∞ ([0,R]) ∇ν L ∞ + C I ( r L ∞ ([0,T ];L 1 ) ) .
Proof. Below we denote q(ρ) = ρ v(ρ), so that f (t, x, ρ) = q(ρ) ν(x) + I r(t) (x) .

The existence of a solution follows from Lemma 2.2. The rest of the proof is obtained through the following steps.

Estimates in L 1 and L ∞ . Since ρ is a solution to a conservation law, (5.1) is immediately satisfied.

Besides, ρ ≡ 0 and ρ ≡ R are solutions to (2.2) associated respectively to the constant initial conditions ρ 0 ≡ 0 and ρ 0 ≡ R. Hence, we can apply the comparison Theorem [24, Theorem 3] and obtain that if the initial condition takes value in [0, R], then the solution takes values in

[0, R], so that ρ ∈ L ∞ R + ; L 1 (R N ; [0, R]) .
Continuity in time. Thanks to [START_REF] Colombo | Stability and total variation estimates on general scalar balance laws[END_REF]Remark 2.4] or [START_REF] Lécureux-Mercier | Improved stability estimates on general scalar balance laws[END_REF]Corollary 2.4], ρ is continuous in time if, for any T > 0,

q L ∞ ([0,R]) T 0 R N div ν(x) + I r(t) (x) dx dt < ∞ . (5.3)
We have, thanks to (I.1), for all t 0 div ν + I r(t)

L 1 div ν L 1 + C I ( r(t) L 1 ) .
Hence, condition (5.3) is satisfied under hypotheses (v)-(ν)-(I.1) and we also obtain

ρ ∈ C 0 (R + ; L 1 (R N ; [0, R])).
Estimate in BV. To prove the bound on the TV norm we use [START_REF] Lécureux-Mercier | Improved stability estimates on general scalar balance laws[END_REF]Theorem 2.2]. To this aim, we have to check that for any T > 0,

∇∂ ρ f ∈ L ∞ ([0, T ]×R N ×[0, R]; R N ×N ) and T 0 R N ∇ div f (t, x, •) L ∞ ([0,R]) dx dt < ∞ .
Note first that ∇∂ ρ f (t, x, ρ) = q ′ (ρ) ∇ ν + I r(t) . Hence, thanks to (v)-(ν) and (I.1), we have

∇∂ ρ f (t, x, ρ) L ∞ q ′ L ∞ ([0,R]) ∇ν L ∞ + C I ( r(t) L 1 ) ,
Thanks to (v)-(ν) and (I.2), we have, for all t 0 ∇ div ν + I r(t)

L 1 ∇ div ν L 1 + C I ( r(t) L 1 )
.

Applying [25, Theorem 2.2], we obtain

κ * 0 = (2N + 1) ∇∂ ρ f L ∞ (2N + 1) q ′ L ∞ ([0,R]) ∇ν L ∞ + C I ( r L ∞ ([0,T ];L 1 ) ) , TV (ρ(t)) TV (ρ 0 ) e κ * 0 t +te κ * 0 t N W N q L ∞ ([0,R]) ∇ div ν L 1 + C I ( r L ∞ ([0,T ];L 1 ) ) ,
completing the proof.

Proof of Theorem 2.3. Thanks to Lemma 5.1, under hypotheses (v)-(ν) and (I), for any r ∈ C 0 [0, T ]; L 1 (R N ; [0, R]) there exists a unique solution ρ ∈ C 0 R + ; L 1 (R N ; [0, R]) to (2.2) satisfying (5.1)-(5.2). For any time T > 0, introduce the space

X(T ) = r ∈ C 0 [0, T ]; L 1 (R N ; [0, R]) : ∀ t ∈ [0, T ], r(t) L 1 ρ 0 L 1 .
Note that X(T ) is a Banach space with respect to the norm r X(T ) = sup t∈[0,T ] r(t) L 1 .

Let Q be the map that associates to r ∈ X(T ) the solution ρ ∈ X(T ) to (2.2). We want to find a condition on T so that Q is a contraction.

Let r 1 , r 2 ∈ X(T ). To apply [25, Theorem 2.6], we have to check that, for all t ∈ [0, T ],

∂ ρ (f 1 -f 2 ) ∈ L ∞ ([0, T ] × R N × [0, R]; R N ) , and T 0 R N div (f 1 -f 2 )(t, x, •) L ∞ ([0,R]) dx dt < +∞ .
Thanks to (I.3) we have

∂ ρ (f 1 -f 2 ) L ∞ = I(r 1 ) -I(r 2 ) q ′ L ∞ K I q ′ L ∞ ([0,R]) r 1 -r 2 L ∞ ([0,T ];L 1 ) , T 0 R N div (f 1 -f 2 )(t, x, •) L ∞ ([0,R]) dx dt T K I q L ∞ ([0,R]) r 1 -r 2 L ∞ ([0,T ];L 1 ) .
Hence, we get for all t 0

(ρ 1 -ρ 2 )(t) L 1 t K I q W 1,∞ ([0,R]) r 1 -r 2 L ∞ ([0,T ];L 1 ) × 1 + e κ * 0 t TV (ρ 0 ) + N W N 2 t q L ∞ ([0,R]) ∇ div ν L 1 + C I ( r L ∞ ([0,T ];L 1 )
) .

Denoting

k = (2N + 1) q ′ L ∞ ([0,R]) ∇ν L ∞ + C I ( ρ 0 L 1 ) , (5.4)
C 1 = K I q W 1,∞ ([0,R]) , C 2 = N W N 2 q L ∞ ([0,R]) ∇ div ν L 1 + C I ( ρ 0 L 1 ) .
The above estimate can be written

(ρ 1 -ρ 2 )(t) L 1 C 1 t r 1 -r 2 L ∞ ([0,T ];L 1 ) 1 + e kt TV (ρ 0 ) + C 2 t ,
We choose now T so that

C 1 T 1 + e kT TV (ρ 0 ) + C 2 T = 1 2
and, applying Banach Fixed Point Theorem, we obtain a unique fixed point for Q on X(T ). We now iterate the procedure above. To this aim, we use the total variation estimate (5.2)

TV ρ(t)

TV (ρ 0 ) + 2 C 2 t e kt ,
and, given T n , we recursively define T n+1 so that

C 1 (T n+1 -T n ) 1 + e kT n+1 TV (ρ 0 ) + C 2 e k(T n+1 -Tn) (2T n e kTn + T n+1 -T n ) = 1 2
and the above procedure ensures the existence of a fixed point on the interval [T n , T n+1 ].

The sequence T n is strictly increasing. If it is bounded, then the latter relation above yields 0 = 1/2. Hence, T n → ∞ when n → ∞, completing the proof.

Stability

Proof of Theorem 2.4. Here, we apply [START_REF] Lécureux-Mercier | Improved stability estimates on general scalar balance laws[END_REF]Corollary 2.8] to compare the solutions ρ 1 and ρ 2 of (2.4). Let k be as in (5.4). As in the proof of Theorem 2.3, (v)-(ν)-(I) ensure that the set of hypotheses [25, (H1*)-(H2*)-(H3*)] are satisfied. Thus, for any t 0, we obtain

ρ 1 (t) -ρ 2 (t) L 1 ρ 0,1 -ρ 0,2 L 1 + e kt TV (ρ 0,1 ) + N W N q 1 L ∞ ([0,R]) t 0 e k(t-τ ) R N ∇ div (ν 1 (x) + I 1 (ρ 1 )) dx dτ × t 0 q ′ 1 (ρ) ν 1 (x) + I 1 (ρ 1 ) -q ′ 2 (ρ) ν 2 (x) + I 2 (ρ 2 ) L ∞ (R N ) dτ + t 0 R N q 1 (ρ) div ν 1 (x) + I 1 (ρ 1 ) -q 2 (ρ) div ν 2 (x) + I 2 (ρ 2 ) L ∞ ([0,R]) dx dτ ρ 0,1 -ρ 0,2 L 1 + TV (ρ 0,1 ) + t N W N q 1 L ∞ ([0,R]) ∇ div ν 1 L 1 + C I ( ρ 0,1 L 1 ) × t q ′ 1 -q ′ 2 L ∞ ([0,R]) ν 1 + I 1 (ρ 1 ) L ∞ +t q ′ 2 L ∞ ([0,R]) ν 1 -ν 2 + I 1 (ρ 1 ) -I 2 (ρ 1 ) L ∞ + q ′ 2 L ∞ ([0,R]) t 0 I 2 ρ 1 (τ ) -I 2 ρ 2 (τ ) L ∞ dτ +t q 2 L ∞ ([0,R]) div ν 1 -ν 2 + I 1 (ρ 1 ) -I 2 (ρ 1 ) L 1 +t q 1 -q 2 L ∞ ([0,R]) div (ν 1 + I 1 (ρ 1 )) L 1 + q 2 L ∞ ([0,R]) t 0 div I 2 (ρ 1 ) -I 2 (ρ 2 ) dτ .
Hence, denoting

a(t) = TV (ρ 0,1 ) + t N W N q 1 L ∞ ([0,R]) ∇ div ν 1 L 1 + C I ( ρ 0,1 L 1 ) b(t) = K I q ′ 2 L ∞ ([0,R]) e kt a(t) + q 2 L ∞ ([0,R]) γ 1 = q ′ 1 -q ′ 2 L ∞ ([0,R]) ν 1 L ∞ + C I ( ρ 0,1 L 1 ) + q ′ 2 L ∞ ([0,R]) ν 1 -ν 2 L ∞ + q ′ 2 L ∞ ([0,R]) sup ρ∈L 1 (R N ;[0,R]) I 1 (ρ) -I 2 (ρ) L ∞ γ 2 = q 1 -q 2 L ∞ ([0,R]) div ν 1 L 1 + C I ( ρ 0,1 L 1 ) + q 2 L ∞ ([0,R]) div (ν 1 -ν 2 ) L 1 + sup ρ∈L 1 (R N ;[0,R]) div I 1 (ρ) -div I 2 (ρ) L 1 we get ρ 1 -ρ 2 (t) L 1 ρ 0,1 -ρ 0,2 L 1 + γ 1 t e kt a(t) + t γ 2 + b(t) t 0 ρ 1 -ρ 2 (τ ) dτ .
By Gronwall Lemma, we obtain

ρ 1 -ρ 2 (t) L 1 ρ 0,1 -ρ 0,2 L 1 + γ 1 t e kt a(t) + t γ 2 +b(t)e t 0 b(u)du t 0 e -τ 0 b(u)du ( ρ 0,1 -ρ 0,2 L 1 + γ 1 τ e kτ a(τ ) + τ γ 2 ) dτ .
That is to say, setting H ϑ (t) = {(y, z) : w(t, y + ϑz) • ñ(y) 0} and ∂Ω M = ∂Ω ∩ B(0, M )

R N ρ(0, x) -ρ(T, x) ψ(x) dx = T 0 ∂K M B(0,1)
q ρ(t, y + ϑz) w(t, y + ϑz) • ñ(y) µ( z ) dz dσ(y) dt

T 0 ∂K M B(0,1) χ H ϑ (t) (y, z) q ρ(t, y + ϑz) w(t, y + ϑz) • ñ(y) µ( z ) dz dσ(y) dt = T 0 ∂Ω M B(0,1) χ H ϑ (t) (y, z) q ρ(t, y + ϑz) w(t, y + ϑz) • ñ(y) µ( z ) dz dσ(y) dt + T 0 ∂K M \∂Ω M B(0,1) χ H ϑ (t) (y, z) q ρ(t, y + ϑz) w(t, y + ϑz) • ñ(y) µ( z ) dz dσ(y) dt T 0 ∂Ω M B(0,1) χ H ϑ (t) (y, z) q L ∞ w(t, y + ϑz) • ñ(y) µ( z ) dz dσ(y) dt + m
where m can be arbitrarily small provided M is sufficiently large, since ρ ∈ L 1 ([0, T ] × R N ; [0, R]). Letting ϑ → 0, M → ∞ and using (3.2), by the Dominated Convergence Theorem we get

0 c Ω ρ(0, x) -ρ(T, x) dx = - c Ω ρ(T, x) dx ,
since by hypothesis spt ρ 0 ⊂ Ω. Besides, the positivity of ρ gives usc Ω ρ(t, x) dx 0. Finally, we have c Ω ρ(T, x) = 0 and spt ρ(T ) ⊂ Ω.

Proof of Lemma 3.2. Let ρ ∈ L 1 (R N ; [0, R]). Note that I can be rewritten as I(ρ) = -ε ρ * (∇η)

1 + ρ * ∇η 2 . By the assumptions on η and the properties of the convolution product, we have

I(ρ) L ∞ ε ρ L 1 ∇η L ∞ , I(ρ) L 1 ε ρ L 1 ∇η L 1 , ∇I(ρ) L ∞ ε ρ L 1 ∇η W 1,∞ 1 + R 2 ∇η L 1 ∇ 2 η L 1 , ∇I(ρ) L 1 ε ρ L 1 ∇η W 1,∞ 1 + R 2 ∇η L 1 ∇ 2 η L 1
.

Hence (I.1) is satisfied. Let us check now (I.2)

∇ 2 I(ρ) L 1 ε R 2 ρ L 1 ∇η W 2,1 ∇η 2 W 1,1 (1 + 4 + 3R 2 ∇η 2 W 1,1 ) . Let r 1 , r 2 ∈ L 1 (R N ; [0, R]).
We have:

I(r 1 ) -I(r 2 ) L ∞ ε 1 + R 2 ∇η 2 L 1 ∇η L ∞ r 1 -r 2 L 1 , I(r 1 ) -I(r 2 ) L 1 ε 1 + R 2 ∇η 2 L 1 ∇η L 1 r 1 -r 2 L 1 , div (I(r 1 ) -I(r 2 )) L 1 ε r 1 -r 2 L 1 ∇η W 1,1 (1 + 8R 2 ∇η 2 W 1,1 + 3R 4 ∇η 4 W 1,1
) , completing the proof.

Proof of Proposition 3.3. The fact that (ν) holds follows through simple computations from (Ω), (g), (η), (α) and the assumptions on ∂Ω or α, d ∂Ω . Condition (I) follows from Lemma 3.2.

To prove the invariance property (P), we verify that g(x) + δ(x) + I(ρ) (x) • n(x) 0 for all x ∈ ∂Ω. By (g), (3.5) and (1.2), it is sufficient to prove that

λ ε sup x∈∂Ω (ρ * ∇η)(x) • n(x) 1 + (ρ * ∇η)(x) 2 . Note that sup x∈∂Ω (ρ * ∇η)(x) • n(x) 1 + (ρ * ∇η)(x) 2 sup x∈∂Ω (ρ * ∇η)(x) • n(x) ρ L ∞ ∇η L 1 n L ∞ R ∇η L 1
completing the proof, by Proposition 3.1.

Proof of Lemma 3.4. Note that I as

I(ρ) = ε R N ρ(y) ∇ η(x -y) g (y -x) • ν(x) dy. Let r η be such that B(0, r η ) ⊃ spt η. Fix ρ ∈ L 1 (R N ; [0, R]).
Let us check (I.1):

I(ρ) L ∞ ε ρ L 1 η W 1,∞ g W 1,∞ 1 + (1 + r η ) ν W 1,∞ , I(ρ) L 1 ε ρ L 1 η W 1,1 g W 1,∞ 1 + (1 + r η ) ν W 1,∞ , ∇I(ρ) L ∞ ε ρ L 1 η W 2,∞ g W 2,∞ × 1 + 2(r η + 1) ν W 1,∞ + (1 + r η ) 2 ν 2 W 1,∞ + (2 + r η ) ∇ν W 1,∞ , ∇I(ρ) L 1 ε ρ L 1 η W 2,1 g W 2,∞ × 1 + 2(r η + 1) ν W 1,∞ + (1 + r η ) 2 ν 2 W 1,∞ + (2 + r η ) ∇ν W 1,∞ .
Passing to (I.2):

∇ 2 I(ρ) L 1 ε ρ L 1 η W 3,1 g W 3,∞ × 1 + 3(1 + r η ) ν W 1,∞ + 3(1 + r η ) 2 ν 2 W 1,∞ + 3(2 + r η ) ∇ν W 1,∞ +(1 + r η ) 3 ν 3 W 1,∞ + 3(2 + r η ) 2 ν 2 W 2,∞ + (3 + r η ) ∇ 2 ν W 1,∞ . Let r 1 , r 2 ∈ L 1 (R N ; [0, R]). The operator I is linear in ρ, hence I(r 1 ) -I(r 2 ) L ∞ ε r 1 -r 2 L 1 η W 1,∞ g W 1,∞ (1 + (1 + r η ) ν W 1,∞ ) , I(r 1 ) -I(r 2 ) L ∞ ε r 1 -r 2 L 1 η W 1,1 g W 1,∞ (1 + (1 + r η ) ν W 1,∞ ) , div (I(r 1 ) -I(r 2 )) L ∞ ε r 1 -r 2 L 1 η W 2,1 g W 2,∞ × 1 + 2(1 + r η ) ν W 1,∞ + (1 + r η ) 2 ν 2 W 1,∞ +(2 + r η ) ∇ν W 1,∞
showing that (I.3) is satisfied.

Proof of Proposition 3.5. Fix ρ ∈ L 1 (R N ; [0, R]). Let us check (I.1) and (I.2):

I(ρ) L ∞ ε ρ L 1 η W 1,∞ ϕ W 1,∞ 1 + (1 + ℓ) g W 1,∞ , I(ρ) L 1 ε ρ L 1 η W 1,1 ϕ W 1,∞ 1 + (1 + ℓ) g W 1,∞ , ∇I(ρ) L ∞ ε ρ L 1 η W 2,∞ ϕ W 2,∞ × 1 + 2(ℓ + 1) g W 1,∞ + (1 + ℓ) 2 g 2 W 1,∞ + (2 + ℓ) ∇g W 1,∞ × 1 + R 2 η 2 W 1,1 ϕ 2 W 1,∞ 1 + (1 + ℓ) g W 1,∞ 2 , ∇I(ρ) L 1 ε ρ L 1 η W 2,1 ϕ W 2,∞ × 1 + 2(ℓ + 1) g W 1,∞ + (1 + ℓ) 2 g 2 W 1,∞ + (2 + ℓ) ∇g W 1,∞ × 1 + R 2 η 2 W 1,1 ϕ 2 W 1,∞ 1 + (1 + ℓ) g W 1,∞ 2 .
Passing to (I.2):

∇ 2 I(ρ) L 1 ε ρ L 1 η W 3,1 ϕ W 3,∞ × 1 + 3(1 + ℓ) g W 1,∞ + 3(1 + ℓ) 2 g 2 W 1,∞ + 3(2 + ℓ) ∇g W 1,∞ +(1 + ℓ) 3 g 3 W 1,∞ + 3(2 + ℓ) 2 g 2 W 2,∞ + (3 + ℓ) ∇ 2 g W 1,∞ × 1 + R 2 η 2 W 1,1 ϕ 2 W 1,∞ 1 + (1 + ℓ) g W 1,∞ 2 +3 ε R η W 1,1 ϕ W 1,∞ 1 + (1 + ℓ) g W 1,∞ × ρ L 1 η W 2,1 ϕ W 2,∞ × 1 + 2(ℓ + 1) g W 1,∞ + (1 + ℓ) 2 g 2 W 1,∞ + (2 + ℓ) ∇g W 1,∞ × 1 + R 2 η 2 W 1,1 ϕ 2 W 1,∞ 1 + (1 + ℓ) g W 1,∞ 2 ε ρ L 1 η W 3,1 ϕ W 3,∞ × 1 + 3(1 + ℓ) g W 1,∞ + 3(1 + ℓ) 2 g 2 W 1,∞ + 3(2 + ℓ) ∇g W 1,∞ +(1 + ℓ) 3 g 3 W 1,∞ + 3(2 + ℓ) 2 g 2 W 2,∞ + (3 + ℓ) ∇ 2 g W 1,∞ +3 R η W 1,1 ϕ W 1,∞ 1 + (1 + ℓ) g W 1,∞ × 1 + 2(ℓ + 1) g W 1,∞ + (1 + ℓ) 2 g 2 W 1,∞ + (2 + ℓ) ∇g W 1,∞ × 1 + R 2 η 2 W 1,1 ϕ 2 W 1,∞ 1 + (1 + ℓ) g W 1,∞ 2 .
Finally, (I.3) is proved exactly as in Lemma 3.2.

To prove the invariance property (P), following the same procedure as in the proof of Proposition 3.3, we check that ν(x) + I(ρ) (x) • n(x) 0 for all x ∈ ∂Ω. In fact, if

λ(x) -ε R ϕ L ∞ ∇η L 1 -ε R ϕ ′ L ∞ Ω η(x -y) (y -x) ∇g(x) n(x) dy . λ(x) -ε R ϕ L ∞ ∇η L 1 + ϕ ′ L ∞ η L 1 ℓ ∇g L ∞
completing the proof, by Proposition 3.1.

Proof of Proposition 3.6. Denote ℓ = diam spt η. Straightforward computations give:

I(ρ) L ∞ ε ρ L 1 ∇η L ∞ ϕ L ∞ , ∇I(ρ) L ∞ ε ρ L 1 ∇ 2 η L ∞ ϕ L ∞ + ∇η L ∞ ϕ ′ L ∞ g L ∞ + ℓ ∇g L ∞ , I(ρ) L 1 ε ρ L 1 ∇η L 1 ϕ L ∞ , ∇I(ρ) L 1 ε ρ L 1 ∇ 2 η L 1 ϕ L ∞ + ∇η L 1 ϕ ′ L ∞ g L ∞ + ℓ ∇g L ∞ , ∇ 2 I(ρ) L 1 ε ρ L 1 ∇η W 2,1 ϕ W 2,∞ 1 + 3 g L ∞ + ∇g L ∞ (3ℓ + 2) + ℓ ∇ 2 g L ∞
giving (I.1) and (I.2). The proof of (I.3) is immediate by the linearity of I in ρ.

A Appendix: Geometrical Issues Related to Ω

The framework presented in Section 3 can be adapted to various real situations. For instance, the "walls" ∂Ω may stop the visibility of the individuals. Then, from a modeling point of view, it can be reasonable to introduce the set The above relation means that the individual at x evaluates an average of the densities ρ(y) at all values y that are visible from x. With these choices, the validity of condition (I) essentially

The choice (1.2)-(A.2) is appealing from the modeling point of view, but not easily tractable from both the numerical and the analytical points of view, without major restrictions on the geometry of Ω and on η. Therefore, we consider also the following choice. where α is as in (α).

The proof is entirely similar to that of Proposition A.1.

Proposition 3 . 1 .

 31 Let ν, I and Ω satisfy (ν), (I) and ( Ω). If for all x ∈ ∂Ω and ρ ∈

Proposition 3 . 3 .

 33 Let ε > 0, (v), ( Ω), (η), (g) and (α) hold. Define ν by (3.4), δ by (3.5) and I by (1.

Proposition 3 . 5 .

 35 3 and 2.4 to (1.1)-(3.8) and the validity of property (P) is ensured by the following proposition. Let ε > 0, (v), ( Ω), (η) and (g) hold, with moreover g ∈

  4), δ by (3.5) and I by (3.6). Then, (1.1)-(3.6) satisfies the assumptions of Theorem 2.3 and Theorem 2.4.

Figure 1 :

 1 Figure 1: Solution to (1.1)-(1.2)-(3.4)-(4.1) at times t = 0, 2.529, 5.043, 7.557, 10.071, 15.014. First 3 lanes are formed, then the middle lane bifurcates forming the fourth lane.

Figure 2 :

 2 Figure 2: Solution to (1.1)-(1.2)-(3.4)-(4.1) computed at time t = 9.014 and with spt η with radius r = 0.5, 0.6, 0.8, 0.9, 1.0, 1.4. Note that as r increases, the number of lanes diminishes.

Figure 3 :

 3 Figure 3: Solution to (1.1)-(1.2)-(3.4)-(4.1) with different initial data at time t = 0, 5.014, 10.014 and above with r = 0.6, below with r = 0.9. Note that above 4 lanes form and below 5, similarly to what obtained in Figure 1.

Figure 4 :

 4 Figure 4: Initial datum and room geometry considered in § 4.2.

Figure 5 :

 5 Figure 5: Solution to (1.1)-(1.2)-(3.4)-(4.2) with different geometries, computed at time t = 2.521, 5.043 and 7.563.

Figure 6 :

 6 Figure 6: Solution to (1.1)-(1.2)-(3.4)-(4.2) with different geometries, computed at time t = 2.521, 5.043 and 7.563.

Figure 7 :

 7 Figure 7: Left, the initial datum and, right, the solution to (1.1)-(1.2)-(4.3) at time t = 1. Note the formation of vertical faces in the initially smooth distribution.

x

  ∈ ∂Ω, ν(x) + I(ρ) (x) • n(x) ) + λ n(x) -ε ∇ Ω ρ(y) η(x -y) ϕ (y -x) • g(x) dy 1 + ∇ Ω ρ(y) η(x -y) ϕ (y -x) • g(x) dy ) + g(x) • n(x) -ε ϕ L ∞ Ω ρ(y) ∇ η(x -y) • n(x) dy -ε Ω ρ(y) η(x -y) ϕ ′ (y -x) • g(x) (y -x) ∇g(x) n(x) -g(x)• n(x) dy

Ω

  x = y ∈ Ω : x + σ(y -x) ∈ Ω , ∀ σ ∈ [0, 1] (A.1)of the points in Ω visible from x. Correspondingly, the nonlocal operator I in (1.2) can be modified intending the convolution as follows:(ρ * η)(x) = Ωx ρ(y) η(x -y) dy . (A.2)

Proposition A. 2 .

 2 Assume that ( Ω) and (η) hold. Let I be as in(1.2), with the usual convolution(3.3). Ifδ(x) = ε R Hx ∇η(x -y) dy for x ∈ ∂Ω then δ(x) + I(ρ)(x) • n(x) 0 for all x ∈ ∂Ω and ρ ∈ L 1 (R N ; [0, R]), so that (P) holds.With this choice, a way to define δ on all Ω could beδ(x) = ε R α(x)Hx ∇η(x -y) dy ,

  is a weak entropy solution to (1.1) if it is a Kružkov solution to the Cauchy problem for the scalar conservation law

Dipartimento di Matematica, Università degli studi di Brescia, Italia. rinaldo@ing.unibs.it

Di.S.T.A., Università del Piemonte Orientale, Italia, mauro.garavello@mfn.unipmn.it

Université d'Orléans, UFR Sciences, Bâtiment de mathématiques -Rue de Chartres B.P. 6759 -45067 Orléans cedex 2, France, magali.lecureux-mercier@univ-orleans.fr

we get ρ 1 -ρ 2 (t) L 1 ρ 0,1 -ρ 0,2 L 1 + γ 1 t e kt a(t) + t γ 2 b(t) b(0) e tb (t) completing the proof.

Driving examples

Proof of Proposition 3.1. Let ρ solve (1.1) in the sense of Definition 2.1. Define w(t, x) = ν(x)

) is also a weak solution to (2.1).

Then, set as above q(ρ) = ρ v(ρ) and, for any

ϑ N µ x /ϑ and define ψ = χ K M * µ ϑ , the convolution being over all of R N .

Note that ∇ψ = 0 on R N \ B(∂K M , ϑ). If x ∈ B(∂K M , ϑ), we have

where ñ(x) is the outer normal to K M at any x ∈ ∂K M . Then, we choose

Letting α → 0, we obtain

q ρ(t, y + ϑz) w(t, y + ϑz) • ñ(y) µ( z ) dz dσ(y) dt .

depends on the geometry of Ω. In particular, if Ω is convex, then Ω x = Ω and (I) holds by Proposition 3.3.

Here we only show how to choose a discomfort so that (P) holds for (1.1)-(1.2)-(A.2). The case of the nonlocal operator (1.3) is entirely similar.

To this aim, introduce the set

Proposition A.1. Assume that ( Ω) and (η) hold. Let I be as in (1.2), where the convolution is intended as in (A.2). If In the latter expression, the supremum is attained for ρ ⊇ Ω x .

The latter inclusion holds since Ω x is a convex set contained in Ω and containing x.


Yves Lechevallier
email: yves.lechevallier@inria.fr

Francisco De A T De Carvalho

Thierry Despeyroux
email: thierry.despeyroux@inria.fr

Filipe M De Melo

Clustering of Multiple Dissimilarity Data Tables for Documents Categorization

Keywords: Clustering Analysis, Relational Data, Documents Categorization

This paper introduces a clustering algorithm that is able to partition objects taking into account simultaneously their relational descriptions given by multiple dissimilarity matrices. These matrices could have been generated using different sets of variables and a fixed dissimilarity function, using a fixed set of variables and different dissimilarity functions or using different sets of variables and dissimilarity functions. This method, which is based on the dynamic hard clustering algorithm for relational data, is designed to provided a partition and a prototype for each cluster as well as to learn a relevance weight for each dissimilarity matrix by optimizing an adequacy criterion that measures the fit between clusters and their representatives. These relevance weights change at each algorithm iteration and are different from one cluster to another. Experiments aiming at obtaining a categorization of a document data base demonstrate the usefulness of this partitional clustering method.

Introduction

Clustering is a popular task in knowledge discovering and it is applied in various fields including data mining, pattern recognition, computer vision, etc [START_REF] Gordon | Classification. Chapman and Hall/CRC[END_REF], [START_REF] Jain | Data clustering: A review[END_REF]). Clustering methods aims at organizing a set of objects into clusters such that items within a given cluster have a high degree of similarity, while items belonging to different clusters have a high degree of dissimilarity. The most popular clustering techniques are hierarchical and partitioning methods. Partitioning methods seek to obtain a single partition of the input data into a fixed number of clusters. Such methods often look for a partition that optimizes (locally) an adequacy criterion function.

There are two common representations of the objects upon which clustering can be based : (usual or symbolic) feature data and relational data.

When each object is described by a vector of quantitative or qualitative values the set of vectors describing the objects is called a feature data. When each (complex) object is described by a vector of sets of categories, intervals or weight histograms, the set of vectors describing the objects is called a symbolic feature data. Symbolic data has been mainly studied in Symbolic Data Analysis (SDA) [START_REF] Bock | Analysis of Symbolic Data[END_REF]). Alternatively, when each pair of objects is represented by a relationship, then we have relational data. The most common case of relational data is when we have (a matrix of) dissimilarity data, say R = [r il ], where r il is the pairwise dissimilarity (often a distance) between objects i and l. This paper gives a clustering algorithm that is able to partition objects taking simultaneously into account their relational descriptions given by multiple dissimilaritity matrices. The main idea is to obtain a collaborative role of the different dissimilarity matrices [START_REF] Pedrycz | Collaborative fuzzy clustering[END_REF]) in order to obtain a final consensus partition [START_REF] Leclerc | Concensus en classification : une revue bibliographique[END_REF]). These dissimilarity matrices could have been generated using different sets of variables and a fixed dissimilarity function (the final partition gives a consensus between different views (sets of variables) describing the objects), using a fixed set of variables and different dissimilarity functions (the final partition gives the consensus between different dissimilarity functions) or using different sets of variables and dissimilarity functions. Moreover, the influence of the different dissimilarity matrices is not equally important in the definition of the clusters in the final consensus partition. Thus, in order to obtain a meaningful partition from all dissimilarity matrices, it is necessary to learn cluster-dependent relevance weights for each dissimilarity matrix. [START_REF] Frigui | Clustering and aggregation of relational data with applications to image database categorization[END_REF] proposed CARD, a clustering algorithm that is able to partition objects taking into account multiple dissimilaritity matrices and that learns a relevance weight for each dissimilarity matrix in each cluster. CARD is mainly based on the well know fuzzy clustering algorithms for relational data RFCM [START_REF] Hathaway | Relational duals of the c-means algorithms[END_REF]) and FANNY [START_REF] Kaufman | Finding Groups in Data[END_REF]).

The clustering algorithm given in this paper is designed to give a partition and a prototype for each cluster as well as to learn a relevance weight for each dissimilarity matrix by optimizing an adequacy criterion that measures the fit between clusters and their representatives. These relevance weights change at each algorithm iteration and are different from one cluster to another. However, this method is based on the dynamic hard clustering algorithm for relational data [START_REF] Lechevallier | Optimisation de quelques critères en classification automatique et application a l'étude des modifications des protéines sériques en pathologie clinique[END_REF][START_REF] De | Clustering methods in symbolic data analysis[END_REF], De Carvalho et al (2009)) as well as on the dynamic clustering method based on adaptive distances (Diday andGovaert (1977), De Carvalho andLechevallier(2009)). The adaptive dynamic clustering method gives a partition as well as a prototype for each cluster and learns a relevance weight for each variable in each cluster. In order to demonstrate the usefulness of this clustering algorithm, experiments were designed in order to obtain a categorization of a document data base.

This paper is organized as follows. Section 2 presents a partitioning clustering algorithm based on multiple dissimilarity matrices. In order to illustrate the usefulness of this clustering method, section 3 shows the application of this algorithm in order to obtain a categorization of a document data base. Finally, section 4 presents the conclusions.

A Dynamic Clustering Algorithm Based on Multiple Dissimilarity Matrices

In this section, we introduce an extension of the dynamic clustering algorithm for relational data [START_REF] De | Clustering methods in symbolic data analysis[END_REF]) which is able to partition objects taking simultaneously into account their relational descriptions given by multiple dissimilaritity matrices.

Let E = {e 1 , . . . , e n } be a set of n examples and let p dissimilarity n × n matrices (D 1 , . . . , D j , . . . , D p ) where D j [i, l] = d j (e i , e l ) gives the dissimilarity between objects e i and e l on dissimilarity matrix D j . Assume that the prototype g k of cluster C k belongs to the set of examples E, i.e., g k ∈ E ∀k = 1, . . . , K.

The dynamic hard clustering algorithm with relevance weight for each dissimilarity matrix looks for a partition P = (C 1 , . . . , C K ) of E into K clusters and the corresponding prototype g k ∈ E representing the cluster C k in P such that an adequacy criterion (objective function) measuring the fit between the clusters and their prototypes is locally optimized. The adequacy criterion is defined as

J = K k=1 ei∈C k d (k) (e i , g k ) = K k=1 ei∈C k p j=1 λ j k d j (e i , g k ) (1)
in which

d (k) (e i , g k ) = p j=1 λ j k d j (e i , g k ) (2)
is the dissimilarity between an example e i ∈ C k and the cluster prototype g k ∈ E parameterized by relevance weight vector λ k = (λ 1 k , . . . , λ j k , . . . , λ p k ) where λ j k is the weight between the dissimilarity matrix D j and the clusters C k , and d j (e i , g k ) is the local dissimilarity d j between an example e i ∈ C k and the cluster prototype g k ∈ E.

The relevance weight matrix λ composed by K relevance weight vectors λ k = (λ 1 k , . . . , λ j k , . . . , λ p k ) changes at each iteration, i.e., they are not determined absolutely, and are different from one cluster to another. Our clustering algorithm alternates the three following steps:

Step 1: Definition of the Best Prototypes

In this step, the partition P = (C 1 , . . . , C K ) of E into K clusters and the relevance weight matrix λ are fixed.

Proposition 1. The prototype g k = e l ∈ E of cluster C k , which minimizes the clustering criterion J, is computed according to:

l = arg min 1≤h≤n ei∈C k p j=1 λ j h d j (e i , e h ) (3)
Step 2: Definition of the Best Relevance Weight Matrix In this step, the partition P = (C 1 , . . . , C K ) of E and the vector of prototypes g = (g 1 , . . . , g K ) are fixed.

Proposition 2. The element j of the relevance weight vector λ k = (λ 1 k , . . . , λ p k ), which minimizes the clustering criterion J under λ j k > 0 and p j=1 λ j k = 1, is calculated by the following expression:

λ j k = p h=1 ei∈C k d h (e i , g k ) 1 p ei∈C k d j (e i , g k ) (4)
Step 3: Definition of the Best Partition In this step, the vector of prototypes g = (g 1 , . . . , g K ) and the relevance weight matrix λ are fixed.

Proposition 3. The cluster C k , which minimize the criterion J, is updated according to the following allocation rule:

C k = {e i ∈ E : d (k) (e i , g k ) < d (h) (e i , g h )∀h = k } (5)
If the minimum is not unique, e i is assigned to the class having the smallest index

It's easy to demonstrate that each preview step decreases the criterion J. The dynamic hard clustering algorithm with relevance weight for each dissimilarity matrix sets an initial partition and alternates three steps until convergence, when the criterion J(P, λ, g) reaches a stationary value representing a local minimum. This algorithm is summarized below.

The Dynamic Hard Clustering Algorithm with Relevance Weight Matrix

Initialization.

Fix the number K of clusters; Randomly select K distinct objects g k ∈ E; Set the Relevance Weight Matrix λ where λ k = (λ 1 k , . . . , λ p k ) = (1, . . . , 1); Assign each object e i to the closest prototype in order to obtain the partition P = (C 1 , . . . , C K ) where C k is constructed by the rule (5).

Step 1: definition of the best prototypes.

The partition P = (C 1 , . . . , C K ) and the relevance weight matrix λ are fixed.

Compute the prototype g k ∈ E of cluster C k according to equation (3) 3.

Step 2: definition of the best relevance weight matrix.

The vector of prototypes g and the partition P = (C 1 , . . . , C K ) are fixed.

For each k compute the component of the weight vector λ k according to equation (4) 4.

Step 3: definition of the best partition.

The vector of prototypes g and the relevance weight matrix λ are fixed.

Construct the new partition P ′ = (C ′ 1 , . . . , C ′ K ) with the rule given by ( 5) and control the convergence by: test ← 0; for i = 1 to n do e i belonged to the class C m and belongs to the winning cluster C ′ k if k = m then test ← 1; P ← P ′ ; 5. Stopping criterion. If test = 0 then STOP, otherwise go to 2 (Step 1).

Application: document data base categorization

To illustrate the usefulness of the proposed clustering algorithm, we use it to categorize a document data base. The document data base is a collection of reports produced by every Inria (The French National Institute for Research in Computer Science and Control) research team in 2007. Research teams are grouped into scientific themes that do not correspond to an organizational structure (such as departments or divisions), but act as a virtual structure for the purpose of presentation, communication and evaluation. Choice of themes and team allocation are mostly related to strategic objectives and scientific closeness between existing teams, but also take in account some geographical constraints, such as the desire for a theme to be representative of most Inria centers. Our aim is to compare the categorization given automatically by the clustering algorithm introduced in this paper with the a priori expert categorization given by INRIA.

These reports are written in English. The sources are LaTeX documents, and are automatically translated into XML, then to HTML to be published on the Web. In the rest of the paper we implicitly refer to the XML version of the Activity Report. The logical structure of the RA is defined by an XML DTD with a few mandatory sections and some optional parts.

In this application we considered activity reports from 164 INRIA research teams in 2007. On each activity report, 4 sections have been selected to describe a research team: overall objectives, scientific foundations, dissemination and new results. The overall objectives part defines the research objectives and scientific foundations provides the scientific background followed From these activity reports we initially obtained 4 feature data tables, each table with 164 individuals (INRIA research team) described by the frequent words (categories) present in one of 4 variables The number of frequent words in overall objectives section is 220, 210 for scientific foundations, 404 for dissemination and 547 for new results sections. Each cell on a data table gives the frequency of a word for the considered activity report section and research team.

Then, 4 relational data tables have been obtained from the 4 feature data tables through a dissimilarity measure derived from the affinity coefficient (Barcelar-Nicolau (2000)). We assume that each individual is described by one set-valued variable ("presentation", etc.) which has m j modalities (or categories) {1, . . . , m}. An individual e i is described by x i = (n i1 , . . . , n im ) where n ij is the frequency of modality j. The dissimilarity between a pair of individuals e i and e i ′ is given by:

d(x i , x ′ i ) = 1 - m j=1 n ij n i• n i ′ j n i ′ • where n i• = m j=1 n ij .
All these relational data tables were normalized according to their overall dispersion [START_REF] Chavent | Normalized k-means clustering of hyper-rectangles[END_REF]) to have the same dispersion. This means that each dissimilarity d(x i , x ′ i ) in a given relation data table has been normalized as

d(xi,x ′ i ) T
where T = n i=1 d(e i , g) is the overall dispersion and g = e l ∈ E = {e 1 , . . . , e n } is the overall prototype, which is computed according to l = argmin 1≤h≤n n i=1 d(e i , e h ).

Results

The clustering algorithm has been performed simultaneously on these 4 relational data tables ("presentation", "foundation", "dissemination" and "bibliography") in order to obtain a partition in K ∈ {1, . . . , 15}. For a fixed number of clusters K, the clustering algorithm is run 100 times and the best result according to the adequacy criterion is selected.

In order to determine the number of clusters, we used the approach described by Da Silva ( 2009), which consists on the choice of the peaks on the graph of the "second order differences" of the clustering criterion (equation ( 1)): J (K-1) + J (K+1) -2J (K) , K = 2, . . . , 14. According to this approach, we fixed the number of clusters in 4 and 9.

The 4-cluster and the 9-cluster partitions obtained with this clustering algorithm were compared with the INRIA research team categorization 5-class partition known a priori. The 5-class a priori categorization is as follows: "Applied Mathematics, Computation and Simulation (M)", "Algorithmics, Programming, Software and Architecture (A)", "Networks, Systems and Services, Distributed Computing (N) ", "Perception, Cognition, Interaction (P)" and "Computational Sciences for Biology, Medicine and the Environment (C)". These 5 categories are themselves divided into several sub-categories. In many points we retrieve in the 9-cluster partition the categorization done a priori by INRIA. For example the sub-category "Networks and Telecommunications" of N fits exactly in one cluster. The two sub-categories "Distributed Systems and Services" and "Distributed and High Performance Computing" are merged into a unique cluster, indicating that from the language used point of view the distinction between these two categories is artificial. Some teams have also migrate. For example it seems that the language used in Cryptography (that is part of A in the a priori categorization) is closer to the language used in math (M). Looking at the 4-cluster partition, some migrations are also clearly detected, which have a political sense, in particular when the concerned team is found in a cluster corresponding to the "right" category in a former categorization done by INRIA. Finally, is seems that teams in the C category share the same language as teams in M or in P stressing the fact that in the activity report the weight of the scientific foundations is important, and this fact showing up in both partitions is however clearer in the 4-cluster partition than in the 9-cluster one.

Concluding Remarks

This paper introduced a clustering algorithm that is able to partition objects taking into account simultaneously their relational descriptions given by multiple dissimilaritity matrices. These matrices could have been generated using different sets of variables and dissimilarity functions. This algorithm provides a partition and a prototype for each cluster as well as a relevance weight for each dissimilarity matrix by optimizing an adequacy criterion that measures the fit between clusters and their representatives. These relevance weights change at each algorithm iteration and are different from one cluster to another. The usefulness of this algorithm was illustrated comparing the categorization of INRIA research teams given by the clustering algorithm with the a priori expert categorization given by INRIA. The clustering algorithm was able to retrieve the a priori categorization, the observed minor divergences being explained by political choices of INRIA.

Fig. 1 .

 1 Fig. 1. INRIA research categorization and example of the Activity Report summary