Modélisation du mouvement de foules denses :
phénoménologie et couplage de modèles

Observables et validation des modèles

Les études expérimentales dans le domaine des mouvements de foules se multiplient ces dernières années pour identifier et caractériser les différents phénomènes pouvant émerger des comportements individuels. Ces différentes données peuvent servir à nourrir la modélisation, en donnant des valeurs à des paramètres utilisés pour la taille des individus ou la vitesse de marche libre par exemple. Beaucoup d’expériences se concentrent en revanche sur l’identification et la caractérisation de phénomènes émergents propres aux foules en mouvement, comme le Faster is Slower mentionné ci-après. Ces phénomènes peuvent servir à construire de nouveaux modèles ou à comparer des calculs à la réalité dans un cadre contrôlé. La validation de modèles. Il serait tentant de chercher à valider les modèles, en comparant des résultats de simulation à chaque phénomène observé. Il n’existe pas à notre connaissance de cas test globalement accepté par la communauté pour chaque situation. Un document fréquemment utilisé est un guide publié par l’International Maritime Organisation (IMO), qui fournit un ensemble de tests et de données d’entrée pour un logiciel de simulation de mouvement de foules ([IMO07]). 26 Figure 2.5 – Quelques exemples de géométries dans les tests IMO ([IMO07]) Ces tests permettent de vérifier un certain bon fonctionnement du logiciel (en vérifiant qu’on est capable de fixer correctement les paramètres d’entrée, ou qu’il n’y a pas de pénétration de murs solides, etc). Assez peu de critères de validation précis sont cependant donnés : le test n°9 par exemple demande de simuler l’évacuation d’une salle rectangulaire donnée en figure 2.5 et de vérifier qu’en fermant deux des quatre portes le temps d’évacuation double approximativement. Dans le cadre de mouvements de foules cette démarche de validation n’est en général pas possible, les phénomènes observés étant très souvent par nature difficiles à observer systématiquement. Les comportements humains possèdent une part inhérente d’aléatoire qui peut être délicate à estimer et certains processus de pensées complexes peuvent être difficiles à étudier en laboratoire : lors d’un exercice d’évacuation en cas d’incendie en laboratoire, à quel point le comportement observé diffère d’une évacuation réelle ? On ne peut ainsi rejeter ou accepter des modèles selon leur capacité à reproduire des résultats expérimentaux, mais plutôt catégoriser les modèles selon les phénomènes qu’ils peuvent reproduire et décider de la démarche de modélisation la plus adaptée. Vitesse de déplacement et densité. En cas de marche libre des individus on observe des vitesses de marche comprises entre 1.2 m/s et 1.5 m/s. D’après [Wei92], la vitesse de marche libre des individus dans le monde suit une loi normale de moyenne 1.34 m/s, avec un écart-type de 0.26m/s. En cas de congestion en revanche il est connu que la vitesse individuelle diminue. Beaucoup d’études ont été menées afin de déterminer la relation entre vitesse et densité locale. Des expériences dans des couloirs circulaires sont proposées dans [JALP12] et [MGM+12]. Dans le cas de mouvements de panique, les vidéos d’une bousculade à la Hadj en 2006 ont été analysées dans [HJA07], et montrent des vitesses positives mêmes à des densités bien supérieures à ce qui est généralement considéré dans la littérature. Ces différentes observations sont représentées en Figure 2.6. À partir d’observations de terrain, différentes relations ont été proposées dans [Wei92] ou [PM78]. 27 Figure 2.6 – Relations entre densité et vitesse dans [Wei92 ; PM78 ; HJA07]. Les traits pleins représentent des formules théoriques et les points des données expérimentales. Flux en amont de rétrécissements Dans le cadre d’évacuations de bâtiments, le flux à un goulot d’étranglement est un indicateur significatif de la fluidité et la vitesse de l’évacuation. Les portes par exemple sont des structures où on peut s’attendre à des congestions plus importantes. Les expériences menées dans [Cep09] montrent que, pour des fortes densités en amont d’une porte, le flux passant diminue (effet parfois appelé capacity drop). Plusieurs études essayent de caractériser ces effets de congestion et d’en estimer les causes lorsque cela est possible. Effet zipper. L’effet zipper est un phénomène d’auto-organisation lors du passage de couloirs où les piétons s’organisent spontanément en rangées pour optimiser l’espace disponible. Cet effet a été observé dans [HD05] pour des flux unidirectionnels dans différentes largeurs de couloirs. Ces résultats mettent en évidence que le flux augmente par crans successifs, et non linéairement. Dans [PMLW14] la répartition des plus proches voisins dans des expériences de flux unidirectionnels a été étudiée et montre que cette organisation favorise la formation de groupes de deux piétons, qui crée une alternance au niveau des temps de passage entre deux individus en sortie de couloir (deux piétons côte à côte vont sortir pratiquement en même temps, les deux suivants seront après un intervalle plus grand, etc.). Ces corrélations ont été observées dans d’autres expériences menées sur des évacuations d’urgence à travers une porte [NBK17]. Les auteurs ont observé qu’en cas de fortes bousculades le flux pouvait devenir instable et observent que le passage par la porte se fait par clusters de deux ou trois personnes. Cette similarité phénoménologique les conduit à parler d’effet zipper généralisé, qui pourrait s’appliquer dans des situations plus variées

Couplage de modèles 

Dans ce chapitre nous allons nous intéresser aux méthodes de couplage pour des modèles à une dimension. Ces situations simplifiées nous permettrons d’identifier différentes approches du problème qui serviront ensuite pour le passage à deux dimensions. Le problème du couplage entre micro et macro qui nous intéresse consiste à décrire le passage d’entités lagrangiennes en quantité eulérienne appropriée à partir d’une interface fixée en amont, comme l’illustre la figure 3.1. Lors de ce changement d’échelle la modélisation de l’effet de chacun des modèles sur l’autre n’est pas unique et aura un rôle sur la propagation d’information entre les deux modèles. Dans ce chapitre, nous écrivons des conditions préservant la stabilité du flux à l’interface pour des modèles similaires à échelle microscopique et macroscopique, ainsi que des schémas numériques adaptés. Nous abordons le couplage à l’interface pour deux types de modèles. Pour des modèles basés sur l’anticipation individuelle, nous regarderons les différentes formulations possibles pour transmettre x = 0 f(0, t) ρ(x, t) Xi(t) XL(t) x Figure 3.1 – Couplage microscopique vers macroscopique. Le mouvement se fait vers les x positifs et le changement de modèle à l’interface x = 0. 31 l’information au niveau de l’interface, ainsi que des approches numériques pour ce type de problèmes. Nous nous intéresserons ensuite à des modèles raides de type granulaire, en montrant le lien que l’on peut faire entre ces deux approches lorsque l’on considère une anticipation de plus en plus raide. Nous développerons ensuite les schémas numériques pour le couplage de ces modèles. L’organisation de cette section est résumée dans la figure 3.2 Modèle FTL Modèle LWR Modèle d’inhibition Modèle d’inhibition limite anticipation contact microscopique macroscopique section 3.1 section 3.2 sous-section 3.2.2 microscopique macroscopique couplage couplage Figure 3.2 – Organisation des parties de cette section. 3.1 Couplage de modèles d’anticipation Nous allons utiliser les modèles Follow-the-Leader (FTL) et Lighthill-Whithams-Richards (LWR), décrits dans le chapitre 2. Nous rappelons ici les équations principales pour faciliter la lecture du manuscrit. Le modèle Follow-The-Leader, ou FTL, décrit le mouvement d’individus 1, …, L de positions Xi = Xi(t) et de vitesses Vi = Vi(t) pour tout i tel que 1 ≤ i ≤ L par le système d’équations : Vi = φ(Xi+1 − Xi) ∀ 1 ≤ i < L, (3.1) VL = V (t), (3.2) Xi+1(0) − Xi(0) > 2r ∀ 1 ≤ i < L, (3.3) où φ : R + → R + est une fonction croissante à valeurs dans [0, U], U est la vitesse maximale d’un individu et V (t) est une fonction donnée (typiquement continue et à valeurs dans [0, U]) et r 32 est le rayon d’un individu, supposé identique entre tous les individus. On supposera que φ est nulle sur [0, 2r]. Le nombre d’individus L change au cours du temps, par conséquent dans la suite nous utiliserons L(t) la fonction donnant le nombre d’individus dans le modèle microscopique à l’instant t. Le modèle Lighthill-Whithams-Richards (LWR), est basé sur une équation de conservation du premier ordre avec une vitesse dépendant de la densité locale : ∂tρ + ∂x(f(ρ)) = 0, (3.4) où f(ρ) = ρv(ρ), v 7→ v(ρ) est une fonction donnée. Plusieurs formes de fonction ont été proposées pour v(ρ) dans la littérature (voir section 2.3). On suppose généralement que v est décroissante et que pour un certain ρmax on a v(ρmax) = 0. En choisissant v(ρ) = φ(1/ρ), on peut interpréter le modèle FTL comme une formulation lagrangienne du modèle LWR. En effet, il est possible de prouver qu’en choisissant une limite appropriée du modèle FTL lorsque le nombre de particules tend vers +∞ on retrouve le modèle LWR (voir [DR15] pour une formulation précise de ce résultat). Dans [HR17] il est aussi montré comment le modèle FTL peut être vu comme un schéma numérique lagrangien pour le modèle LWR. Discrétisation et schémas numériques. Le modèle FTL est discrétisé avec un schéma implicite en temps pour l’équation (3.1). Pour un pas de temps ∆t, on note Xn i l’approximation de Xi(t) à l’instant n∆t et on a : X n+1 i = Xn i + ∆tφ(X n+1 i+1 − X n+1 i ) ∀ 1 ≤ i < Ln , (3.5) où L n est le leader à l’instant n∆t. La résolution numérique se fait à l’aide d’une méthode de Newton qui ne pose pas de problème particulier. Ce type de schéma a l’avantage d’être stable, même pour des grands pas de temps. Ce choix permettra d’utiliser un pas de temps commun entre le modèle micro et macro. Pour le modèle LWR, l’espace est discrétisé avec un pas ∆x sur une grille xj = (j+ 1 2 )∆x ∀j ∈ Z. On utilise un schéma numérique de type volumes finis avec un pas de temps ∆t pour toute cellule j : ρ n+1 j − ρ n j ∆t + F n j+1/2 − F n j−1/2 ∆x = 0, (3.6) où ρ n j est l’approximation de ρ(n∆t, xj ) et F n j+1 le flux à l’interface entre les cellules xj et xj+1. Pour notre schéma numérique nous prenons le flux introduit dans [FFL+16] : F n j+1/2 = F(ρ n i , ρn j+1) = ρ n i v(ρ n j+1). (3.7) Ce schéma repose sur le fait que la matière se propage vers l’aval, alors que l’information se propage vers l’amont. On a F(ρ, ρ) = ρv(ρ) donc ce schéma est consistant, et il est monotone sous la condition CFL : ∆t < ∆x/(U + ρmax max ρ (v ′ (ρ))). (3.8) On peut le montrer en réécrivant le schéma sous la forme : ρ n+1 j = G(ρ n j−1 , ρn j , ρn j+1), 33 la monotonicité est vérifiée lorsque G(ρ1, ρ2, ρ3) est croissante par rapport à chacune de ses variables. On a : G(ρ1, ρ2, ρ3) = ρ2 − ∆t ∆x (ρ2v(ρ3) − ρ1v(ρ2)) v est décroissante, donc G est croissante par rapport à ρ1 et ρ3. En dérivant par rapport à son second argument on a : ∂ρ2G(ρ1, ρ2, ρ3) = 1 − ∆t ∆x (v(ρ3) − ρ1v ′ (ρ2)). La condition CFL assure ainsi la monotonicité du schéma.

Couplage microscopique vers macroscopique : premiers schémas heuristiques

Dans le cadre du couplage microscopique vers macroscopique, nous devons attribuer d’une part une vitesse au leader de la partie microscopique et un flux pour le modèle macroscopique. On supposera par la suite que l’interface entre les deux modèles se situe en x = 0. Pour la vitesse du leader un premier choix de modélisation possible est de prendre la vitesse macroscopique à l’interface : VL(t) = v(ρ(t, 0 +)), (3.9) où L(t) donne l’indice du leader du modèle microscopique à l’instant t et VL(t) = X˙ L(t) . Le couplage est complété en attribuant un flux à l’interface pour le modèle macroscopique. Le leader ayant une vitesse donnée, en lui attribuant une densité ρb(t) observée par le modèle macroscopique à l’interface le flux pourra être exprimé : f(t, 0 +) = ρb(t)VL(t) (3.10) Une première possibilité de modélisation consiste à écrire : ρb(t) = ( 1/2r si − 2r ≤ XL(t) ≤ 0 0 sinon, (3.11) La densité vue au bord est ainsi concentrée autour du leader qui s’apprête à passer dans le modèle macro. En rassemblant toutes les équations concernées, le système d’équations pour le problème de couplage micro-macro est donné par : ∂tρ + ∂x(f(ρ)) = 0 ∀(x, t) ∈ R + × R + (3.12) Vi = φ(Xi+1 − Xi) ∀ 1 ≤ i < L(t), (3.13) VL(t) = v(ρ(t, 0 +)), (3.14) L(t) = sup {i, Xi < 0} , (3.15) ρb(t) = ( 1/2r si − 2r ≤ XL(t) ≤ 0 0 sinon, (3.16) f(t, 0 +) = ρb(t)VL(t). (3.17) 34 Schémas numériques et cas test. Nous reprenons ici la discrétisation des équations (3.1) et (3.4) introduites avec les schémas numériques (3.5) et (3.6). Pour le modèle LWR on a un maillage xj = (j + 1 2 )∆x ∀j ∈ N. Pour la vitesse du leader on a : X n+1 Ln − Xn Ln ∆t = v(ρ n 0 ). (3.18) Une difficulté réside ici dans l’approximation du flux F n −1/2 à l’interface. On peut exprimer ce flux en utilisant la conservation de la masse entre les deux modèles. À chaque pas de temps on peut déterminer la quantité intégrée de piéton de chaque côté de l’interface. On note m[X] la masse intégrée à gauche de l’interface telle que : m[X] = Z 0 −∞ ρL(t) [X](x)dx. (3.19) Pour avoir une conservation de la masse entre les deux modèles, il est nécessaire que la perte de masse dans le modèle lagrangien soit compensée par un ajout équivalent dans le modèle macroscopique. La perte de masse dans le modèle microscopique est donnée par : F n −1/2 = m [Xn ] − m  Xn+1 ∆t , (3.20) où Xn = (Xn i )i . Si les positions X(n∆t) sont connues, m(t) peut être calculée directement sans approximation supplémentaire. Entre deux pas de temps on peut ainsi avoir un ou plusieurs piétons qui passent l’interface. Le schéma numérique global devient alors : X n+1 i = Xn i + ∆tφ(X n+1 i+1 − X n+1 i ) ∀ 1 ≤ i < Ln , (3.21) L n = sup {i, Xn i < 0} , (3.22) ρ n+1 j − ρ n j ∆t + F n j+1/2 − F n j−1/2 ∆x = 0 ∀j ∈ N, (3.23) F n j+1/2 = ρ n j v(ρ n j+1), (3.24) X n+1 Ln − Xn Ln ∆t = v(ρ n 0 ), (3.25) F n −1/2 = m [Xn ] − m  Xn+1 ∆t . (3.26) Cas test. Nous présentons ici à deux cas tests : • le premier en figure 3.3 montre une portion de route de 100 m avec 90 individus à gauche répartis uniformément sur 50 m et à droite et une densité constante à 1.8 pers/m. À 100 m la densité au bord est fixée à ρmax. • le second cas test en figure 3.4 avec 40 individus à gauche et une densité à 0.9 pers/m à droite. Des oscillations apparaissent au niveau de l’interface entre les deux modèles. Ces oscillations sont d’autant plus marquées à basse densité par rapport aux hautes densités. Le cas test en figure 3.3 montre en revanche la capacité du couplage à faire remonter l’information à l’interface. Ces oscillations sont dues au choix de représentation pour la densité au bord vue par le modèle macroscopique (3.11) : en choisissant une fonction aussi locale, le flux à l’interface en est impacté et oscille fortement.

Modification de l’absorption

Dans la littérature ([DR15] par exemple) on peut voir que les outils de convergence du modèle micro vers le modèle macroscopique reposent sur la définition d’une densité pour chaque individu telle que : ρi(t, x) = ( 1 Xi+1−Xi si Xi ≤ x < Xi+1 0 sinon. (3.27) Pour le leader cette définition n’est pas applicable directement : celui-ci ne perçoit aucun individu devant lui. Une première solution serait de choisir pour la densité ρb(t) : ρb(t) = 1  XL(t)   + 2r (3.28) qui correspond à la densité que voit le leader si un piéton imaginaire se tient au niveau de l’interface. Avec cette définition cependant la densité augmente toujours au fur et à mesure qu’un individu s’approche de l’interface, créant de nouvelles oscillations. 36 paramètre valeur v(ρ) U(1 − ρ ρmax ) U 1.2m/s ρmax 2.5m−1 ∆x 6.25m ∆t 1.0s Table 3.1 – Paramètres pour le cas test considéré en Figures 3.3-3.4. On choisit ainsi une autre fonction pour la densité du leader. On introduit la fonction τ donnant le temps d’intersection du précédent leader avec l’interface, c’est-à-dire : τ (t) = ( 0 si L(t) = L(0) sup  t, XL(t)+1 < 0 sinon. (3.29) Dans cette définition on prend comme convention qu’à l’instant 0 un piéton imaginaire se situe à l’interface. On écrit enfin la densité vue par le modèle macroscopique comme étant : ρb(t) = 1  XL(τ) (τ )   . (3.30) On fixe ainsi la valeur de la densité au moment où il devient leader. En approchant de l’interface on aura alors une injection continue à la vitesse v(ρ(t, 0 +)). Notre nouvelle version du couplage s’obtient donc en remplaçant l’équation (3.16) par l’équation (3.30). Modification du schéma numérique. Le schéma numérique doit prendre en compte les temps τ où le calcul de la vitesse du leader change, et qui peuvent se trouver entre deux itérations. Nous modifions ainsi le schéma numérique en introduisant une détection des instants τ où un piéton traverse l’interface. Nous allons ajouter entre deux itérations une étape de vérification pour savoir si le leader traverse l’interface, et si c’est le cas on fait avancer le temps jusqu’au moment où le leader traverse l’interface. On obtient ainsi la densité au bord pour le leader suivant et on continue l’itération, comme illustré en figure 3.5. On aura ainsi autant de sous-pas de temps qu’il y aura de leaders qui franchissent l’interface pendant ∆t. Ce sous-pas de temps est tel que : ∆tk = min    ∆t − X j¯ 0 de l’interface et de piéton uniformément répartis avec la même densité à gauche de l’interface à un pas n, k quelconque. En supposant toujours que φ(1/ρ) = v(ρ), la vitesse des agents microscopiques sera égale à la vitesse locale de la densité macroscopique. À l’interface le flux vaut ρbv(¯ρ) = ¯ρv(¯ρ) si la densité ρb a été initialisée correctement, donc le leader avancera aussi à vitesse constante et le flux à l’interface sera constant lui aussi. Lorsque le leader rencontre l’interface, la nouvelle densité ρb devenant la condition de bord sera aussi ρ¯ grâce à la stabilité du modèle FTL.