Capucine-Marin Dubroca-Voisin

TOWARDS A GLOBAL OBSERVATORY OF PEDESTRIAN FLOW MANAGEMENT LEVERS FOR THE RAILWAY SYSTEM

Crowd disasters are more and more known and documented, but the ways of preventing them don't follow closely that trend. Train and station operators face an increasing pressure with pedestrian flow management, as ridership often grows more quickly than resources, and cannot rely on a pedestrian flow overview. To help tackling these issues, we launch a global observatory of the current practices and research about pedestrian flow management in stations. This observatory collects, analyses, and synthetises the available information from train and station operators, existing databases, academic and grey literature.

In the literature, some papers and books are already explaining very interesting pedestrian flow management techniques, but no generalized analysis has been made. Several websites are listing crowd disasters and near-miss cases, but they generally don't analyse the causes of the accidents, and they don't offer a panorama of the levers that could be used for preventing disasters or in daily use. Moreover, flow management is not limited to safety and should also enable comfort, accessibility, or railway operation efficiency, which is not broadly covered by these sources. Finally, there is no consensual definition of the capacity of a station, making it difficult to assess the necessity of expansion works.

To address these issues, the observatory collects and shares well-sourced information about crowd and pedestrian flow management, focusing on the case of train stations. That observatory in run is an open-source and collaborative manner, involving railway companies, station managers, and researchers. It aims at providing a complete list of crowd accidents, and a near-complete list of flow managements levers that can be used in the context of train stations. This paper describes these three tasks of the Observatory: know the risks, collect the solutions, assess the saturation. It also explains how to contribute to this collective effort.

RISKS OF TRANSPORT-RELATED CROWD CRUSHES ARE BADLY KNOWN

One of the main reasons to improve pedestrian movement in train stations is safety. In this domain, the bigger risk seems to be the occurrence of a crowd crush (also incorrectly called stampede). A crowd crush can lead to dozens, or even hundreds or thousands of deaths, and a similar number of injuries, during a single event. Crowd crushes are not the only crowd-related accidents that can happen in a train station (we can specifically think to the case where a person falls on the track due to the density on the platform), but they might be an important danger. We chose to focus on the crowd crushes in the early stages of our work, as there were already lists of such events.

1.1

Knowing the history of the crowd crushes In a first work that we documented in (Dubroca-Voisin, 2023a), we compiled the crowd crushes listed in six different sources [START_REF] Asgary | World Crowd Disasters Web App V1[END_REF][START_REF] Feliciani | Introduction to Crowd Management Managing Crowds in the Digital Era : Theory and Practice[END_REF], Working with Crowds 2022;[START_REF] Still | Crowd Disasters | Prof. Dr. G. Keith Still[END_REF]Wikidata contributors, 2023;Wikipedia contributors, 2023). We showed that none of these sources were complete; we then may have a partial view of crowd crushes and miss some important events than happened in transport context. Even partial, our compilation might be the most complete list of crowd crushes in history. We were then able to distinguish crowd crushes that were transport-related.

1.2

A small number of events with no common causes In a second paper (Dubroca-Voisin, 2023c), we listed the 16 crowd crushes than were transport-related. Most of these events happened in train stations; some of them were linked to war events or public gatherings, but some others happened during normal train service operations. The five more recent events all happened in Indian stations (that could be linked to the growing pedestrian flows in the country but also to its better anglophone press coverage).

For each of these events, we analysed the potential causes. This work was based on press coverage, which can be very incomplete, therefore this analysis cannot be a scientific approach to these events. However, we can observe that it is very difficult to find common causes to these events: overcrowding generally played a role but the succession of events leading to a transport-related crowd crush is very variable. Simultaneity of trains also seems to have a growing role in the occurrence of such events.

1.3

Future work Our work showed that crowd crushes did happen in the (recent) past in train stations: this confirms they are definitely not risk-exempt. Therefore, we think it's needed to keep investigating those cases: run extensive analyses of past events but also gather information on near-misses cases. In the meantime, we plan to assess both the safety impact of situations other than crowd crushes, and the non-safety aspect of crowdrelated accidents (economical and image impacts, for instance). To improve the reliability of the crowd crushes list in the future, we also propose the collective use of Wikidata. The website organises data in a semantic way and also offers powerful exports tools, making easy to visualise data (Figure 1).

THERE IS A HUGE OPPORTUNITY IN BETTER KNOWING THE LEVERS

By pedestrian flow lever, we understand a management measure that can be implemented with reasonable means and aiming at a substantial effect. Using these levers, station managers can improve how the pedestrian flow works in the station without burning down all their available resources. They might be powerful tools, but they're insufficiently known. Only a few research papers and books, such as [START_REF] Campanella | Managing Large Flows in Metro Stations : The New Year Celebration in Copacabana[END_REF][START_REF] Feliciani | Introduction to Crowd Management Managing Crowds in the Digital Era : Theory and Practice[END_REF][START_REF] King | Using Simulation to Analyze Crowd Congestion and Mitigation at Canadian Subway Interchanges : Case of Bloor-Yonge Station[END_REF][START_REF] Moroni | La régulation ordinaire des flux de voyageurs en gare : L'infrastructure à l'épreuve permanente[END_REF][START_REF] Muñoz | Alleviating a subway bottleneck through a platform gate[END_REF], describe the levers that can be used in a transport-related context. The second part of the Observatory is then to collect as much precise information as possible on the existing implementations. To do so, we propose to use a basic questionnaire but also a detailed analysis grid. Results are made publicly available on the Observatory website (temporarily at the address mind-the-cap.github.io).

2.1

How to analyse a lever To fully understand how a lever works, we propose to use a detailed template (Dubroca-Voisin, 2023b). We divided it in several parts:

-Description of the lever (situation, operation and effects, activation and classification), to be able to understand how the lever is operated (by who and in which conditions, notably) and what place it takes in the station but also in the traveller's journey. Classification offers a first occasion to compare the levers by assessing them on different axes: physical vs virtual, influence vs coercion, easiness of activation, cost, and permanence. -Description if the station itself, to better understand the transport, mobility and urban context of the lever, but also the stakeholders and their interactions.

-Analysis of the lever (cost, gain, measure of performance, impact on risks). This section can include an analysis that has been made specifically when filling the template if the existing measures are not sufficient. -Historic of the lever (decision of implementation, first steps, adaptations and acceptance), to describe the story-behind-the-lever, that can be full of learning experience. -Replicability (similar implementations planned or done, expected ability to implement this lever in a different context). All these topics might need a long analysis time, but we're convinced it is needed to encourage the future dissemination of these levers. However, an alternative is to fill a shorter questionnaire that's available on the Observatory website; it covers mostly the same themes and still helps to understand how a lever works.

2.2

Presentation of levers Once they're obtained, these results are presented on the Observatory website. Figures 2 and3 show the first presentation choices that were; the first on a lever already evaluated by the literature, the second on a lever that was observed directly. This presentation may change to include more details and be clearer.

Future work

The primary goal is to include more and more levers to the collection. We could then compare the efficiency of these levers in different contexts. An additional possibility could be to collect levers that are not used in the railway system but in different spaces (such as museums, stadiums, public space…); these levers could potentially be adapted to the railway context.

ASSESSING THE SATURATION

In a third part of the Observation, we plan to observe the saturation of stations and their different elements. This is currently the less advanced part of the work; a first

Figure 1 :

 1 Figure 1: cartography of crowd crushes listed on Wikidata, September 2023.

FigureFigure 3 :

 3 Figure 2: First version of the presentation

© AET 2023 and contributors


Komlanvi Lampoh
email: komlanvi.lampoh@inrae.fr

Farhang Radjai

Claire Mayer-Laigle

Xavier Rouau

Jean-Yves Delenne

Charging kinetics of a vibrated bed of particles by numerical simulations and experiments

Keywords: triboelectric charge, electrostatic experiments, Discrete Element Method, granular materials PACS: 45.70.Mg, 61.43.Gt 2000 MSC: 76T25, 74F15, 78M25

Through particle dynamics simulations and experiments, we investigate the charging kinetics of a granular bed in a vibrating cell. The numerical method is based on a theoretical model combining an equation for charge transfer at the particle scale with a theoretical model of charge relaxation. Simultaneously, experiments are carried out to calibrate and validate the numerical model using a vibrated cell with glass beads. Specific charge per unit mass induced by vibrations initially increases before leveling off. Saturation specific charge and charging characteristic time decrease with increasing bed mass. Data collapse onto a master curve is observed when normalizing specific charges and times. At the particle scale, we find that charge saturation results from the balance between charge transfer and relaxation. Saturation specific charge correlates with the relative rate of collisions between particles and the vibrating plate. We introduce a simple expression for the average ballistic time based on single-particle dynamics and a restitution coefficient dependent on bed mass. The characteristic time is proportional to the ballistic time and decreases with increasing bed mass due to more inelastic collisions between particles.

Introduction

The handling of powders is crucial for the formulation of many products such as foods, pharmaceuticals, cosmetics, civil engineering materials, and particle embedded composites. In addition to the inherent complexity of granular media, electrostatic charges due to particle-particle or particlesurface contacts generate attractive and repulsive forces that deeply modify the mechanical behavior of powders. These interactions are at the origin of technological issues related to the handling of powders such as agglomeration effects, clogging and adhesion to surfaces. For example, electrostatic forces are critical in space missions since their effects are enhanced by low gravity and lack of moisture [START_REF] Lee | Adhesion and cohesion mechanisms of lunar dust on the moon's surface[END_REF][START_REF] Kawamoto | Practical performance of an electrostatic cleaning system for removal of lunar dust from optical elements utilizing electrostatic traveling wave[END_REF]. Electrostatic discharge can occur due to charge accumulation causing explosions of silos or flour mills [START_REF] Glor | Electrostatic ignition hazards in the process industry[END_REF][START_REF] Bailey | Incendivity of electrostatic discharges in dust clouds: the minimum ignition energy problem[END_REF]. On the contrary, these interactions can also be exploited in processes. Well known examples are electrostatic separation [START_REF] Gupta | Electrostatic separation of powder mixtures based on the work functions of its constituents[END_REF][START_REF] Landauer | Triboelectric separation of a starchprotein mixture -impact of electric field strength and flow rate[END_REF], coating [START_REF] Yang | Modeling of inter-tablet coating uniformity of electrostatic dry powder coating by discrete element method[END_REF][START_REF] Vanamu | An overview on dry powder coating in advancement to electrostatic dry powder coating used in pharmaceutical industry[END_REF], and powder inks for laser printers [START_REF] Hays | Paper documents via the electrostatic control of particles[END_REF].

Although all particle flows can give rise to triboelectric charges, the exchanged charge level strongly depends on the materials used and their relative position in the triboelectric series [START_REF] Park | Triboelectric series and charging properties of plastics using the designed verticalreciprocation charger[END_REF]. The charge migration process in the system is driven by the particle collision dynamics. For this reason, highly dynamic processes such as vibrated or fluidized beds, which are widely used in industrial processes for drying, mixing, and sorting, are strongly affected by electrostatic interactions [START_REF] Zhu | On the electrostatics of pneumatic conveying of granular materials using electrical capacitance tomography[END_REF][START_REF] Mehrani | Electrostatic charge generation in gas-solid fluidized beds[END_REF][START_REF] Liao | The effect of vibrating conditions on the electrostatic charge in a vertical vibrating granular bed[END_REF].

Triboelectric charging can occur by electron or ionic transfer or due to a transport of material. Several theoretical models have been introduced based on these different types of charge transfer. However, the underlying physical phenomenon remains poorly understood [START_REF] Bailey | Incendivity of electrostatic discharges in dust clouds: the minimum ignition energy problem[END_REF]. The most detailed approaches such as quantum simulations [START_REF] Yoshida | Experimental and theoretical approaches to charging behavior of polymer particles[END_REF] or models based on the distribution of electrons according to Fermi energy [START_REF] Cross | Electrostatics, principles, problems and applications[END_REF] directly consider the atomic scale but are hardly usable at the grain scale. Indeed, charging may also depend on grain scale characteristics such as their size, surface roughness, and mechanical loading [START_REF] Shi | Numerical study of the effect of temperature and h2o concentration on the electrostatic precipitator characteristics at high temperatures[END_REF][START_REF] Tilmatine | Effect of ambient air relative humidity on the triboelectric properties of polypropylene and polyvinyl chloride slabs[END_REF]. Furthermore, environmental factors such as humidity and temperature influence the characteristic time of return to electrostatic neutrality [START_REF] Kolehmainen | Effect of humidity on triboelectric charging in a vertically vibrated granular bed: Experiments and modeling[END_REF].

Several grain-scale phenomenological models have been proposed in the past [START_REF] Cheng | Investigation on electrostatic charging and its effect on mixing of binary particles in a vibrating bed[END_REF][START_REF] Grosshans | A model for the non-uniform contact charging of particles[END_REF]. These models involve the so-called work function providing the minimal energy necessary to gain or lose an electron, and were introduced for the first time in the seminal paper of Gallo et al. [START_REF] Gallo | Some charge exchange phenomena explained by a classical model of the work function[END_REF]. This work function should be calibrated experimentally. For example, the charging of a single sphere as a function of the number of collisions was investigated using various materials [START_REF] Matsuyama | Charge relaxation process dominates contact charging of a particle in atmospheric conditions[END_REF] and relative humidity and temperature [START_REF] Greason | Investigation of a test methodology for triboelectrification[END_REF].

To simulate the charging of a collection of particles, particle dynamics methods based on Discrete Element Models (DEM) can be efficiently used with the advantage of accounting for particle motions, interactions between particles, and collective effects such as arching, wall effects, force chains, and inelastic cooling. However, although most DEM-based simulations of electrostatic charging include a description of triboelectric charging at the particle scale based on the work function [START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF], they usually neglect the charge relaxation occurring during particle collisions. Various engineering applications were studied using DEM such as, for example, charge generation during mixing in vibrated granular beds [START_REF] Kolehmainen | Effect of humidity on triboelectric charging in a vertically vibrated granular bed: Experiments and modeling[END_REF][START_REF] Lu | Mixing in vibrated granular beds with the effect of electrostatic force[END_REF], triboelectric charging of particles rolling on inclined planes [START_REF] Hogue | Discrete element modeling (dem) of triboelectrically charged particles: revised experiments[END_REF], and the pneumatic conveying [START_REF] Watano | Numerical simulation of electrostatic charge in powder pneumatic conveying process[END_REF].

In this paper, we focus on the modeling of triboelectric charging of a bed of particles subjected to vertical vibrations. This configuration is particularly important for the characterization of charging efficiency of a powder bed in contact with a plate of given material. We introduce a model of pairwise charge exchange between particles, which relies both on the equation of Laurentie et al. [START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF] for charge transfer and on the theoretical framework of Matsuyama and Yamamoto [START_REF] Matsuyama | Charge transfer between a polymer particle and a metal plate due to impact[END_REF] for charge relaxation. We also calibrate and validate our numerical setting with experiments in bed charging under vertical vibrations. As we shall see, charge relaxation plays a key role in triboelectric charging. We apply this model to investigate the evolution of the specific charge for increasing bed mass and we analyze the origins of the saturation charge and and charging rate at the particle scale.

In the following, we first briefly describe in Section 2 the principles of the implementation of our DEM-based model for electric field and triboelectric charging. Then, we present in Section 3 our experimental electrostatic charging setup and tests performed on granular beds of different masses in order to calibrate the physical parameters of the numerical model. Section 3 describes the results of experiments performed to validate the charging kinetics for different masses of granular beds. A detailed numerical parametric study is presented with a broad range of values of granular bed mass. The calibrated numerical model is then used in Section 4 to investigate the charging kinetics for increasing bed mass. In Section 5, we introduce a scaling of the charging curves and a simple model based on single-particle assumption to clarify the role of particle-particle and particle-plate collision dynamics on the charging regimes. We conclude with a summary of the main findings of this work and its possible extensions.

Modeling of triboelectric charging of particles

In this section, we describe our numerical model of particle charging. We first detail the theoretical model of charge transfer at the contact scale. Then, we propose an implementation of this model in the framework of the code that we developed for the simulations based on the Discrete Element Method (DEM).

Charge transfer at contact scale

Matsuyama and Yamamoto [START_REF] Matsuyama | Charge transfer between a polymer particle and a metal plate due to impact[END_REF] proposed one of the first models of particle electrification due to contact and collisions. In this model, it is assumed that charge transfer between two particles i and j occurs only when they undergo a mechanical contact and it depends on a characteristic length δ c called tunneling distance, which defines their maximum interaction distance. The particle pair behaves as a capacitor in which the electric charges are transferred by tunneling effect from one surface to the other. For simple shapes such as spheres, the corresponding surface of transfer S is explicitly known through the contact force law. To account for material properties such as chemical composition, particle size, and surface state, Matsuyama and Yamamoto used the so-called work function introduced in the earlier seminal paper of Gallo et al. [START_REF] Gallo | Some charge exchange phenomena explained by a classical model of the work function[END_REF]. The work function φ is the minimal energy necessary to extract one electron from an initially chargeless surface. It was theoretically evaluated for spherical particles [START_REF] Gallo | Some charge exchange phenomena explained by a classical model of the work function[END_REF]. In experiments, φ can be determined by photoemission [START_REF] Lang | Theory of local-work-function determination by photoemission from rare-gas adsorbates[END_REF] and Atomic Force Microscopy [START_REF] Kwon | Probing of surface potential using atomic force microscopy[END_REF].

In the absence of an external electric field, for two particles i and j with different materials, the increment of charge gained or lost per unit time during a collision is [START_REF] Matsusaka | Triboelectric charging of powders: A review[END_REF][START_REF] Matsusaka | Electrification of an elastic sphere by repeated impacts on a metal plate[END_REF] ∆q

ij ∆t = C ij (φ i -φ j ), (1)
with [START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF]

C ij = r 0 4πδ c ∆S ij ∆t , (2)
where ∆S ij /∆t is the contact area change rate, e is the absolute value of the electron charge, 0 is the dielectric constant of void, and δ c is the tunneling distance between two particles for charge transfer. This model makes it possible to predict the evolution of the electric potential as a function of particle positions and to take into account the effect of their initial charge.

N P b N j u 1 N c Q t a f 8 w q t 3 i D g = " > A A A C y 3 i c j V H L S s N A F D 3 G V 6 2 v q k s 3 w S K 4 K o k U d V l 0 4 0 a o Y B / Q l j K Z T u v Q N A n J R K j V p T / g V v 9 L / A P 9 C + + M K a h F d E K S M + e e c 2 f u v V 7 k y 0 Q 5 z u u c N b + w u L S c W 8 m v r q 1 v b B a 2 t u t J m M Z c 1 H j o h 3 H T Y 4 n w Z S B q S i p f N K N Y s J H n i 4 Y 3 P N P x x o 2 I E x k G V 2 o c i c 6 I D Q L Z l 5 w p o p r t n v A V 6 / J u o e i U H L P s W e B m o I h s V c P C C 9 r o I Q R H i h E E A i j C P h g S e l p w 4 S A i r o M J c T E h a e I C 9 8 i T N y W V I A U j d k j f A e 1 a G R v Q X u d M j J v T K T 6 9 M T l t 7 J M n J F 1 M W J 9 m m 3 h q M m v 2 t 9 w T k 1 P f b U x / L 8 s 1 I l b h m t i / f F P l f 3 2 6 F o U + T k w N k m q K D K O r 4 1 m W 1 H R F 3 9 z + U p W i D B F x G v c o H h P m x j n t s 2 0 8 i a l d 9 5 a Z + J t R a l b v e a Z N 8 a 5 v S Q N 2 f 4 5 z F t Q P S + 5 R q X x Z L l Z O s 1 H n s I s 9 H N A 8 j 1 H B O a q o m T k + 4 g n P 1 o W V W L f W 3 a f U m
In the presence of an electric field E c due to particle charges or an external potential, following the work of Laurenti et al [START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF], the model of charge transfer can be written as:

∆q t ij = ∆S ij ε 0 δce (∆φ ij -E c • n ij δ c e) if ∆S ij > 0, ∆q ij = 0 if ∆S ij ≤ 0, (3)
in which n ij is the normal unit vector from particle j to particle i. For spheres, the contact surface area can be approximated as [START_REF] Hertz | Ueber die berührung fester elastischer körper[END_REF]:

S ij = πR δ n ij (4)
where R = 2 1/R i +1/R j . Hence, the variation of contact area is proportional to that of δ n .

In addition to charge transfer due to potential difference, we need a model for charge relaxation due to micro-discharges (electric arcs). In the literature, this effect has been modeled in different ways. For example, Chen et al. [START_REF] Chen | Tribocharging properties of wheat bran fragments in air-solid pipe flow[END_REF] considered that the electric charge undergoes an exponential decrease with time q(t)/q ∞ = 1 -e -(1/τ +β)t where q ∞ is the maximal sustainable charge, and τ and β are relaxation parameters describing the charge transfer and the effect of environmental conditions, respectively. Mastsuyama and Yamamoto [START_REF] Matsuyama | Charge relaxation process dominates contact charging of a particle in atmospheric conditions[END_REF] assumed that micro-discharges occur when the gap between two surfaces is close to the tunneling distance δ c . This process is described by Paschen's curve, which gives the maximum potential that can sustain the surrounding gas before a discharge occurs, as shown in Fig. 1. In this figure, we assume that when a particle impacts a wall (or two particles collide), as illustrated in Fig. 2, different steps should be distinguished:

1. A particle of charge q p = q 0 (where q 0 is the initial charge) is approaching a wall. As the difference of potential is proportional to the difference of charges between the particle and the wall, ∆V is independent of δ n at large δ n . However, as the gap declines, the potential decreases and eventually vanishes when the particle touches the plane [START_REF] Go | A mathematical model of the modified paschen's curve for breakdown in microscale gaps[END_REF]. 2. During collision, charges are exchanged between the two materials, leading to a sudden increases of |q p |. After collision, the particle moves away and the potential ∆V between the particle and the plane increases. 3. When ∆V reaches the Paschen's limit curve, the charge relaxes due to micro-discharges. 4. The particle moves away with a residual charge, which can be above the initial charge (Fig. 2). Note that the minimum of the Paschen curve occurs at the distance δ c .

The charge relaxation occurs over the characteristic distance δ c . In the model of [START_REF] Matsuyama | Charge relaxation process dominates contact charging of a particle in atmospheric conditions[END_REF] the relaxation due to charge transfer through gas molecules is assumed to be negligible compared to electric arc discharges. In this work, we assume that charge transfer due to relaxation takes place only during the separation phase; see Fig. 2c. We use a linear relaxation rule with the parameter K r defined by ∆q r ij = -K r (q i -q j ).

(

) 5
Although experimentally the charge density at the surface of a particle may be inhomogeneous [START_REF] Grosshans | A model for the non-uniform contact charging of particles[END_REF], considering a nonuniform charge distribution at the surface of each particle would lead to prohibitive computation time in 3D. For this reason, in this work we assume that particles are uniformly charged. For a particle i with contacts j ∈ C i , where C i is the set of particles in electric contact with i, the charge of the particle at time t + dt is given by

q i (t + dt) = q i (t) + j∈C i (∆q t ij + ∆q r ij ) (6)

DEM implementation

For the simulations, we used the Discrete Element Method (DEM) [START_REF] Radjai | Discrete-element modeling of granular materials[END_REF][START_REF] Luding | Introduction to discrete element methods[END_REF][START_REF] O'sullivan | Particulate Discrete Element Modelling: A Geomechanics Perspective[END_REF] by adding the electrostatic forces to frictional contact interactions between particles. For time integration of the equations of motion of the particles the velocity-Verlet scheme was used [START_REF] Allen | Computer simulation of liquids[END_REF]. At each time step, the forces and moments applied to all particles are calculated from particle positions and force laws. Then, using Newton's second law, the translational and angular accelerations of each particle are evaluated. Finally, using these accelerations and current velocities, the positions and velocities are updated.

The mechanical contact force is

F c = F c n n + F c t t,
where n and t are unit vectors defining the local reference frame at the contact between the two particles and F c n and F c t are its normal and tangential components; see Fig.

The normal component is given by the Hertz model

F c n = k n δ 3/2 n + η n v ij • n, (7)
where

δ n is the overlap distance, v ij = v i -v j , k n = 4 3 E R 1/2
with R and E defined as the harmonic means of the radii and the moduli of the two particles, respectively, and η n is the normal damping coefficient.

The tangential component is computed using a regularized Coulomb law

F c t = min(k t δ t -η t || v s ||, µ|| F c n ||), (8)
where δ t is the cumulative tangential displacement, η t is the tangential damping coefficient, and

v s = ( v ij -( v ij • n ij ) n ij ) + (R i Ω i + R j Ω j ) ∧ n (9)
is the sliding velocity and

k t = 8G R 1/2
ij is the tangential stiffness with G defined as the effective shear modulus computed using an harmonic mean of the shear moduli of the two particles. In addition to the mechanical contact forces acting on each particle i, the latter is subjected to gravity and electrostatic forces. The electrostatic force F e i = E i q i is calculated using the electric field at the center of mass of particle i:

E i = 1 4π 0 j =i q j r ij 2 n. (10)
The damping coefficients η n and η t can be scaled by the characteristic viscosities √ m k n and √ m k t , respectively [START_REF] Richefeu | Approche par éléments discrets 3d du comportement de matériaux granulaires cohésifs faiblement contraints[END_REF]. Therefore, we set

η n = 2α n m k n ; η t = 2α t m k t (11)
The values of α n and α t are related to the normal and tangential restitution coefficients e n and e t , respectively [START_REF] Ting | Computational laboratory for discrete element geomechanics[END_REF]:

α n = ln(e n ) ln 2 (e n ) + π 2 ; α t = ln(e t ) ln 2 (e t ) + π 2 (12)
Finally, in all simulations we set the time step ∆t = 10 -6 s. With this value, the maximum gap δ n in the system is below 3 × 10 -3 d and each collision involves at least 15 time steps. The CPU time for the simulation of charging of several hundred particles up to the maximum charge Q ∞ lasts from a few hours to four weeks on a moderate single-CPU station, depending on the number of particles. The simulation time is considerably longer here than in purely mechanical simulations (without electrostatic forces) mainly due to the calculation of the electrostatic field from particle charges in each time step. It is noteworthy that due to the strong inhomogeneity of particle charges during charging, it is not physically possible to cut off the electrostatic field below a few particle diameters.

Charging of a bed of particles

In this section, we present the principles of well-controlled experiments that we performed to calibrate and validate our triboelectric model. We first describe an in-house setup used to measure the evolution of the charge of a granular bed fluidized by means of a vibrating plate. These experiments together with DEM simulations will allow us to evaluate the tunneling distance δ c and the relaxation coefficient K r .

Experimental setup

The experiments are similar to those carried out by Liao et al. [START_REF] Liao | The effect of vibrating conditions on the electrostatic charge in a vertical vibrating granular bed[END_REF]. It consists of a bed of mass m b of single-size glass beads placed on a a vibrating plate of given material, as schematized in Fig. 5(a). The particles are confined by a cylindrical box which vibrates with the plate. In this work, for all experiments we used polypropylene for the plate material. As compared to aluminium and PVC that we tested, polypropylene provides an efficient charge transfer since polypropylene and glass are well-separated in the triboelectric series based on the work function [START_REF] Pan | Fundamental theories and basic principles of triboelectric effect: A review[END_REF][START_REF] Zou | Quantifying the triboelectric series[END_REF]; see Fig. 4. Table 1 summarizes the material parameters of the glass beads and the polypropylene container.

The friction coefficient and normal restitution coefficient between the beads and with the container were determined experimentally. For friction a tripod of three beads glued on a piece of steal was used. At the beginning of the experiment, the tripod is placed on a glass or polypropilene plane. Then, the plane is gradually tilted until the tripod slides. From this angle the friction coefficient was found to be µ b = 0.17 ± 0.02 for bead-bead and µ p = 0.23 ± 0.03 bead-polypropilene contacts. This value of bead-bead friction coefficient is close to that (µ b = 0.16) measured elsewhere [START_REF] Lorenz | Measurements of impact properties of small, nearly spherical particles[END_REF][START_REF] Tang | Measurement of restitution and friction coefficients for granular particles and discrete element simulation for the tests of glass beads[END_REF]. Regarding the restitution coefficient, we used the classical method of free fall of a bead on a surface. We recorded the drop of particles using a highframe camera and we determined the drop and rebound heights. We found a restitution coefficient of e b = 0.91 ± 0.02 for glass bead on glass surface and e p = 0.6 ± 0.02 for glass bead on polypropilene surface.

Let us note that all parameters required for DEM simulations using the electrostatic charging model are provided experimentally except the tunneling distance and relaxation coefficient, which are generally chosen arbitrarily in the literature [START_REF] Naik | Triboelectrification: A review of experimental and mechanistic modeling approaches with a special focus on pharmaceutical powders[END_REF][START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF]. As we shall see below, we use the experimental data as compared with numerical simulations to adjust these parameters. Glass beads [START_REF] Lorenz | Measurements of impact properties of small, nearly spherical particles[END_REF][START_REF] Tang | Measurement of restitution and friction coefficients for granular particles and discrete element simulation for the tests of glass beads[END_REF] The box is shaken using a Brünel & Kjaer 4810 mini-shaker at a frequency ν =20 Hz with an amplitude A =1.5 mm. For these values, the relative acceleration is A(2πν) 2 /g 2.4, leading to the fluidization of the granular bed across all its thickness. The vibration of the box is controlled by an amplified sinusoidal signal. The vibration amplitude was calibrated by means of a laser position sensor.

The total mass m b of the beads and the vibration time t v were varied and the bed charge was determined using and an electrometer and a Faraday cage; see Fig 5b). For each value of m b and t v , at least three independent measurements were performed. All experiments were carried out during a short period of time of the order of a few hours to avoid potential alteration of the results due to variations of humidity and temperature. The atmospheric conditions during experiments were a temperature of 26.5 ± 0.1 • C and a humidity of 60 ± 0.5%.

Figure 6 shows a typical example of the evolution of specific charge Q (charge per unit mass) as a function of time for m b = 1 g. We see that Q first increases almost linearly at a rate Q0 , and then slows down and tends to a steady-state value Q ∞ in a characteristic time τ . This steady-state specific charge is related to the maximum charge Q max that a spherical particle can carry before the electrostatic field around the particle reaches the "breakdown field". It corresponds to the yield potential beyond which the air molecules are ionized and the particle charge is removed by the air charged ions. This theoretical limit can be approximated as Q max 2.64 10 -5 4πR 2 p, where p = 3 r /( r + 2) [START_REF] Cross | Electrostatics, principles, problems and applications[END_REF]. Using this semi-empirical relation, we find Q max 68 nC g -1 which is of the same order of magnitude as Q ∞ 21 nC g -1 that we measured for a layer with m b = 1 g of beads of diameter of 2 mm. In the following, we will use the theoretical value Q max = 68 nC g -1 to normalize the specific charge.

To compare the experimental measurements with the numerical data, we apply numerically the same boundary conditions and vibrations in both cases. Regarding the tunneling distance and relaxation coefficient, a specific calibration procedure was necessary and it is described in the following section.

Determination of tunneling distance and relaxation coefficient

Little information can be found on the determination of the tunneling distance δ c and relaxation coefficient K r . The measurement of δ c has been shown to be sensitive to surface properties such as roughness [START_REF] Naik | A combined experimental and numerical approach to explore tribocharging of pharmaceutical excipients in a hopper chute assembly[END_REF][START_REF] Duff | Particle dynamics simulations of triboelectric charging in granular insulator systems[END_REF]. Arbitrary values have been used for different powders, such as δ c = 500 nm for millimeter-sized polyamide particles [START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF] and δ c = 250 nm for a pharmaceutical powder [START_REF] Naik | A combined experimental and numerical approach to explore tribocharging of pharmaceutical excipients in a hopper chute assembly[END_REF]. Regarding K r , to our best knowledge, no data are available since this parameter is not generally taken into account.

Figure 7 displays the evolution of the specific charge as a function of time for 1 g of particles and different values of δ c and K r . For a fixed value K r = 10 -6 s -1 we varied gradually δ c from 10 nm to 400 nm for 8 different values of δ c . Ten values of K r were also tested in the range from 10 -7 s -1 to 5 × 10 -6 s -1 for δ c = 2 × 10 -7 m. Fig. 7(a) indicates that both the charging rate Q0 and saturation charge Q ∞ are strongly influenced by δ c . On the contrary, Fig. 7(b) shows that K r only controls the saturation charge, which decreases with increasing values of K r .

Figure 8a shows the evolution of charging rate Q0 as a function of 1/δ c . Interestingly, Q0 is proportional to 1/δ c . This is consistent with Eq. ( 3) for weak values of E c at the beginning of the charging process. Hence, in the (a) following we will assume

Q0 (δ c ) 1.50 10 -7 /δ c , (13)
independently of K r . From this equation and the experimental data of Fig. 6 (established for a granular bed of the same mass as as in simulations), we get δ c 2 10 -7 m. For this value of δ c , Fig. 8(b) shows that Q0 declines only very slightly with K r so that its value can be considered to be practically constant.

We also consider the evolution of saturation charge Q ∞ (δ c , K r ) as a function of δ c shown in Fig. 9. We find a clear power-law dependence for both Q ∞ (δ c , K r = 10 -6 ) and Q ∞ (δ c = 2 10 -7 , K r ). Using the experimental value of Q ∞ = 21.4 nC for m b = 1 g (Fig. 6) and the data of Fig. 9b, we get K r 6 10 -7 s -1 .

Effect of the bed mass

With the parameters calibrated in the case of m b = 1 g, it is possible to simulate the charging of any amount of particles. Fig. 10(a) shows that our numerical simulations with the calibrated values of all parameters compare very well with the experimental data for bed masses of 2 g, 3 g and 4 g. Additional experiments were performed to evaluate the saturation charge Q Exp to a small systematic error, which may be attributed to the measurement method of the mass or saturation charge in experiments. Finally, it is worth noting that, contrary to the work of Jari [START_REF] Kolehmainen | Triboelectric charging of monodisperse particles in fluidized beds[END_REF], who scaled the effective work function with mass, in our model the work function is an intrinsic parameter. However, its dependence on moister and temperature will have to be clarified in future works.

Charging kinetics

In this section, we use the calibrated numerical model described in the previous sections, to perform a detailed analysis of charging kinetics mainly as a function of the granular bed mass. Note that the influence of the bed mass has been considered in several previous studies, but no scaling has been so far proposed [START_REF] Kolehmainen | Effect of humidity on triboelectric charging in a vertically vibrated granular bed: Experiments and modeling[END_REF][START_REF] Laurentie | Discrete element modeling of triboelectric charging of insulating materials in vibrated granular beds[END_REF]. We performed a parametric study in which we varied the bed mass from 0.1 g to 8 g for a total number of 13 simulations. The number of particles varies from 100 to 1000. In each case, we simulated the charging kinetics up to 90 s. We first analyze below the charging process by distinguishing different mechanisms. Then, in the next section, we propose a general scaling of the charging kinetics. Figure 11 displays the charge ratio Q/Q max as a function of time for all our simulations. Some of the same data were partially used in Fig. 10 for comparison with experiments. We see that, as shown previously, the total charge in the granular bed increases almost linearly with time at the beginning of the test and then tends to a constant asymptotic value Q ∞ . As observed in Fig. 12(a), Q ∞ decreases slowly and in a nonlinear fashion with increasing bed mass m b . Note that, while the specific mass declines with increasing mass m b , the total saturation charge Q ∞ m b is an increasing function of m b . This nonlinear dependence of Q ∞ on m b seems to hold only for thin beds considered in this study since Q ∞ in Fig. 12 This evolution of the charging characteristics reflects the dynamics of collisions with the vibrating plate, which is at the origin of particle charging due to charge transfer. At the same time, particle-particle collisions tend to redistribute charges among particles by charge relaxation mechanism, but they have no effect on the total charge. The overall behavior depends therefore on the particle bed mass, which affects the balance of the two types of collisions and the amount of charge transfer and charge relaxation during each collision. To quantify these processes, we focus here on the collision rates and the transferred charges for each type of collision.

Figure 13(a) shows the evolution of the rates N pp and N pw of particleparticle and particle-plate collisions, respectively. Both N pp and N pw increase almost quadratically with m b due to the increase of the total number of particles. Since the collision rates depend on bed mass while the specific electric charge acquired by the bed is described per unit it is more convenient to consider the 'relative rates' defined by

K pp = N pp N pp + N pw K pw = N pw N pp + N pw . (14)
Figure 13(b) displays K pp and K pw as function m b . We see that particleplate collisions prevail at low bed but K pw declines as m increases and beyond m b = 5 g, K pp slightly exceeds K . This transition is consistent with that observed in Fig. 12(b). It is physically plausible since the number of particles that can interact with the vibrating plate per unit time is increasingly inhibited by collisions between particles as the total number of particles increases. Interestingly, the decrease of K pw with increasing m b is also linearly correlated with Q ∞ as displayed in Fig. 14:

Q ∞ Q max 0.33K pw (15)
In other words, lower values of the saturation specific charge correspond to lower values of K pw for larger values of m b . The total saturation charge varies therefore as

Q ∞ m b = 0.33Q max K pw m b .
It is also noteworthy that the average relative velocity v p between particles and between particles and the vibrating plate declines with increasing m b due to enhanced dissipation by inelastic as shown in Fig. 15. Hence, we have two contradictory mechanisms that affect the charge transfer during collisions. The increase of N pw implies more charge transfer whereas the decrease of v p leads to lower of charge at each collision with the plate.

Figure 16(a) shows the average value ∆Q t pw b the transferred in a collision between a particle and the vibrating plate (bottom wall) as a function of time for different values of bed mass m b . The average charge per collision is calculated over a mobile time interval of 3 seconds. We see that ∆Q t pw b reaches a constant value in a short time interval and its steadystate value declines as m b increases. This is consistent with the decrease of the relative velocity v r with increasing m b due to increased dissipation as observed in Fig. 15. This gain of charge is reduced on average by an amount equal to ∆Q r pw b due to relaxation. Figure 16(b) shows ∆Q r pw b as a function of time. The variation of ∆Q r pw b in absolute value with time is similar to that of ∆Q t pw b . The net amount of charge gained by a particle varies therefore as ∆Q pw b = ∆Q t pw b + ∆Q r pw b . In the same way, the collisions with the lateral walls lead to an average gain of charge per collision ∆Q pw h = ∆Q t pw h + ∆Q r pw h by particle due to charge transfer ∆Q t pw h and loss of charge ∆Q r pw h by relaxation. Figure 16(c) displays ∆Q t pw h as a function time. It shows that charge transfer at the lateral walls is an order of magnitude below that at the vibrating plate due to the much lower values of relative velocity at the lateral walls. Furthermore, this transfer occurs only during the transient (at most 20 s) and then drops to zero. Hence, the total average charge gained by a particle per collision ∆Q can be approximated by

∆Q ∆Q pw b ∆Q t pw b + ∆Q r pw b . (16)
Figure 16(d) shows ∆Q as a function time. We see that ∆Q is initially high, but falls off and vanishes in finite time in the steady state as a result of balance between gain by collisions with the vibrating plate and loss by relaxation. This falloff is exponentially fast as shown by the fits to the data for smallest and largest bed masses m b in Fig. 16(d). The faster decrease of ∆Q with increasing m g is consistent with the observed increase of Q0 in Fig. 12

Hence, the faster increase of N pw /m b overcompensates the decrease of ∆Q with increasing bed mass and leads to the increase of Q0 with m b .

Scaling of charging kinetics

The evolution of specific charge Q with time in Fig. 11 shows that both the charging time τ and the saturation specific charge Q ∞ depend on the bed mass. In the last section, it was shown that this evolution is mainly controlled by the collisional dynamics of the granular bed. This suggests that the underlying charging process is similar for all values of the bed mass. This implies that the data must be scalable by normalizing the specific charge by Q ∞ and the time t by τ . Figure 17 shows that the normalized data do collapse on a master curve. This scaling does not fully capture the crossover from the nearly linear regime of rising charge to the saturation charge, which occurs at time t τ , but, given the low number of simulated particles, it is remarkably good. According to Eq. ( 15), the decrease of Q ∞ towards a constant value with increasing bed mass is correlated with that of the relative rate of collisions K pw between the particles and the vibrating plate. However, the physical origin of the rise time τ is not obvious. We introduce here a simple model to quantify the time scales involved in the vibrated bed. Fig. 18 displays a snapshot of particles and their color-encoded charges after 20 s of vibration. We see that, although the particles collide both with the bottom plate and with other particles, the packing fraction is low. Hence, we assume that the dynamics is essentially collisional (for the values of vibration frequency and amplitude used in our simulations). We consider therefore the dynamics of a single particle and we account for collisions with other particles by introducing an effective restitution coefficient.

For a plate vibrating vertically with amplitude A and frequency ω = 2πν, the contact of the particle with the plate opens only if aω 2 > g. When this condition is fulfilled, the velocity of the particle initially placed on the vibrating plate at the separation time is v 0 = Aω. The particle rises while decelerating, reaches its highest position ∆H, and falls back into the plate with a velocity v 0 . The collision with the plate may occur while the plate is moving upward or downward. However, the largest relative impact velocity 2v 0 occurs if the particle hits the plate while the latter is at its maximum upward velocity v 0 . We assume that this is the case in order to calculate the maximum possible height of the particle. Given the restitution coefficient e p between the particle and the plate, the relative velocity after the impact is 2v 0 e p and the velocity in the laboratory frame is 2v 0 e p + v 0 . More generally, the velocity v k (in the laboratory frame) of the particle at the separation time after k impacts is given by

v k = (v 0 + v k-1 )e p + v 0 , (18)
where v k-1 is the particle velocity at impact k -1. By iterating this relation, we obtain the following expression for the particle velocity:

v k = 1 + e p -2e k+1 p 1 -e p v 0 . (19)
Hence, since e k+1 p tends to zero with k, the largest steady-state velocity is

v max = Aω 1 + e p 1 -e p . ( 20
)
As a result, the corresponding maximum height reached by the particle is

∆H = v 2 max 2g = A 2 ω 2 2g 1 + e p 1 -e p 2 . (21)
Figure 19 shows the increase ∆H = H -H 0 of the bed height H from its initial height H 0 as a function of bed mass m b . For the lowest mass m b = 0.1 g, where the bed contains only a single dilute layer of particles, ∆H is nearly equal to the value ( 0.029 m) given by Eq. ( 21) for the numerical values of A = 1.5 × 10 -3 m, ω = 2πν with ν = 20 Hz, g, and e p = 0.6 used in the simulations. This indicates that single-particle assumption is relevant in the limit of low bed mass. However, ∆H declines as m b increases. This decrease is obviously due to energy dissipation by multiple particle collisions. A simple way to account for this additional dissipation (in addition to dissipation due to incomplete restitution at the vibrating plate) consists in keeping the expression [START_REF] Gallo | Some charge exchange phenomena explained by a classical model of the work function[END_REF]

where m 0 is a characteristic mass and α > 0 is a parameter depending on the dynamics. This collective restitution coefficient is equal to e p when m b m 0 and tends to zero when m b m 0 . Hence, the bed expansion can be expressed as

∆H(m b ) = A 2 ω 2 2g 1 + e p (1 + m b /m 0 ) -α 1 -e p (1 + m b /m 0 ) -α 2 . ( 23
)
This expression is used in Fig. 19 to fit the simulation data. We see that the fit closely follows the data points with m 0 = 5 g and α = 0.75. With this 'calibration' of the expression of ∆H(m b ) from simulation data, we also can express the ballistic time t b , i.e. the time taken by the particle to fall a distance ∆H(m b ):

t b = v max g = 2∆H g 1/2 = Aω g 1 + e c 1 -e c . (24)
Since the onset of fluidization implies the condition Aω 2 > g, we have Aω/g > 1/ω. Hence, the minimum value of t b occurs for e c = 0 (infinite mass):

t bmin = 1 ω . (25
Figure 20 shows the characteristic time τ versus t * with α = 0.75 and m 0 = 5 g for all values of the bed mass. We see that τ is proportional to t * with a prefactor 580. This means that, for all values of m b , the time required for the charging of the particles up to the saturation value corresponds to an average number of 290 collisions of each particle with the vibrating plate. Note that the total ballistic time of a particle (rise and fall) is 2t * . Indeed, we expect that this number is equal to τ N pw /N p , where N p is the number of particles. For example, for m b = 8 g, this ratio (with τ 10 s, N pw = 3 × 10 4 s -1 , and N p = 1000) is 300, which is very close to the above value calculated from the proportionality between τ and t * . The scaling of τ by t * suggests that the decrease of charging time with bed mass reflects the decrease of ballistic time. This assertion needs, however, further validation by means of simulations of charging for different values of vibration amplitude and frequency.

Conclusion

In this study, we employed a combination of particle dynamics simulations and experiments to explore the charging kinetics of a granular bed within a vibrating cell. Our approach involved developing a triboelectric charging model, incorporating a relaxation term, and conducting extensive charging simulations for granular beds with varying mass. To validate our numerical model, we calibrated its parameters through charging experiments using glass beads on a vibrating polypropylene plate. Remarkably, our simulation results for specific charge over time aligned well with experimental data across multiple bed masses. We used the simulations to analyze the charging process by considering the collision rates, relative velocities, and the transfer and relaxation of charges per particle.

A key observation was that the specific charge per unit mass increased with time before stabilizing at a constant value. This saturation was attributed to a balance between charge transfer from collisions with the vibrating plate and charge relaxation. Interestingly, the data showed also that saturation specific charge decreased with increasing bed mass due to a reduction in the relative collision rate between particles and the plate compared to inter-particle collisions. Another important observation was that the charging characteristic time decreased with higher bed masses. We introduced a simple model based on single-particle dynamics, accurately predicting bed expansion with a restitution coefficient dependent on bed mass. This model further forecasted a diminishing ballistic time in relation to particle mass, with the charging characteristic time effectively scaled by this ballistic time within the studied bed mass range. Due to the underlying dynamics of the charging process, our simulation data for specific charge over time remarkably collapsed onto a single master curve when specific charges were normalized by the saturation charge, and times were normalized by the characteristic time.

The numerical model developed in this study has the potential for application in exploring charging kinetics for larger bed masses, varying vibration amplitudes and frequencies, and different particle-plate restitution coefficients. Future investigations should assess the model's validity under different vibration conditions and explore its efficiency for simulating larger particle numbers. This expansion could provide insights into how particle dynamics and charging physics depend on system size.

  t e x i t s h a 1 _ b a s e 6 4 = " G D b U D q

s s 8 OFigure 1 :

 81 Figure 1: Paschen's limit curve describing the electric potential V necessary to initiate charge relaxation as a function of the gap δ n . The bold line illustrates variation of electric potential for a pair of particles before (a) and after (b-d) the creation of mechanical contact leading to a variation of charge. The dashed line represents the variation of electric potential after charge relaxation.

Figure 2 :

 2 Figure 2: Illustration of a particle (a) coming into contact with a wall and gaining electric charge during deformation (b). Micro-discharges occurs in the vicinity of contact for distances ranging from a few nanometers to micrometers (c) and move apart with a residual charge (d).

Figure 3 :

 3 Figure 3: The geometry of a contact between two particles.

Figure 4 :

 4 Figure 4: Typical example of triboelectric series adapted from [44].

Figure 5 :

 5 Figure 5: Schematic diagram of experimental setup with the vibrating system.

Figure 6 :

 6 Figure 6: Evolution of specific charge Q as a function of time, illustrating the measured quantities: charging rate Q0 , characteristic time τ , and asymptotic charge Q ∞ .

6 Figure 7 :

 67 Figure 7: Evolution of the specific charge Q for a particle bed of mass equal to 1 g as a function of time for different values of the parameter δ c (a) with K r = 10 -6 s -1 , and K r (b) for δ c = 2 10 -7 m.

∞(Figure 8 :Figure 9 :

 89 Figure 8: Charging rate for a granular bed of mass of 1 g (a) as a function of the inverse of tunneling distance 1/δ c with K r = 10 -6 s -1 and its linear approximation Q0 = 1.50 10 -7 /δ c ; (b) as a function of the relaxation parameter K r for δ c = 2 10 -7 m with its linear fit Q0 = -1.36K r + 0.75 nC/s.

Figure 10 :

 10 Figure 10: Comparison between experiments and simulations: (a) Specific charge (charge per unit mass) as a function of time for different values of bed mass; (b) Saturation charge as function of the bed mass m b . The error bars represent standard deviation over three independent measurements in experiments; (c) Saturation charge obtained from numerical simulations versus their experimental values.

Figure 11 :

 11 Figure 11: Specific charge Q normalized by Q max as a function of time for different values of bed mass.

Figure 12 :

 12 Figure11displays the charge ratio Q/Q max as a function of time for all our simulations. Some of the same data were partially used in Fig.10for comparison with experiments. We see that, as shown previously, the total charge in the granular bed increases almost linearly with time at the beginning of the test and then tends to a constant asymptotic value Q ∞ . As observed in Fig.12(a), Q ∞ decreases slowly and in a nonlinear fashion with increasing bed mass m b . Note that, while the specific mass declines with increasing mass m b , the total saturation charge Q ∞ m b is an increasing function of m b . This nonlinear dependence of Q ∞ on m b seems to hold only for thin beds considered in this study since Q ∞ in Fig.12(a) tends to a constant value with increasing m b . Figure 11 also shows that the initial charging rate Q0 per unit mass increases with m b . This is what we observe in Fig. 12(b), which displays the initial charging rate as a function of m b . The rate increases, but a transition

Figure 13 :Figure

 13 Figure 13: (a) Rates of collisions N pp and N pw between particles and with the vibrating plate, respectively, averaged over the whole period of vibration as a function of bed mass m b ; (b) Average relative rates K pp and K pw of particle-particle and particle-wall collisions as a function of bed mass m b .

Figure 15 :

 15 Figure 15: Relative velocity v p averaged over the whole vibration period as a function of m b .

Figure 16 :

 16 Figure 16: (a) Average charge ∆Q t pw b transferred per collision between particles and the vibrating plate, (b) Average charge loss ∆Q r pw b per collision between the particles and the vibrating plate due to charge relaxation, (c) average charge ∆Q t pw h transferred per collision between particles and the lateral walls, and (d) total transferred charge per collision between particles for different values of bed mass m b as a function of time. The solid lines in (d) are exponential fitting functions to the data points for the largest and smallest values of m b .

= 1 m b τ τ 0 N

 0 Figure16(d) shows ∆Q as a function time. We see that ∆Q is initially high, but falls off and vanishes in finite time in the steady state as a result of balance between gain by collisions with the vibrating plate and loss by relaxation. This falloff is exponentially fast as shown by the fits to the data for smallest and largest bed masses m b in Fig.16(d). The faster decrease of ∆Q with increasing m g is consistent with the observed increase of Q0 in Fig.12as a function of m b . During the transient time τ , we have Q0 = 1 m b τ

Figure 17 :

 17 Figure 17: Specific charge Q normalized by saturation specific charge Q ∞ as a function of time t normalized by characteristic time τ for all values of bed mass m b .

Figure 18 :

 18 Figure 18: Snapshot of a triboelectrically charged layer of particles after a vibration time t v = 20 s. Color level is proportional to particle charge.

Figure 19 :

 19 Figure 19: Bed height expansion ∆H as a function of bed mass m b from simulations and model.

Figure 20 :

 20 Figure19shows the increase ∆H = H -H 0 of the bed height H from its initial height H 0 as a function of bed mass m b . For the lowest mass m b = 0.1 g, where the bed contains only a single dilute layer of particles, ∆H is nearly equal to the value ( 0.029 m) given by Eq. (21) for the numerical values of A = 1.5 × 10 -3 m, ω = 2πν with ν = 20 Hz, g, and e p = 0.6 used in the simulations. This indicates that single-particle assumption is relevant in the limit of low bed mass. However, ∆H declines as m b increases. This decrease is obviously due to energy dissipation by multiple particle collisions. A simple way to account for this additional dissipation (in addition to dissipation due to incomplete restitution at the vibrating plate) consists in keeping the expression[START_REF] Gallo | Some charge exchange phenomena explained by a classical model of the work function[END_REF] but replacing the restitution coefficient e p by a 'collective'

)

  We can therefore define a time scale t * of fluidization by subtracting t bmin from t b :t * (m b ) = Aω g 1 + e p (1 + m b /m 0 ) -α 1 -e p (1 + m b /m 0 ) -α -1 ω .

Table 1 :

 1 Material parameters for glass beads and polypropilene box.

	Polypropilene box [50]

Acknowledgment

This work has been realized with the support of MESO@LR-Platform at the University of Montpellier. The authors also thank the PLANET facility (DOI: 10.15454/1.5572338990609338E12) run by the IATE joint research unit for providing process experiment supports.