PhD Philippe Pecol
email: philippe.pecol@enpc.fr

Stefano Dal Pont
email: dalpont@lcpc.fr

Silvano Erlicher
email: erlicher@enpc.fr

Joanna Bodgi
email: joanna.bodgi@fs.usj.edu.lb

Professor Pierre Argoul
email: pierre.argoul@enpc.fr

A 2D discrete model for crowd-structure interaction

Keywords: granular assembly, contact, crowd movement, crowd-structure interaction, synchronization, vibrations, resonance. 1

An emerging problem for civil engineering in the field of dynamics is the modelling of dynamic effects due to crowdstructure interaction. We are particularly concerned with structures such as footbridges, which oscillate due to the crossing of a group of pedestrians. The objective of this study is double: to model the movement of pedestrians with consideration of pedestrian-pedestrian and pedestrian-obstacle interactions and to include a pedestrian-structure coupling in the proposed model. This 2D discrete model would be able to study the problem of synchronization between pedestrians and structure on footbridges. Our idea is to extend the modelling of particles movements to study the movements of a crowd. The non-smooth discrete model of Frémond initially proposed to simulate a granular assembly, is chosen. It applies a rigorous thermodynamic framework in which the local interactions between particles are managed by the use of dissipation pseudo-potentials. Social forces as well as a desired direction/ velocity are introduced in order to simulate the behaviour of pedestrians. Concerning the pedestrian-structure coupling, a differential equation of Kuramoto allows one to manage the evolution of the pedestrians' phase. Two cases are studied depending on the sensitivity of pedestrians to the footbridge's oscillations. Numerical simulations on the Millenium Bridge are performed and discussed.

Introduction

Several footbridges, built recently in accordance with the architectural criteria of lightness and elegance, such as the Millennium Footbridge in London and the Solferino footbridge in Paris have proved sensitive to the excitation induced by the crossing of pedestrians. Several experimental measurement campaigns have allowed a better understanding of this phenomenon: the crowd walking on a bridge imposes a lateral excitation on the structure which has a frequency of about 1Hz. When a lateral mode of vibration of the bridge, often the first, has a frequency close to this value, a resonance phenomenon is activated. It follows that the amplitude of oscillations of the bridge increases. If the number of pedestrians is not high, the amplitude of oscillations stays small and pedestrians continue to walk as they would on a rigid floor. If instead the number is high enough (beyond a certain critical number), the oscillations' amplitude becomes larger, enough so that pedestrians are encouraged to change their way of walking, in particular their frequency of walking, until there is a phase synchronization between pedestrians and structure. The behaviour of a pedestrian is also influenced by the crowd around him: if the density of the crowd is very low, the walking is "free", i.e. each individual walks as if there was nobody close to him. Conversely, if the crowd is dense, the pedestrian is forced to walk "at the same speed" as the others.

In this communication, the first part deals with the crowd movement microscopic model. To describe the acceleration of each pedestrian along his trajectory, we use the non-smooth discrete model of Frémond [START_REF] Fremond | Rigid bodies collisions[END_REF][START_REF] Fremond | Collisions[END_REF], proposed to simulate a granular assembly, in order to manage the local pedestrian-pedestrian and pedestrian-obstacle interactions. The adaptation of the Frémond's model to the crowd gives a desired trajectory to each pedestrian. To take into account the pedestrian's oscillations around his trajectory due to his walking and his action on the bridge, an alternating (sinusoidal) sideways force is used. This force allows one to define the acceleration of each pedestrian's oscillations around his trajectory. In the second part, we give the crowd-footbridge interaction model in the case of lateral oscillations of a bridge. The synchronization of the walking frequency of each pedestrian with the lateral oscillations frequency of the bridge can be managed thanks to a differential equation of Kuramoto [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF][START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF], which allows one to govern the evolution of the total phase of the walking force generated by each pedestrian on the bridge. In the last part, the whole model is applied to the case of the North span of the Millennium Bridge and results are presented and discussed.

Representation of the crowd

We need to model the dynamic loading of a single pedestrian on a vibrating structure like footbridges in order to study the phenomenon of crowd-structure synchronization. The human body is a very complex mechanical system made of several parts in mutual interaction. It can be modeled according to the desired degree of refinement and objectives to be achieved. In [START_REF]Stability, scaling and chaos in passive-dynamics gait models[END_REF], each part of the human body's skeleton is represented by a rigid body connected to the others by means of springs, dampers, and torques which impose a relative movement between rigid bodies in order to reproduce the walking. In [START_REF] Nigg | The effect of the muscle stiffness and damping on simulated impact force peaks during running[END_REF][START_REF] Zadpoor | A model-based parametric study of impact force during running[END_REF], the human body is represented with a set of masses, springs and dampers, in order to reproduce the vertical force of foot's impact on the ground during a footstep.

The human body can also be represented by a simple mass interacting with the ground. The expression of the restoring force between the pedestrian and the ground allows one to distinguish between models. Erlicher and al. [START_REF] Erlicher | Modeling the lateral pedestrian force on a rigid floor by a self-sustained oscillator[END_REF] suggest a self-sustained single degree of freedom oscillator, able to accurately predict the lateral walking force. In [START_REF] Abrams | Two coupled oscillator models: the Millenium bridge and the chimera state[END_REF][START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF], the action of a pedestrian on a bridge is represented by a sinusoidal lateral force with amplitude of 30N and the evolution of its total phase is managed by a differential Kuramoto equation. Bodgi [START_REF] Bodgi | Lateral vibration of footbridges under crowd-loading: continuous crowd modelling approach[END_REF][START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF] was inspired by this equation to achieve an original macroscopic model of pedestrian-structure coupling.

The knowledge of the trajectory of each pedestrian and of his action on the ground is sufficient for us to conduct our study on the synchronization. Hence, we represent a pedestrian i by a mass i m , supposed lumped at the gravity centre of the pedestrian, and subjected to the inertial, gravitational and restoring forces. In the Galilean reference frame ( )

3 2 1 e e e , , , 0
, the fundamental principle of dynamics applied to each pedestrian i gives:

g q F q u i i i abs i i i m t t m = + ) , ( ) , ( & & (1)
where the dot indicates the derivative function of time; g is the vector of the gravity's acceleration; abs i u & & is the absolute acceleration of i m ; i q is the position of the pedestrian i in the Galilean reference frame; and ( )

t , i i q F
is the exchanged force between the pedestrian i and the ground.

The absolute acceleration is expressed as: , taking into account local pedestrian-pedestrian and pedestrian-obstacle interactions. The microscopic model of crowd presented in [START_REF] Pecol | Modelling crowd-structure interaction[END_REF][START_REF] Pecol | Discrete approaches for crowd movement modelling[END_REF] is used to obtain

) ( ) ( ) ( ) , ( ) ( ) ( ) ,
) (t tr i u & &
. This model has been first proposed by Frémond [START_REF] Fremond | Rigid bodies collisions[END_REF][START_REF] Fremond | Collisions[END_REF] to model the rigid particles' movement, and is then adapted to the crowd movement.

The model of Frémond, presented in the first part of this section, is based on the theory of rigid bodies' collisions and follows the lead of the works of Moreau [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF]. The numerical aspects were later developed in [START_REF] Dal | A theory for multiple collisions of rigid solids and numerical simulation of granular flow[END_REF][START_REF] Dal | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF]. In this model, the set of N particles is considered as a deformable system composed of rigid solids. 3) is applied almost everywhere except at the instant of the collision where it is replaced by Equation ( 4). When contact is detected, the velocities of colliding particles are discontinuous and we introduce in Equation ( 4) the percussions int p and ext p , interior and exterior to the system respectively. By definition, percussions have the dimension of a force multiplied by a time. The int p percussions are unknown; they take into account the dissipative interactions between the colliding particles (dissipative percussions), and the reaction forces in order to avoid overlapping among particles (reactive percussions). Frémond [START_REF] Fremond | Rigid bodies collisions[END_REF][START_REF] Fremond | Collisions[END_REF] showed that these interior percussions are defined in duality with the velocity of deformation at the moment of impact ( )

p p u u M f f u M + - = - + - = (3) (4) Equation (
) ( ) ( - + + tr, tr, u ∆ u ∆ & & 2
, in the sense of the work of internal forces. int p depends on (

)

2 ) ( ) ( - + + tr, tr, u ∆ u ∆ & &
, and Frémond used a pseudopotential of dissipation Φ (a convex function [START_REF] Moreau | Sur les lois du frottement, de la viscosité et de la plasticité[END_REF]) defined as:

reac diss Φ + Φ = Φ
, where diss Φ and reac Φ are two pseudopotentials which define the dissipative and reactive interior percussions respectively, to express these quantities:

        + Φ ∂ ∈ - + 2 ) ( ) ( tr, tr, int u ∆ u ∆ p & & (5)
where the operator ∂ is the subdifferential which generalizes the derivative for convex functions [START_REF] Fremond | Collisions[END_REF]. The pseudopotential diss Φ is chosen to be quadratic:

( ) ( ) ( ) 2 2 ) ( ) ( 2 - + + = Φ tr, tr, u ∆ u ∆ & & K diss
where K is a coefficient of dissipation. This choice allows one to find the classical results of the coefficient of restitution. Others choices of diss Φ allow one to obtain a large variety of behaviours after impact [START_REF] Cholet | Chocs de solides rigides[END_REF][START_REF] Fremond | Collisions[END_REF].

In Equation ( 4), the problem is to find the velocity after collision of grains

( ) ( ) [ ] w M p M u w ∆ w M w w v ext 1 tr, t t - -+ - Φ + ℜ ∈ = 2 ) ( min arg 2 N (6)
where the solution is

( ) 2 - + + = tr, tr, u u v & & & .
The existence and uniqueness of the solution are proved in [START_REF] Dal | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF][START_REF] Fremond | Rigid bodies collisions[END_REF][START_REF] Fremond | Collisions[END_REF] and the minimization problem ( 6) is numerically solved using a time-stepping scheme and the classical Uzawa algorithm [START_REF] Dal | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF].

So, when there is no contact, Equation (3) gives directly the real value of

) (t tr u & &
. When a contact is detected, Equation (4) allows to get an approximated value of

) (t tr u & & after the determination of + tr, u & : h t - + - ≈ tr, tr, tr u u u & & & & ) ( (7)
where h is the chosen time step of the numerical scheme.

In the second part of this section, we show how the model of Frémond initially applied to rigid grains' movement can be adapted to the crowd movement. The first step is to give a desired trajectory to each particle. The strategy of the shortest path to get from one point to another [START_REF] Kimmel R | Fast marching methods for computing distance maps and shortest paths[END_REF] is implemented through a Fast Marching algorithm and is used to obtain the desired direction, tr i d, e , of an individual i . This direction depends on the evolution space (obstacles, etc.), the time and also the characteristics of the individual (gender, age, hurried steps or not, etc.). It is defined by: of the pedestrian i and on the lateral angular frequency i ω of his free desired walking [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF]:

i i step L ω π , 0 , = tr i d, u & (8)
In a second step, we introduce the desired velocity of each pedestrian i in the model of grains' movement to obtain the one of crowd movement. Let

) ( ) ( t t ext a ext f f =
, where the external acceleration force [START_REF] Helbing | Social force model for pedestrian dynamics[END_REF] allows one to give a desired direction and amplitude of the velocity to each pedestrian. Each component where i τ is a relaxation time, allowing the pedestrian to recover the desired velocity after a contact. Smaller values of i τ let the pedestrians walk more aggressively [START_REF] Helbing | Simulating dynamic features of escape panic[END_REF][START_REF] Pecol | Discrete approaches for crowd movement modelling[END_REF]. The pedestrians' behavior can be enriched with the adding of other external social forces [START_REF] Helbing | Simulating dynamic features of escape panic[END_REF][START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF] in

) (t ext a f
) (
) (t ext f
, in order to become more realistic (socio-psychological force, attractive force, group force, etc.).

Expression of

) (t osc i u & &
When a pedestrian is walking, its center of gravity oscillates around its mean trajectory. Thus, to define the acceleration of the pedestrian's oscillations around his trajectory osc i u & & , we consider the case where the pedestrian is walking in place.

We assume that the forces osc i F applied by the two legs of the pedestrian i on the floor are identical. In the frame ( )

2
1 e e , , 0

, Equation (1) becomes:

0 F u osc i osc i = + ) ( ) ( t t m i & & (10)
As each pedestrian is walking in a plane, the force applied by the pedestrian i on the floor has a normal osc i N, F and a tangential osc i T, F component to the pedestrian's motion. We assume that each pedestrian imparts an alternating (sinusoidal) sideways force to the bridge [START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF]. As the frequency of osc i T, F is usually taken to be twice the frequency of osc i N, F [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF], the expression of osc i F can be written as follows: A combination of [START_REF] Fremond | Collisions[END_REF] and [START_REF]Stability, scaling and chaos in passive-dynamics gait models[END_REF] allows us to determine

( ) ( ) ) ( ) ( sin ) ( ) ( 2 sin ) ( ) ( ) ( ) ( ) ( , , t t N t t T t t F t t F t i i i i osc i N osc i T ⊥ ⊥ + = + =
) (t osc i u & & : ( ) ( ) (
)

) ( ) ( sin ) ( ) ( 2 sin 1 ) ( t t N t t T m t i i i i i ⊥ + - =

Crowd-footbridge interaction

We consider N pedestrians walking on a bridge of length L and width l . Some assumptions are done: the floor is flat and horizontal; the main axis (longitudinal) of the bridge is straight, along 1 e ; the floor's movement is governed by its lateral oscillations, i.e. the floor's oscillations take place in the horizontal plane while the longitudinal oscillations are negligible; the displacements of the floor are supposed constant along 2 e , unit direction vector of the lateral footbridge's oscillations; and the lateral displacements of the floor are small. Thus, we define the following Cartesian coordinate system with x-axis which is parallel to the main axis, the z-axis vertical and the y-axis oriented accordingly; the structure has a geometrically linear behaviour; and the displacement of the footbridge is given by:

2 e q U ) , ( ) , ( t x U t y = .
The response of a linear structure can be expanded into a sum of modal responses. For the lateral component, limiting to the m N first modes, we find:

∑ = = m N j j y j y t U x t x U 1 , ) ( ) ( ) , ( ψ (13)
where ) (x j ψ is the th j modal shape and ) ( , t U j y is the th j modal coordinate of the bridge's displacement.

As the lateral excitation of the pedestrian due to his walking has a frequency closed to 0.9-1 Hz, the footbridge's lateral mode having the nearest frequency to the pedestrian's excitation one, is sufficient to model the dynamics of the footbridge: The equation describing the dynamic behaviour of the footbridge projected onto its first mode shape is: The expression of the ground's acceleration (in Equation ( 2)) is: 16) states that the bridge moves in the lateral direction only (

) ( ) ( ) ,
) , ( ) ( ) ( ) ( ) ( , 1 1 t x F x t U K t U C t U M i i
2 ent e u ) ( ) ( ) , ( 1 t U x t x y & & & & ψ = (16) Equation (
) , ( t x ent u & & is parallel to 2 e
), with a spatial distribution proportional to the lateral mode shape

) ( 1 x
ψ . This assumption is acceptable for numerous practical situations. More complexes bridge motions could be considered, but this is beyond the purposes of this paper. Synchronization can be seen as an adaptation of the frequency of the force generated by a pedestrian to the frequency of the "crowd-structure" system. We choose to study two differential Kuramoto equations [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF][START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF] to govern the evolution of the total phase ) , ( t x i i φ (see Equation ( 11)) of the walking force generated by the pedestrian i on the bridge. Both differential equations allow the instantaneous walking angular frequency of the pedestrian to converge to that of the footbridge.

In the first one, pedestrians are sensitive to the amplitude of the footbridge's displacement:

( ) α φ ψ ψ ε ω φ + - + = ∂ ∂ ) , ( ) ( sin ) ( ) ( ) , ( 1 t x t x t A t t x i dis i str i y dis i i i dis i (17)
as in [START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF] but their instantaneous frequency also depends on the footbridge's modal shape. In Equation ( 17), dis i ε quantifies pedestrians' sensitivity to the bridge's oscillations, its value was estimated by Strogatz [START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF] by comparing its results (without modal shape) with the experimental results achieved on the Millennium bridge [START_REF] Dallard | The millennium bridge Londonproblems and solutions[END_REF]. Strogatz chose 2 π α = to have the lateral force exerted by pedestrian i on the ground in phase with the lateral velocity of the bridge when synchronization occurs.

In the second one, pedestrians are sensitive to the amplitude of the footbridge's acceleration [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF]: where acc i ε quantifies pedestrians' sensitivity to the bridge's oscillations, its value was estimated by Bodgi [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF] by comparing its results (with modal shape) with the experimental results achieved on the Millennium bridge [START_REF] Dallard | The millennium bridge Londonproblems and solutions[END_REF];

( ) ( ) α φ ψ ψ ψ ε ω φ + - + = ∂ ∂ ) , ( ) ( sin ) ( ) ( ) ( 2 ) , ( 2 1
) ( , t A m y is the maximum value of ) (t A y
during the last 5 seconds of simulation;

( ) 2 1 , ) ( ) ( ) ( ) , ( t x t A t x a str i m y i foo ψ ψ & =
is the footbridge acceleration; and Bodgi has found that

2 π α = too.
For a walking on a rigid floor, both of previous equations are reduced to:

i i t ω φ = ) ( & .
Hence, the expression of

) , ( , t x F i i y becomes: ( ) ( ) ( ) (
) To take into account the oscillations of the structure in the crowd movement model, the amplitude of the desired velocity of the pedestrian i in Equation ( 8

) ( cos ) , ( sin ) ( sin ) , ( 2 sin )
( ) ( ) ( ) , ( , 1 , t t x N t t x T t u m t U x m t x F i i i i i i i i tr i y i y i i i i y θ φ θ φ ψ + + - - = & & & & (19)
i i i step i ∂ ∂ = ) , ( ) , ( ) , ( , mod φ π tr
i d, u & (20)
Three expressions of the instantaneous step length can be proposed according to the model we chose to apply.

If we don't take into account the deceleration of the pedestrian, then:

i step i i step L t x L , 0 , , ) , ( = (21
) If we take into account the deceleration and Equation ( 17), so as in [START_REF] Pecol | Modelling crowd-structure interaction[END_REF]:

        - = i i y i step i dis i step A x t A L t x L max, 1 , 0 , , ) ( ) ( 1 , 0 max ) , ( ψ (22)
where i A max, is the maximum amplitude of the lateral displacement of the footbridge that the pedestrian i can stand. When

) ( ) ( 1 i y x t A ψ reaches i A max,
, the amplitude is so large that the pedestrian stops walking.

If we take into account the deceleration and Equation [START_REF] Nigg | The effect of the muscle stiffness and damping on simulated impact force peaks during running[END_REF], so as in [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF]: where min a (resp. max a ) is the amplitude of the lateral footbridge's acceleration from which pedestrians feel the oscillations (resp. stop to walk).

       ≤ < < - - ≤ = ) , ( 0 ) , ( ) , ( ) , ( 1 ) , ( max

Simulations

The experimental study made by Arup [START_REF] Dallard | The millennium bridge Londonproblems and solutions[END_REF] on the North span of the Millennium Bridge showed that when the average walking frequency of pedestrians is equal to the modal frequency of the structure (1.03 Hz), the critical number of pedestrians Nc, above which the frequency synchronization between a portion of the crowd and the structure is triggered, is: Nc = 166 pedestrians. We have tried to find this result with our proposed model of interaction.

In this section, we consider the case where pedestrians cross the North span of the Millennium Bridge in London and we present the results obtained. For each simulation, the number of pedestrians on the structure is constant and all the pedestrians move in the same direction. To keep the number of pedestrians walking on the bridge constant, once one pedestrian exits the bridge on one side, one similar pedestrian enters the bridge with a random position along the y-axis on the other side. The parameters used in simulations are listed in 2) without modal shape : 1 We performed simulations of 10 minutes long considering the following numbers of pedestrians: N = 100, 150, 170, 180, 210, 250, 300 and 400. For each N , we study the influence of the modal shape and of the pedestrians' deceleration (see Equations [START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF][START_REF] Zadpoor | A model-based parametric study of impact force during running[END_REF](23)) on the triggering of synchronization. In order to study the influence of the modal shape on the results obtained by Equations ( 17) and ( 18), we consider the case where the modal shape is not taken into account ( Figure 1 shows examples of non-synchronization and synchronization results, obtained in the following case: with Equation ( 17), without modal shape and without deceleration. Figures 1 a. andb. present a case of non synchronization, obtained with 100 pedestrians crossing the structure. The lateral displacement of the bridge is very low and pedestrians do not synchronize their walking frequency with that of the "pedestrian-bridge" system. Figures 1 c. andd. show a case of synchronization, obtained with 300 pedestrians crossing the structure. The lateral displacement of the bridge reaches 0.15 m, and pedestrians synchronize their walking frequency with that of the "pedestrian-bridge" system. Table 2 shows the preliminary results obtained on the phenomenon of synchronization. It gives the number of simulations where the phenomenon of synchronization is detected compared with the number of achieved simulations. We find that, when the modal shape is not taken into account ( 1) ( 1 = x ψ ), i.e. the pedestrian close to the bounds of the span is feeling the same oscillations as that at mid-span, synchronization is highlighted. The number of pedestrians being synchronized is higher, and therefore the amplitude of the force generated by pedestrians too.

We consider first the case without modal shape ( 1) ( 1 = x ψ ), and with or without deceleration. The use of Equation [START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF], where pedestrians are sensitive to the amplitude of the footbridge's displacement, gives the critical number obtained by Abrams [START_REF] Abrams | Two coupled oscillator models: the Millenium bridge and the chimera state[END_REF], around 150 pedestrians. Using Equation [START_REF] Nigg | The effect of the muscle stiffness and damping on simulated impact force peaks during running[END_REF], where pedestrians are sensitive to the amplitude of the footbridge's acceleration, gives a critical number less than 100 pedestrians. The estimated value of acc i ε given in [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF] with modal shape and with deceleration is not correct for the case without modal shape.

We consider then the case with modal shape ( ( )

L x x π ψ sin ) ( 1 =
). With deceleration, the use of Equation [START_REF] Nigg | The effect of the muscle stiffness and damping on simulated impact force peaks during running[END_REF] gives the critical number obtained by Bodgi [START_REF] Bodgi | Synchronisation piétons-structure: Application aux vibrations des passerelles souples[END_REF], around 170 pedestrians. With or without deceleration, using Equation [START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF] gives a critical number around 250 pedestrians. The estimated value of dis i ε given in [START_REF] Strogatz | Theoretical mechanics: Crowd synchrony on the Millenium bridge[END_REF] without modal shape and without deceleration is not correct for the case with modal shape.

This application to the interaction model on the North span of the Millennium Bridge shows that the value of the parameter which quantify pedestrians' sensitivity to the bridge's vibrations, has a significant influence on the critical number. The value of this parameter has to be chosen according to the differential Kuramoto equation and to the choice of taking into account the modal shape or not.

Conclusion

This paper presents a 2D discrete model for crowd movement on rigid ground obtained from the adaptation of a granular model based on the theory of rigid bodies'collisions developed by Frémond. This crowd movement model is then used in the crowd-structure interaction model. A differential Kuramoto equation allows to manage the evolution of the total phase of the walking force generated by each pedestrian. Depending on the sensitivity of the pedestrians, two differential equations of Kuramoto have been studied.

Numerical simulations on the North span of the Millennium Bridge show that the value of the parameter which quantifies pedestrians' sensitivity to the bridge's oscillations, has to be estimated according to the differential Kuramoto equation and to the modal shape. Other numerical simulations are in progress to fully validate the model from other experimental data. Two studies are planned: first, determining the influence on the phenomenon of synchronization, of the phase shift value between the lateral force exerted by pedestrians on the ground and the lateral velocity of the bridge, and second, comparing the stationary amplitude values of the lateral bridge's displacement when synchronization occurs, with those found in literature.

  is the acceleration of the pedestrian i along is trajectory on the ground; and osc i u & & is the acceleration of the pedestrian's oscillations around his trajectory.The Coriolis acceleration for each pedestrian i is neglected: which are the components of the relative acceleration of the pedestrian i are given. In this purpose, the ground is considered rigid and the reference frame becomes ( ) of each pedestrian i along his mean trajectory on the ground, we control the movements of the N pedestrians in the plane ( )

  we have to solve the following constrained minimization problem:

  is the desired velocity of the pedestrian i . The amplitude tr i d, u & of the desired velocity represents the speed at which the pedestrian i wants to move on the considered structure, it can be influenced by the nervousness of the pedestrian. It depends on the desired step length

  the unit direction vector of the pedestrian's motion (resp. normal to the pedestrian's motion), phase of the walking force generated by the pedestrian i on the bridge, such that is not influenced by external excitations such as the bridge's oscillations.

  C and str K are the modal (generalized) mass, damping and stiffness of the lateral first mode of the footbridge, respectively and where

1 e

 1 of the longitudinal axis of the bridge.

  of a simply supported beam at both ends)

1 =.

 1 For each studied case, several simulations are achieved with different set of parameters. When the final seconds of simulation correspond to the stationnary state, synchronization has been reached.

Fig. 1 :

 1 Fig. 1: Phenomenon of synchronization for the model using Equation (17), without modal shape and without deceleration. Figures a. and b. (resp. c and d) show the results obtained with 100 pedestrians (resp. 300 pedestrians) -a. and c. represent the time evolution of the lateral displacement of the bridge; b. and d. represent the time evolution of the walking frequency of pedestrians on the structure. The red curve in bold represents the instantaneous frequency of the "pedestrian-bridge" system.

.

  Several cases are considered according to N , the type of the differential Kuramoto equation governing the evolution of ) , ( t x i i φ , the expressions of the modal shape and of the instantaneous step length number of simulations: achieved / synchronized pedestrians are sensitive to the amplitude of the footbridge's displacement (Equation (17))pedestrians are sensitive to the amplitude of the footbridge's acceleration (Equation (18))

cor i osc i tr i i ent i cor i rel i i ent i i abs i u u u q u

	u	( q	,	t	)	=	u	( q	t	+	u	t	+	u	t
						=			t	+		t	+		t	+	t

ext int tr, i tr, i ext int tr i

  Each particle i is rigid, circular, of radius

	i r and of mass i m . We introduce : M , the dimension N 2 , ),..., N 2 2 × ( N )) ( ), ( ( ) ( t t t t tr N t tr 2 t tr 1 t tr t u u u )) ( ( t tr ij u ∆ & , the relative deformation velocity between the particles i and j , mass matrix of the particles; ) (t tr u & , the velocity vector of ) ( ) ( )) ( ( t t t tr j tr i tr ij u u u ∆ & & & -= . The motion equations of the system are: everywhere ) ( ) ( )) ( ) ( ( t t t t u applied to the system; and everywhere almost ) ( ) ( ) ( t t t

  t

																f	ext i a,	t (	)	=	i m	u &	tr d,	i	e	tr d,	i	i t τ ) (	-	u &	tr i	t (	)	(9)
	dimension N 2 :	t	f	ext a	( t	)	=	(	t	f	ext a,1	( t	),	t	f	ext a,2	( t	)...,	t	f	ext N a,	( t	)	)	f	ext i a,	of the force vector of

, is associated with pedestrian i and can be expressed as: th INTERNATIONAL CONFERENCE

  (

	) normalized by: (t str ψ and an instantaneous amplitude	) A y , and (t	ψ	( 1 x	)	is the first lateral modal shape of the footbridge,
	where	U	y	( t	)	=	A y	( t	)	sin( ψ	str	( t	))	1 is the lateral modal displacement of the footbridge having an instantaneous phase t U x t x U y ψ ≈ (14) y

  ,

	acc i	t	x i	t	i	acc i	A y	m	t	x i	&	str	t	str	t	acc i	x i	t	(18)

Table 1 :

 1

		Pedestrian's parameters		Footbridge's parameters
	symbol	value	symbol	value

Table 1 :

 1 Parameters used in simulations

Table 2 :

 2 Results obtained on the phenomenon of synchronization

& & & & & & & & & & & & & & & &

& & & & =; the indices -and + referring to values before and after collisions, respectively; ext f (resp. int f ), the exterior forces vector (resp. interior forces vector) of dimension N 2 ,

& & & &

& & & & & & &


Ferdaouss Mattoussi

Vincent Roca
email: vincent.roca@inria.fr

Bessem Sayadi
email: bessem.sayadi@alcatel-lucent.com

Design of Small Rate, Close to Ideal, GLDPC-Staircase AL-FEC Codes for the Erasure Channel

This work introduces the Generalized Low Density

Parity Check (GLDPC)-Staircase codes for the erasure channel, that are constructed by extending LDPC-Staircase codes through Reed Solomon (RS) codes based on "quasi" Hankel matrices. This construction has several key benefits: in addition to the LDPC-Staircase repair symbols, it adds extra-repair symbols that can be produced on demand and in large quantities, which provides small rate capabilities. Additionally, with selecting the best internal parameters of GLDPC graph and under hybrid Iterative/Reed-Solomon/Maximum Likelihood decoding, the GLDPC-Staircase codes feature a very small decoding overhead and a low error floor. These excellent erasure capabilities, close to that of ideal, MDS codes, are obtained both with large and very small objects, whereas, as a matter of comparison, LDPC codes are known to be asymptotically good. Therefore, these properties make GLDPC-Staircase codes an excellent AL-FEC solution for many situations that require erasure protection such as media streaming.

I. INTRODUCTION

Low Density Parity Check (LDPC) codes have been intensively studied during the last decade due to their near-Shannon limit performance under iterative Belief-Propagation (BP) decoding [START_REF] Gallager | Low Density Parity Check Codes[END_REF][START_REF] Mackay | Near Shannon limit performance of low density parity check codes[END_REF][START_REF] Tanner | A recursive approach to low complexity codes[END_REF]. A (N, K) LDPC code, where N is the code length and K is its dimension, can be graphically represented as a bipartite graph with N "variable nodes" (VNs) and M = N -K "check nodes" (CNs). Equivalently, LDPC codes can be represented through their H L parity check matrix translating the connection between (VNs) and (CNs). The degree of a VN or a CN is defined as the number of edges connected to it. A VN of degree n can be interpreted as a "Length Repetition Code" (n, 1), i.e. as a linear block code repeating n times its single information bit towards the CN set. Similarly, a CN of degree n can be interpreted as a Single Parity Check (SPC) code (n, n -1), i.e. as a linear block code associated with one parity equation.

To improve error floor, minimal distance and decoding complexity performances, a generalization of these codes was suggested by Tanner in [START_REF] Tanner | A recursive approach to low complexity codes[END_REF], for which subsets of the variable nodes obey a more complex constraint than an SPC constraint. The SPC check nodes in a GLDPC structure are replaced with a generic linear block codes (n, k) referred to as sub-codes or component codes while the sparse graph representation is kept unchanged. More powerful decoders at the check nodes have been investigated by several researchers in recent years after the work of Boutros et al. [START_REF] Boutros | Generalized Low Density (Tanner) Codes[END_REF] and Lentmaier and Zigangirov [START_REF] Lentmaier | On Generalized Low-Density Parity-Check Codes Based on Hamming Component Codes[END_REF] where BCH codes and Hamming codes were proposed as component codes respectively. Later several works, on several channels, have been carried out in order to afford very large minimum distance and exhibit performance approaching Shannon's limit. Each construction differs to other by modifying the linear block codes (components codes) on the check nodes such as [START_REF] Yue | Generalized Low-Density Parity-Check Codes Based on Hadamard Constraints[END_REF][START_REF] Chen | A Hybrid Coding Scheme for the Gilbert Elliott Channel[END_REF][START_REF] Miladinovic | Generalized LDPC Codes with Reed-Solomon and BCH Codes as Component Codes for Binary Channels[END_REF][START_REF] Djordjevic | Generalized Low-Density Parity-Check Codes for Optical Communication Systems[END_REF][START_REF] Wang | Doubly Generalized LDPC Codes[END_REF][START_REF] Paolini | Analysis of Doubly-Generalized LDPC Codes with Random Component Codes for the Binary Erasure Channel[END_REF], or/and the distribution of the structure of GLDPC codes [START_REF] Yue | Generalized Low-Density Parity-Check Codes Based on Hadamard Constraints[END_REF] to offer a good compromise between waterfall performance and error floor under iterative decoding.

Recently, a construction of GLDPC using LDPC Staircase code as base code and Reed-Solomon (RS) codes as component codes has been proposed in [START_REF] Cunche | Lowrate coding using incremental redundancy for GLDPC codes[END_REF]. This construction allows each component code to produce a potentially large number of repair symbols in terms of RS codes (named extrarepair symbols) on demand, a feature that is well suited to situations where the channel conditions can be worse than expected, or to fountain like content distribution applications. The production of these extra-repair symbols allow to extend the initial LDPC Staircase code (based code) to a GLDPC code and very small rates are easily achievable. In this work the performance of the GLDPC-Staircase codes are assessed through the use of an Iterative IT/RS decoding. However it is shown that the extra-repair symbols don't totally solve the problem of small stopping sets that stuck the IT decoding process (i.e. IT remains sub-optimal).

In this paper, we revisit the [START_REF] Cunche | Lowrate coding using incremental redundancy for GLDPC codes[END_REF] proposal by improving the way of generating the extra-repair symbols and also by analyzing its performance in a finite length behavior. First, the extra-repair symbols are now generated through the use of "quasi" Hankel Matrix. The use of "quasi" Hankel matrices is well adapted to situations where the RS code parameters are dynamic determined. Thanks to the "quasi" Hankel matrix structure, the generator matrix creation complexity is significantly reduced without compromising its MDS characteristics. Second, the decoding algorithm is improved through the use of Maximum Likelihood (ML) decoder. In order to reduce the decoding complexity, we introduce an IT/RS/ML decoding scheme, called here hybrid decoding. Through it, the configuration parameters of the GLDPC codes such as base code rate, source variable nodes degree, and the extra-repair symbols distribution are discussed.

Simulation results show that our GLDPC construction features excellent decoding performance (i.e. good waterfall region and small error floor) and approach channel capacity even for small objects, both in terms of average decoding overhead and error floor. Additionally, the proposed GLDPC codes can easily be tuned to behave either like predefined rate LDPC-Staircase codes at one extreme, or like a single Reed Solomon code at another extremism, or like small rate codes. This flexibility makes it possible to tune GLDPC-Staircase codes to better match each use-case requirements, making GLDPC-Staircase codes an almost universal Application-Level FEC (AL-FEC) solution.

The paper is organized as follows. Following the introductory section, section II focuses on the design of GLDPC Staircase codes based on RS codes using "quasi" Hankel matrices. Section III provides a performances evaluation of these codes. Finally, we conclude.

II. PROPOSED LOW RATE CODING SCHEMA

In this section we first introduce the GLDPC-Staircase code design, as well as their encoding and decoding.

A. GLDPC-Staircase code construction

The GLDPC-Staircase (N G , K) codes studied in this paper can be represented by a Tanner graph (Figure 1) with the following meaning: For the reasons detailed in section III-B1 (i.e. improved performance under ML decoding), we assume that all the check nodes have the same number, E, of extra-repair symbols. So for a fixed LDPC-Staircase base code rate, r L , the code rate of the GLDPC-Staircase code is given by:

r G = r L 1 + (1 -r L )E (1)
Let N L and K be the LDPC-Staircase code length and dimension, and N G be the length of the GLDPC-Staircase code (which is also of dimension K). Then N G = N L +n extra , where n extra is the total number of extra-repair RS symbols.

Let H L be the binary parity-check matrix of the LDPC-Staircase code, of size M L = N L -K rows and N L columns.

H L has the form (H 1 |H 2 ). H 1 is the M L × K left-hand side part (information part
) and each column is of degree N 1 (number of "1s" per column). H 2 is the M L × M L righthand side part (redundancy part) and features a staircase (i.e. double diagonal) structure. H 1 is created in a fully regular way, in order to have constant column and row degrees. More precisely, each column of H 1 is of degree N 1 (number of "1s" per column), which is an input parameter during the LDPC-Staircase code creation [START_REF] Roca | Low density parity check (ldpc) staircase and triangle forward error correction (fec) schemes[END_REF]. Each row m of H 1 is of degree

N1 1 r L -1
, and because of the staircase structure of H 2 , depending on whether m = 1 or m > 1, a row m of H L is of degree

d 1 = N1 1 r L -1 + 1 and d m>1 = N1 1 r L -1 + 2.
The E extra-repair symbols associated to the m th row of

H L are generated by RS(n m , k m ) encoding over GF (2 8 ).
Here n m and k m are respectively the RS code length and dimension, and they are related to the other parameters of row m as follows. For row m > 1, the various source symbols (i.e. from the user point of view) that are involved in this row plus the previous repair symbol are considered as source symbols from the RS point of view. The new LDPC-Staircase repair symbol for this row plus the E extra-repair symbols are considered as repair symbols from the RS point of view. For the first row the only difference is the fact there is no previous repair symbol (it's the beginning of the staircase). So:

n m = k m + 1 + E, with k m = d m -1 (no matter the row).

B. Encoding algorithm of GLDPC codes

Two types of repair symbols are produced during encoding: For each row m, the E extra-repair symbols are computed by multiplying the k m symbols of LDPC symbols by the systematic generator matrix G m of an RS (n m , k m ) code based on an "quasi" Hankel matrix. As mentioned in the previous section, the k m symbols are defined by x symbols plus p m-1 (if m = 1). The construction of this systematic matrix is simply obtained by concatenating the identity matrix I km with the k m × (n mk m ) matrix A km,nm-km , where A km,nm-km is a rectangular sub-matrix of the maximal triangular array "quasi" Hankel matrix. For more details on RS codes based on "quasi" Hankel matrices, please refer to [START_REF] Mattoussi | Complexity comparison of the use of Vandermonde versus Hankel matrices to build systematic MDS Reed Solomon codes[END_REF] [START_REF] Roth | On Generator Matrices of MDS Codes[END_REF]. For GLDPC-Staircase codes, this choice has the advantage that p m can be considered indifferently as an LDPC-Staircase or RS repair symbols 1 .

M L LDPC-Staircase repair symbols, (p 1 , • • • , p M L ),
An advantage of this encoding approach is the fact that extra-repair symbols can be produced incrementally, on demand, rather than all at once (unlike LDPC-Staircase repair symbols). Their number can also be rather high since it is only limited by the finite field size, usually GF(2 8 ). Said differently, GLDPC-Staircase codes can easily and dynamically be turned into small rate codes.

C. Decoding algorithm of GLDPC codes

The decoding of GLDPC codes, named hybrid decoding, consists of a joint use of a four decoders XORing two symbols. However it allows reaching the maximum correction capabilities of the code. It is equivalent to the ML decoding over the full system mentioned above, but the system on which it is applied is hopefully simplified by the previous three decoders. Decoding succeeds if one or several of these decoders succeed and recover all the missing source symbols.

This description is detailed in Algorithm 1 where N b ∈ N source or/and LDPC-Staircase or/and extra-repair symbols are received.

III. PERFORMANCE EVALUATION RESULTS

A. Experimental setup

We have developed a GLDPC-Staircase codec based on RS(2 8 ) codes ("quasi" Hankel matrices) under hybrid decoding, in C language, using the OpenFEC.org project (htpp://openfec.org).

Algorithm 1 Hybrid decoding of (N G , K) GLDPC-Staircase codes

1: for s = 1 to N b do 2:
Step1: Decoding with new symbol(s)

3: if (s is LDPC symbol) then 4:
Select IT decoding with received symbol()

5:
if IT decoding successf ul then 6:

Recover symbol 7:

Go to Step1(s ← recovered symbol)

8: else 9:
Select RS decoding with received symbol()

10:
if RS decoding successf ul then 11:

Recover symbols

12:

Go to Step1 (s ← each recovered symbol) Experiments are carried out considering a memoryless channel along with a transmission scheme where all the source and repair symbols are sent in a fully random order. This choice has the benefit to make the performance results independent of the loss model. Instead the target channel loss rate is the only parameter that needs to be considered.

We determined the performances of these codes by testing several LDPC-Staircase matrix that are generated randomly. For each test, a different LDPC-Staircase matrix is used (more precisely we change the PRNG seed used to create the matrix) and then the results, averaged over tests, show the average behavior of GLDPC-Staircase codes.

Performance evaluations make use of the following two metrics:

• Decoding inefficiency ratio and decoding overhead: the decoding inefficiency ratio is defined as the ratio between the number of symbols needed for a successful decoding over the number of source symbols : ineff = nb symbols needed K = 1 + ε where ε is the transmission overhead (ε is often expressed as a percentage). An MDS code, when there is a single block, has a decoding inefficiency ratio equal to 1 (i.e. any K subset of symbols are sufficient for decoding to succeed). With non-MDS codes (e.g. LDPC-Staircase codes) we are usually considering average values, and on average decoding is successful after receiving K(1 + ε) symbols (values for a target decoding success, say 99%, can also be considered instead of the average 50% target we use). • Decoding failure probability: decoding is said to fail when one or more erased source symbols are not recovered after having received all the non-erased encoding symbols and having launched the decoder. This metric defines at a receiver as a function of the number of received symbols, or equivalently the percentage of symbols lost over the network. Since the hybrid decoding scheme we introduced is a Maximum Likelihood (ML) decoder (i.e. provide the same erasure recovery performance), we only mention "ML decoding" in the remaining of this paper. The hybrid approach will only impact decoding complexity and speed, which is not addressed in this work.

B. Preliminary tests and internal parameters set up

GLDPC-Staircase codes depend on three important internal parameters namely the extra-repair symbols distribution across the H L rows, the base code rate r L , and the N1 parameter of the base code. We begin by analyzing the impact of these three parameters on the erasure recovery performance of GLDPC codes, then using the best values of these parameters we evaluate the performance gains made possibly by the use of the decoding scheme, ML versus IT/RS.

1) Performance as a function of the extra-repair symbols distribution:

The distribution of extra-repair symbols per check node is an important parameter. In [START_REF] Cunche | Lowrate coding using incremental redundancy for GLDPC codes[END_REF] it is shown that GLDPC codes with IT/RS decoding perform the best under a uniform (irregular) distribution rule. However in our work we consider an ML decoding scheme and the situation is completely different. Therefore we tested the uniform (irregular) distribution of [START_REF] Cunche | Lowrate coding using incremental redundancy for GLDPC codes[END_REF] against a (regular) constant distribution, where the same number of extra-repair symbols are generated for each parity check matrix row. Figure 2 shows the average decoding inefficiency ratio (i.e, average of 1000 GLDPC codes with r G = 1

2 ) for the proposed GLDPC codes, for different object sizes. It shows that the constant distribution performs significantly better under ML decoding, both with small and large objects. Therefore in the remaining of this work we only focus on this constant distribution, where there are exactly e m = E extra-repair symbols per row.

2) Performance as a function of the base code rate r L : Let us consider a GLDPC-Staircase rate r G . Several values of the base code rate r L , or equivalently of E, enable to achieve r G (see Eq. 1). However choosing a value impacts the performance achieved. In Figure 3, we plot the average decoding inefficiency ratio (i.e, average of 1000 GLDPC codes with r G = 1

3 ) for different object size. This figure provides

(r G = { 1 3 , 1 2 , 2 3 }) AS A FUNCTION OF r L
that increasing rate r L (i.e, increasing the number E), the average decoding inefficiency ratio quickly approaches 1 (i.e. no overhead, decoding is possible with exactly K symbols) as

E = 3 (i.e, r L = 2
3 ), even for very small code dimensions. Based on EXtrinsic Information Transfer (EXIT) functions, asymptotic analysis technique, the impact of this parameter on GLDPC codes performances using the upper bound on the ML decoding threshold δML is studied. For more details on the asymptotic analysis of GLDPC-Staircase codes please refer to [START_REF] Mattoussi | Optimization with EXIT functions of GLDPC-Staircase codes for the BEC[END_REF]. Table I provides, for each r G , as seen in Figure 3 that increasing the base code rate r L (i.e, increasing E) quickly increases the upper bound on the ML threshold until it reaches a stable value (i.e. very close or equal to the Shannon limit δ sh ). For r L = 2 3 (i.e, E = 1 for r G = 1 2 , and E = 3 for r G = 1

3 ) the δML very close to δ sh , and above this base code rate 2 3 the upper bound on the ML threshold improves slightly. During the following tests we chose to set the base code rate equal to r L = 2/3.

3) Performance versus N1: Let us now adjust the third internal parameter of the GLDPC-Staircase codes, N1, which represents the source variable nodes degree of the base code matrix H L .

Table II provides the average decoding inefficiency ratio (i.e, average of 1000 GLDPC codes for each r G ) for different values of N1 under IT/RS and ML decoding of GLDPC codes with K = 1000. Table II shows that, under IT/RS decoding, increasing N1 induces an increase in the average decoding inefficiency ratio because the increase of N1 causes the "densification" of H 1 sub-matrix of base code (i.e, the source sub-matrix) and maybe the decrease of the smallest stopping sets size 2 . This impact leading to the uselessness of extra-repair symbols which are used to cope the problem of the small stopping sets in the base code graph. Even with the increase of E (i.e, reduce of the GLDPC code rate), the problem is not resolved. Therefore, the increase of N1 leads to the deterioration of the ability of correction IT/RS decoding. Table II shows a different behavior for pure ML decoding 3 as the one observed for IT/RS decoding. The average decoding inefficiency ratio decreases by increasing N1 for all the GLDPC code rates. With ML decoding, the most significant performance gains are obtained by switching from N 1 = 3 to 5 and since this value the performance improves slightly (no significant impact of N1). For a given value of N1, the average decoding inefficiency ratio decreases by increasing E per check nodes (i.e, adding extra-repair symbols nodes to the base code graph) due to the role of extra-repair symbols. This value 5 also limits the performance degradation of IT/RS decoding compared to values N 1 > 5. Therefore in the remaining of this work we only focus on the case where N 1 = 5. 4) Performance of IT/RS versus ML decoding: Now that we chose N 1 = 5 as a value that provides a good balance, we can quantify the erasure recovery performance gains made possible by the use of ML decoding 3 as the second step of 2 The residual set of erasures after belief propagation decoding is exactly equal to the maximum size stopping set that is a subset of the originally erased code symbols. The performance of LDPC codes is limited by the size of smallest stopping sets, i.e design of a good LDPC code graph assumes maximizing the size of their minimal stopping sets 3 Here we refer to the ML as Maximum Likelihood decoding but not the hybrid decoding our hybrid decoding scheme. Figure 4 shows the average decoding inefficiency ratio (i.e, average of 1000 GLDPC codes with r G = 1

2 ) versus various object sizes. This figure confirms the major performance gains made possible by the use of ML decoding 3 with GLDPC-Staircase codes for all object sizes. For instance for K = 1000 we obtain a gain equal to 22%. We can remark also based on Table II that for given N1 value, adding extra-repair symbols to the base graph causes the degradation of IT/RS decoding performances unlike the ML decoding. Therefore in the remaining of this work we only focus on ML decoding.

C. Average performance results

Now that the internal parameters are set, let us analyze the GLDPC-Staircase performances in various situations. We show in this section that the average performance achieved is exceptionally good, almost that of ideal codes, no matter the object size. In figure 5(a) we plot the average decoding inefficiency ratio (i.e, average of 1000 GLDPC codes), as a function of the object size, with various curves corresponding to the various code rates (i.e, r G = 1 3 (E = 3), r G = 1 2 (E = 1), and r G = 2 3 (E = 0), and we do the opposite in Figure 5(b). We see in both figures that the GLDPC-Staircase codes with E = 1 or E = 3 exhibit exceptional erasure recovery capabilities, even for tiny objects. On the opposite, codes with E = 0 corresponding to the base codes have performance level significantly smaller for small object sizes. Hence, the addition of extra-repair symbols makes the correction capabilities of GLDPC-Staircase codes under ML decoding close to that of ideal, MDS codes, both for tiny or large objects. These tests also show that GLDPC-Staircase codes can easily achieve very small code rates while keeping exceptionally high erasure recovery performance as shown also in Table II.

D. Decoding failure probability results

We continue the analysis with decoding failure probability, which enables us to more carefully analyze the GLDPC-Staircase codes behavior as a function of the number of received symbols and loss percentage.

Figure 6 shows the average results of 10 7 GLDPC codes with r G = 1 2 for K = 1000, K = 256, and K = 32. The black vertical line corresponds to ideal, MDS code, for which the decoding failure is equal to 0 as long as the experienced loss rate is strictly inferior to 50%. This figure confirms that GLDPC-Staircase codes feature a very small decoding overhead, close to that of ideal, MDS codes, with no visible error floor above 10 -5 decoding failure. This is obvious with objects of size 1000 and K = 256, it remains almost true with the K = 32 case (i.e, error floor started from 8.10 -6 with 42% of loss for K=32, bellow 5.33.10 -6 with 49.45% of loss for K=1000).

Table III gives additional details where we performed 6.10 6 , and 2.10 6 GLDPC codes (r G = 1

2 ) for K = 1000 and K = 32 respectively. It provides the average decoding failure probability as a function of the number of received symbols above K, i.e. as a function of the overhead. It shows that for large or small object sizes, a few symbols above K is sufficient for decoding to success with a high probability. For instance, with K = 1000 received symbols, 1819473 among 6.10 6 GLDPC codes can recover all the erased symbols (i.e, 69.67 % as decoding failure probability) and with two received symbols more than K, 5703206 among 6.10 6 GLDPC codes can decode all the erased symbols (i.e, 4.9 % as decoding failure probability). With K = 32 received symbols, GLDPC codes can decode all the erased symbols with decoding failure probability equal to 3.0584% (i.e, 1938832 among 2.10 6 GLDPC codes that can recover all the erased symbols) and with two received symbols more than K, 1999778 among 2.10 6 GLDPC codes can recover all the erased symbols (i.e, 0.011% as decoding failure probability).

IV. CONCLUSIONS

This paper has introduced the GLDPC-Staircase codes based on RS codes ("quasi" Hankel matrices) under IT/RS/ML decoding. To obtain excellent performances of these codes, we adjust firstly their internal parameters. We have shown that these codes have two key benefits: on the one hand, a significant number of extra-repair symbols can be produced on the fly, which provides "small rate" rate capabilities; on the other hand, these codes feature both a very small decoding overhead and low error floor, not only for large object, but also for very small ones. These properties make GLDPC-Staircase codes an excellent AL-FEC solution for many situations that require erasure protection. If the current paper focuses on practical performance evaluations, in [START_REF] Mattoussi | Optimization with EXIT functions of GLDPC-Staircase codes for the BEC[END_REF] we introduce an asymptotic analysis, in terms of EXtrinsic Information Transfer functions and we derive an upper bound of the ML decoding threshold based on the area theorem. Both papers nicely complement with one another. Future works will focus on encoding and decoding speed and complexity evaluations, a key aspect for practical realizations.

•Fig. 1 .

 1 Fig. 1. GLDPC-Staircase (13, 4) code with two extra repair symbols per check node (E = 2).

  and n e extra-repair symbols, (e 1 , • • • , e ne ). Let S = (s 1 , s 2 , • • • , s K ) be the source symbols. The (p 1 , • • • , p M L ) repair symbols are computed us usual, by "following the stairs" of the H L matrix. Moreover, for each row m in H L , p m is the XOR sum of LDPC source symbols x = (x 1 , • • • , x Km ), where x is a subset of S that correspond to a 1 coefficient in row m of H L , (and the LDPC repair symbol p m-1 if m = 1).

Fig. 3 .

 3 Fig. 3. Performance as a function of base code rate r L with N1= 5 and hybrid decoding.

Fig. 4 .

 4 Fig. 4. Performance for (IT/RS) versus ML decoding schemes, with r G = 1 2 .

Fig. 5 .

 5 Fig. 5. Average decoding inefficiency ratio as a function of the object sizes (a) and code rates (b).

Fig. 6 .

 6 Fig.6. Average decoding failure probability as a function of the channel loss percentage experienced for GLDPC codes (r G = 12 ) with K = 32, K = 256, and K = 1000.

  END Decoding with new symbol(s) 23: end for 24: if (All or some of source symbols are still erased ) then

	13:	end if
	14:	end if
	15:	else if (s is extra repair symbol) then
	16:	Select RS decoding with received symbol()
	17:	if RS decoding successf ul then
	18:	Recover symbols
	19: 20:	Go to Step1 (s ← each recovered symbol) end if
	21:	end if
	22:
	25:	Step2: ML decoding
	26:	Select Binary ML decoding
	27:	if (All or some of source symbols are still erased ) then
	28:	Select Non Binary ML decoding
	29:	if (All or some of source symbols are still erased ) then
	30:	return Decoding is finished with failure
	31:	else
	32:	return Decoding is completed successf ully
	33:	end if
	34:	else
	35:	return Decoding is completed successf ully
	36:	end if
	37: else
	38:	return Decoding is completed successf ully
	39: end if

TABLE II PERFORMANCE

 II AS A FUNCTION OF N1 AND THE DECODING SCHEME (IT/RS VERSUS ML), FOR K = 1000 AND r L = 2 3 .

  TABLE III AVERAGE DECODING FAILURE PROBABILITY OF GLDPC-STAIRCASE CODES (r G = 1 2 ) UNDER ML DECODING FOR K=32 AND K=1000

Technically speaking, this is made possible by the fact that the coefficients on the first column of G i are all equal to 1.

V. ACKNOWLEDGMENTS

This work was supported by the ANR-09-VERS-019-02 grant (Adaptable, Robust, Streaming SOlutions (ARSSO) project) and by the INRIA -Alcatel Lucent Bell Labs joint laboratory.


P Pécol

P Argoul

S Dal Pont

S Erlicher

Vibrations, Shocks and Noise A new crowd movement modeling for pedestrians who hold hands

Keywords: discrete model, crowd movement modeling, oriented pedestrian, subgroup behaviour, pedestrians holding hands, emergency evacuation ;

The aim of this communication is to propose a new crowd movement modeling considering a new particular subgroup behavior: pedestrians holding hands. The first step in this study is to improve our 2D crowd movement modeling [1-3] introducing the orientation of each pedestrian, rarely used in literature. The next step is to describe the particular subgroup behavior of the pedestrians who hold hands. An at-a-distance velocity of deformation allows one to control the interaction between linked pedestrians. The proposed model takes into account the effects of the subgroup as a continuous deformation within the particular subgroup of pedestrians who hold hands while existing methods in literature describe the cohesion of the "normal" subgroup with external forces. Finally, numerical simulations will be presented to study the influence of the pedestrian's orientation on his movement and to compare subgroup behavior models.

Introduction

Crowd movement modeling is currently a challenging problem in various fields, such as in architecture and civil engineering to study flows of people in confined spaces, in the entertainment and visual effects industry to create realistic models of crowds for use in movies, and in psychology/sociology to validate human behavior models under different environmental conditions like panic situations. Over the last fifty years, many studies have been performed to describe the behavior of walking pedestrians. In literature, several discrete models of crowd movements [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF][START_REF] Zheng | Modelling crowd evacuation of a building based on seven methodological approaches[END_REF][START_REF] Venel | Modélisation mathématique des mouvements de foule[END_REF][START_REF] Helbing | Simulating dynamic features of escape panic[END_REF][START_REF] Antonini | Simulation of pedestrian behaviour using a discrete choice model calibrated on actual motion data[END_REF][START_REF] Langston | Crowd dynamics discrete element multi-circle model[END_REF] exist and allow to well reproduce some observed self organization phenomena (counter flow lines, formation of arches, etc.).

Authors have proposed in [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF] a new 2D microscopic approach in which the movement of each pedestrian is represented in time and in space. To control the crowd movement on the ground, three targets have been managed: (i) multiple simultaneous collisions, i.e. to detect and to treat each local interaction, such as pedestrian-pedestrian and pedestrian-obstacle; (ii) the desire of each pedestrian to move in a particular direction with a specific speed at each time; (iii) the possibility to add forces to make the behavior of the pedestrians more realistic (social forces, subgroup forces, etc.).

Pedestrian-pedestrian and pedestrian-obstacle collisions have been modeled using an approach devoted to the movements of particles. The non-smooth approach proposed by Frémond [START_REF] Frémond | Rigid bodies collisions[END_REF][START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF][11] that has mainly been applied to the movement of granular media was selected to treat the contact and to control the type of collisions among particles. This original approach is based on the theory of collisions of rigid bodies, first proposed by Frémond [START_REF] Frémond | Rigid bodies collisions[END_REF][START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF] in a rigorous thermodynamic frame, along the lines of the works of [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF]. Frémond introduces the concept of the coefficient of dissipation to handle the rebound and the velocities of particles after each collision. A further constraint condition is applied to velocities to avoid overlap. The description of the behavior of a collection of discrete bodies is based on the consideration that the global system is deformable even if particles are rigid.

To represent a pedestrian, a willingness is given to each particle to move according to a given possibly time varying target. Social forces as well as a desired direction/velocity are introduced in order to simulate the behavior of pedestrians [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF].

The next step is how to include the behavior of groups or sub-groups in this modeling. A "group" is defined as a physical collection of people following the same route, but who may or may not be part of the same social group and a "subgroup" is defined as people within the same physical "group" who want to stay together [START_REF] Singh | Modelling subgroup behavior in crowd dynamics DEM simulation[END_REF] like those of friends or a family. Several studies have revealed that smaller subgroups constitute the majority of the people in a crowd [START_REF] Singh | Modelling subgroup behavior in crowd dynamics DEM simulation[END_REF][START_REF] Aveni | The not-so-lonely crowd: friendship groups in collective behaviour[END_REF]. But very few studies are able to model the "subgroup" behavior [START_REF] Singh | Modelling subgroup behavior in crowd dynamics DEM simulation[END_REF][START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF]. A particular subgroup concerns pedestrians who hold hands. It is recognized that holding hands is the most effective way of keeping children safe from traffic injury.

The aim of this communication is then to propose a new modeling of the behavior of pedestrians who hold hands. To describe the particular subgroup behavior of the pedestrians who hold hands, we were inspired by the article [START_REF] Caselli | Collision of three balls on a plane[END_REF] describing the interactions between rigid particles as an at-a-distance velocity of deformation. The proposed model takes into account the effects of the subgroup as a continuous deformation within the subgroup of pedestrians who hold hands while existing methods in literature describe the cohesion of the subgroup with external forces [START_REF] Singh | Modelling subgroup behavior in crowd dynamics DEM simulation[END_REF][START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF]. This new subgroup model needs that pedestrians be oriented, so the rotation of the pedestrians on themselves will be integrated in the model to obtain an improvement of the 2D microscopic approach [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF]. It is important to note that the pedestrian's orientation is rarely used in literature [START_REF] Antonini | Simulation of pedestrian behaviour using a discrete choice model calibrated on actual motion data[END_REF][START_REF] Langston | Crowd dynamics discrete element multi-circle model[END_REF].

Numerical simulations will be presented to study the influence of the pedestrian's orientation on his movement and to compare subgroup behavior models.

Crowd model

The proposed crowd model is an improvement of the one presented in [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF], taking into account the rotation of each pedestrian on himself. First, the granular model allowing to manage multiple simultaneous collisions, is presented. Then, the adaptation of this granular model to the crowd one is done. [START_REF] Frémond | Rigid bodies collisions[END_REF][START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF] Let us consider a system of N circular particles moving in a plane, defined by their mass i m , their inertia i I , their radius i r , their center position (i.e. center of gravity i G in this study) described by the vector

A non-smooth granular model based on the theory of Frémond

( ) 2 , R q q q y i x i i t ∈ =
, their velocity of center position denoted by

dt t q d t u i i ) ( ) ( = , their orientation i θ and their rotational velocity dt d t i i θ ω = ) (
. Contacts between particles are assumed to be punctual. Noting

( ) i i t i t q r θ , = , ( ) i i t i t i t u dt t r d t v ω , ) ( ) ( = =
and the generalized displacement vector r of size N 3 , ( )

N t t t t r r r r , ... , , 2 1 =
, the relative deformation velocity between the i th and the j th particles colliding at point

j i A , is: ( ) ( ) ( ) ij j j j ij i i i ij A G u A G u t v × + - × + = ∆ ω ω . ( 1
)
The dynamics equation for each particle can be written as the following differential system:

( ) ( ) ( ) ( ) ( ) ( )      + - = - + - = - + 3 ) ( ) ( 2 ) ( ) ( int int everywhere t p t p t v t v M everywhere almost t f t f t v M ext ext & where M is the N N 3 3 ×
inertial matrix of all the particles; ext f (resp. int f ) is the exterior forces vector (resp. interior forces vector) of dimension N 3 applied to the deformable system. The existence of a solution of the system given by Eqs. ( 2) and ( 3) is proven in [START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF]11,[START_REF] Cholet | Collision d'un point et d'un plan[END_REF]. Equation (2) describes the smooth evolution of the multiparticle system, whereas (3) describes its non-smooth evolution during a collision. Hence, Eq. ( 2) applies almost everywhere, except at the instant of the collision, where it is replaced by Eq. ( 3). When contact is detected, velocities of colliding particles are discontinuous, and so in Eq. ( 3), the percussions int p and ext p , interior and exterior to the system respectively, are introduced. By definition, percussions have the dimension of a force multiplied by a time. The int p percussions are unknown; they take into account the dissipative interactions between the colliding particles (dissipative percussions) and the reaction forces that permit the avoidance of overlapping among particles (reactive percussions). Frémond [START_REF] Frémond | Rigid bodies collisions[END_REF][START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF] defined the velocity of deformation at the moment of impact ( ) ( ) ( )

2 - + ∆ + ∆ v v
and showed that int p is defined in duality with ( ) ( ) ( )

2 - + ∆ + ∆ v v
according to the work of internal forces. He then introduced a pseudopotential of dissipation Φ that allows us to express int p as:

( ) ( )

        ∆ + ∆ Φ ∂ ∈ - + 2 int v v p ( 4
)
where the operator ∂ is the subdifferential that generalizes the derivative for convex functions [START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF]. The convex function Φ [START_REF] Moreau | Sur les lois du frottement, de la viscosité et de la plasticité[END_REF] is defined as the sum of two pseudopotentials:

r d Φ + Φ = Φ
; where d Φ and r Φ characterize the dissipative and reactive interior percussions respectively. The pseudopotential d Φ is chosen to be quadratic:

( )

2 2 y K y d = Φ
; where K is a coefficient of dissipation. This choice allows one to find the classical results when the coefficient of restitution is used. Other choices of d Φ allow one to obtain a large variety of behaviors after impact [START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF][START_REF] Cholet | Collision d'un point et d'un plan[END_REF].

In Eq. ( 3), the problem is to find the velocity + v after particles' collision. To determine + v , the following constrained minimization problem has to be solved:

( ) ( ) ( ) [ ] Y M p M v Y Y M Y R Y X ext t t N 1 3 2 min arg - -+ - ∆ Φ + ∈ = ( 5
)
where the solution 2

- + + = v v X .
Proof of the existence and uniqueness of the velocity + v after the simultaneous collisions of several rigid solids, as well as the dissipativity of the collisions, is presented in [START_REF] Frémond | Rigid bodies collisions[END_REF][START_REF] Frémond | Theoretical approach to instantaneous collisions and numerical simulation of granular media using the A-CD 2 method[END_REF][11].

Adaptation of granular approach to the crowd

A pedestrian can be represented as a circular particle by giving it willingness, i.e. a desire to move in a particular direction with a specific speed at each time. The adaptation [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF] is improved to take into account the rotational degree of freedom.

The first step of the modified approach is to give a desired trajectory to each particle. Several definitions of the desired trajectory of any one pedestrian are possible: either (i) the most comfortable trajectory for him, where he must exert the least effort, e.g. by avoiding the stairs or making the fewest changes in direction, etc.; (ii) the shortest path; or (iii) the fastest path to move from one place to another [START_REF] Hoogendoorn | Microscopic pedestrian wayfinding and dynamics modelling[END_REF]. It is possible to combine two strategies in the same simulation or to change the preferred strategy for any reason during the simulation.

The strategy of the shortest path from one point to another is implemented through a Fast Marching algorithm [START_REF] Kimmel | Fast marching methods for computing distance maps and shortest paths[END_REF] and is used to obtain the desired direction i d e , of an individual i [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF]. This direction depends on the environment in which pedestrians walk (obstacles, etc.), the time of day, and also the characteristics of the individual (gender, age, hurried steps or not, etc.). It is defined by: ( )

i d i d i d u u t e , , ,

=

; where i d u , is the desired velocity of the i-th pedestrian. The amplitude i d u , of the desired velocity represents the speed at which the i-th pedestrian wants to move, and it can be influenced by his nervousness. This velocity is chosen by following a normal distribution with an average of In the second step, the desired velocity of each pedestrian is introduced into the original discrete models to simulate crowd movement. Let ( ) ( )

t h t f a = int , where
( ) t h a
allows one to give a desired direction, amplitude of the velocity and orientation to each pedestrian. Each component

( ) t h a i of the vector of dimension N 3 : ( ) a N t a t a t a t h h h h , ... , , 2 1 =
, is associated with pedestrian i and can be expressed as

( ) (
)

) ( ), ( t l t f t h a i a i t a i t =
, where ( )

t f a i
is the so-called acceleration force [START_REF] Helbing | Simulating dynamic features of escape panic[END_REF] giving the desired direction and amplitude of the velocity to the pedestrian, and ) (t l a i is an angular torque giving him his desired orientation.

( )

t f a i is defined in [1-3, 6] by: ( ) ( ) ( ) i i i d i d i a i t u t e u m t f τ - = , , (6)
where i u is the actual velocity and i τ is a relaxation time, which specifies how long the pedestrian will take to recover his desired velocity either after a contact or after he suddenly changes his walking direction. Helbing [START_REF] Helbing | Simulating dynamic features of escape panic[END_REF] chose s 5 . 0 = τ in his numerical simulations. Smaller values of i τ let the pedestrians walk more aggressively. The influence of the relaxation time parameter τ has been studied in [START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF]. The chosen value of this parameter is less than or equal to s 5 . 0 ; as such, several contacts may occur because the pedestrians walk aggressively. ) (t l a i is the combination of a linear spring (with the desired direction) and linear damper, it is defined by:

( ) i i d i a i c k t l ω θ θ - - - = , ) ( (7)
where i d , θ is the angle between i d e , and a reference direction

[ ] 0 ; 1 = x e
, k is the stiffness and c the damping.

The pedestrians' behavior can be enriched by adding other external social forces [START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF][START_REF] Helbing | Simulating dynamic features of escape panic[END_REF][START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF] so as to become more realistic (socio-psychological force, attractive force, group force, etc.).

Influence of the pedestrian's orientation

Two different simulations are presented in this part to study the influence of the pedestrian's orientation. The first one concerns the interaction between two pedestrians and the second one is about an evacuation exercise.

Three cases are studied in the following simulations, function of the coefficients of dissipation for the tangential components of percussions T K and of the rotational degree of freedom: (i

) T K is neglected, so kg K T 0 = ; (ii) kg K T 100
= and the rotation of each pedestrian on himself is neglected, so the pedestrian is not oriented; and (iii) kg K T 100 = and the rotation of each pedestrian on himself is taken into account, so the pedestrian is oriented.

Pedestrian / pedestrian interaction

Two identical pedestrians are moving in opposite direction. After the collision, they try to recover their desired velocities. The numerical simulations are done depending on the value of T K and on the rotation of each pedestrian on himself. Figure 1 shows the behavior of the pedestrian after a shock is different when an orientation is given to the pedestrian (red curve) or not (other curves). In literature, crowd models rarely take into account this criterion. It seems to be interesting to consider it now.

Evacuation of a square room

The aim of this part is to compare an evacuation situation for the three different cases, considering only the way of treating the local pedestrian-pedestrian contact and the local pedestrian-obstacle contact.

We consider a square room of side 5 m, where 20 pedestrians want to escape by a door of 1 m width (Fig. 2). The parameters used in simulations are given in Table 1.

Table 1. Evacuation of a classroom -parameters used in simulations (* uniformly distributed within their range). The response time is the time needed by pedestrian i to start evacuating after the triggering of the evacuation movement. For each studied case, 50 simulations are performed (Fig. 3). The initial conditions of these runs are the same for each case. Figure 3 a, b andc shows the linear regression of the 50 simulations in each case. The slope allows us to obtain the average flow through the door. The obtained linear regressions can be compared in Figure 3 d. The orientation of each pedestrian which is usually neglected in literature, when kg K T 100 = and when the rotation of each pedestrian on himself is taken into account, influences the results of this evacuation exercise. So, this rotational degree of freedom has to be considered in crowd movement in the future.

Fig. 1 .

 1 Fig. 1. Different trajectories of one pedestrian after a shock with an other one.

Fig. 2 .

 2 Fig. 2. Evacuation of a square room.

Fig. 3 .Fig. 5 .

 35 Fig. 3. Evacuation of a square room -50 simulations and linear regression -a) kg K T 0 = ; b) kg K T 100 = and the rotation of each on himself is neglected; c) kg K T 100 = and the rotation of each pedestrian on himself is taken into account; and d) comparison of the three linear regression.

Pedestrians holding hands

We propose a particular case of subgroup: pedestrians holding hands. While existing methods describe the cohesion of the subgroup with forces [START_REF] Singh | Modelling subgroup behavior in crowd dynamics DEM simulation[END_REF][START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF], the proposed method to describe this original behavior, inspired by [START_REF] Caselli | Collision of three balls on a plane[END_REF] proposing an at a distance interaction between rigid particles, is to model the effects of the subgroup as a continuous deformation of the system for both linked pedestrians.

The at a distance velocity of deformation is the derivative with respect to time of the squared distance between linked shoulders, i.e. points i A and j A for pedestrians i and j (Fig. 4). Taking this example, the at a distance velocity of deformation is expressed as:

In the pseudopotential of dissipation d Φ , as well as quadratic terms taking into account normal N K and tangential T K coefficients of dissipation, an other quadratic term V K is added for the coefficient of viscous dissipation.

On a same example, a comparison has been done between different subgroup approaches: Singh [START_REF] Singh | Modelling subgroup behavior in crowd dynamics DEM simulation[END_REF], Moussaïd [START_REF] Moussaid | The Walking Behaviour of Pedestrian Social Groups and Its Impact on Crowd Dynamics[END_REF], and pedestrians holding hands with

. A subgroups composed of three pedestrians collides with a fourth pedestrians. Figure 5 shows the proposed subgroup approach allows one to re-form the subgroup after a collision faster than the others. Moreover, only the proposed approach allows one to keep the shape of the subgroup throughout the simulation.

Conclusions

An improvement of the crowd movement model presented in [START_REF] Pécol | Modelling crowd-structure interaction[END_REF][START_REF] Pécol | Smooth / non-smooth contact modelling of human crowds movement: numerical aspects and application to emergency evacuations[END_REF][START_REF] Pécol | The non smooth view for contact dynamics by Michel Frémond extended to the modeling of crowd movements[END_REF] is proposed in this communication. First, the orientation of the pedestrians, usually neglected in literature's models, is introduced. Numerical simulations in an evacuation situation show that this criterion influences the pedestrian's behavior and has to be considered in the future.

Second, a particular case of subgroup is proposed: pedestrians holding hands. It's an original approach using the formalism of pseudopotential of dissipation to control the link between two pedestrians. A comparison of this new proposed approach with two others ones of the literature gives satisfactory results.