Simulation of particle resuspension by wind in an urban system

Air pollution caused by particle resuspension is a growing public health problem in many cities. Pollen and anthropogenic pollutants, such as heavy metal particles and micro-plastics debris, settle onto urban ground surfaces. Prolonged urban heat waves are propitious for heavy and continuous deposition. Particles in the submillimeter size range eventually resuspend by urban winds within seconds, may be inhaled, cause allergic reactions and escape the city’s boundaries. Here, the resuspension and subsequent dispersion of generic particles ranging from 10 to 100 μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu$$\end{document}m in size are simulated. The city area “Bayerischer Bahnhof” in Leipzig, Germany, has been chosen as a practical example. To track the resuspended particles, a Lagrangian model is used. Taking advantage of graphics processing unit, turbulent flow simulations at different wind speeds are performed in almost real time. The results show that particle resuspension starts, when the inlet wind speed beyond the canopy, that is at a height of 40 m, exceeds 7 m/s. At wind speed beyond 14 m/s, resuspension occurs in almost all city parts. At moderate wind speed, high-risk areas are identified. The effect of green infrastructures on both the flow field and particle resuspension is also investigated. Extreme strong winds, even short in duration, cause the sudden resuspension of pollutants in the submillimeter size range The resuspension of generic particles ranging from 10 to 100 μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu$$\end{document}m in size is simulated As application, the city quarter “Bayerischer Bahnhof” in Leipzig, Germany, is chosen High-risk areas are located in city parts where the local wind orientation is parallel to narrow street axis In perimeter blocks, almost no resuspension occurs Extreme strong winds, even short in duration, cause the sudden resuspension of pollutants in the submillimeter size range The resuspension of generic particles ranging from 10 to 100 μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu$$\end{document}m in size is simulated As application, the city quarter “Bayerischer Bahnhof” in Leipzig, Germany, is chosen High-risk areas are located in city parts where the local wind orientation is parallel to narrow street axis In perimeter blocks, almost no resuspension occurs


Introduction
Air pollution by resuspending microparticles has become a serious public health problem in many cities and highly populated urban areas around the globe [1]. Heavy metal particles [2], micro-plastics debris such as tire and road wear particles [3], pollen [4], and anthropogenic fine particles are released in urban systems on a daily basis. Their continuous emission results in carcinogenic and cardiovascular diseases, asthma and skin allergies [5]. For these reasons, the investigation of urban flow and pollutant dispersion has received increasing attention in the recent years. Many of these pollutants initially settle on urban ground surfaces. This is particularly true during prolonged urban heat waves. Traffic, construction sites, trees and further pollutant sources exist in all parts of the city and contribute to heavy pollutant deposition [6][7][8]. Within seconds, particles in the submillimeter size range resuspend by strong urban winds and may be inhaled, cause allergic reactions and eventually escape the city's boundaries. It is hence important to understand and be able to predict particle resuspension by wind in an urban system, that is the re-entrainment of previously deposited particles and their subsequent dispersion.
Laboratory experiments can offer excellent insight into particle resuspension by wind. Yet, they are costly and difficult to perform under realistic conditions. To date, most of the experimental work focuses on determining the flow patterns in urban systems and have been performed for some simplified scenarios, often to validate flow simulations [9,10]. In an early work by Oke et al. [11], a theoretical investigation of the flow patterns in a street canyon was presented. By varying the aspect ratio of the building height H to the street width W, the team found out the critical ratio H/W, for which pollutant spread is most propitious. Thanks to the progress in hardware and computational technologies, simulations with Computational Fluid Dynamics (CFD) have gained momentum [12][13][14], especially for scenarios involving pollutant dispersion in street canyons [15]. Most of the work on pollutant dispersion in urban systems is however limited to a portion of an ideal city, often represented by blocks [16]. CFD data on larger city sections do exist, yet only a few studies exist because of the higher computational cost associated with large domains. Tseng [17] used a Large Eddy Simulation (LES) with a dynamic sub-grid model to investigate particle dispersion from a punctual source. The spatial resolution in their simulation was relatively coarse and only a few grid points were used to discretize the buildings with respect to their heights. Xie et al. [18] carried out LES to determine urban flow statistics in central London using a computational domain with a size set to 1.2 × 0.8 × 0.2 km 3 . They found that the mean and fluctuating concentrations in the near-source field was highly dependent on the source location and the layout of the buildings.
They also reported that a one-meter resolution was sufficient to obtain accurate prediction in flow and particle dispersion statistics. Using a dual-grid system, Liu et al. [19] simulated the wind flow field and pollutant dispersion from traffic emission in a portion of the city Macao. Their study highlighted the important contribution of buoyancy to pollutant transport in a street canyon. With the help of LES, Moonen et al. [20] investigated near-field pollutant dispersion from traffic exhaust in an urban street canyon. Their results emphasized the beneficial effect of avenue-like tree planting on public health. More recently, Hertwig et al. [21] simulated the near-surface atmospheric boundary layer for the city center of Hamburg in Germany. They compared eddy statistics and turbulence structures to experimental data obtained in an open-circuit boundary-layer wind tunnel. The scale of the city model in the wind channel was 1:350. Grid resolution was found to have an important effect on turbulence below rooftop. More recently, Lenz et al. [22] used the cumulant lattice Boltzmann method (LBM) for the simulations of urban air flows in Basel. Their work demonstrated the feasibility of real-time large-eddy simulations of flow over urban canopies at the neighborhood scale. All the above simulations mostly focus on either turbulence structures within urban canopies or on the dispersion of pollutant from source point or puff [23]. To the best of the authors' knowledge, there is little if any work on the dispersion of resuspending particles initially at rest on urban ground surfaces. This is exactly where this work comes in.
With the evolving role of extreme weather events, the frequency of low ground winds throughout the city, particularly during prolonged heat waves, will increase in the future [24]. Under such circumstances, particle deposition occurs. Extreme strong winds, even short in duration, eventually cause the sudden resuspension of these pollutants. Following the emerging trends of real-time simulations [22], this work introduces a fast numerical framework for the resuspension and dispersion of particles in the micron range distributed in a city section. The model links all spatial scales that range from the size of the particle to that of a city and incorporates green infrastructures. The contribution of our work hence is the coupling of a LBM solver with a resuspension model inspired from sediment resuspension in riverine flows [25]. Most studies on atmospheric dispersion of pollutants have so far not taken account the wall-pollutant interaction. As an example, we choose the urban quarter "Bayerischer Bahnhof" in Leipzig (Germany). The LBM has several advantages over classical Eulerian/Navier-Stokes solvers used so far. It needs no pressure correction and the locality of the operators allows to take full advantage of recent advances in parallel General Purpose Graphical Processing Units (GPGPU) for fast calculation [26]. The LBM can quickly simulate the turbulent particle-laden flow in almost real time. This tool could for instance be later used in smartphone applications and show high-risk areas associated with large resuspension rate. Here, the effect of green infrastructures, such as bush and grass, on urban flow change and particle resuspension is studied.
This manuscript is organized as follows. In the method section, the particle resuspension, particle-tracking and fluid models are explained. Afterwards, the numerical set up of the urban city part is introduced and the simulation results are discussed.

Urban flow model
An efficient Lattice Boltzmann Method (LBM) implementation is used for the numerical simulation of the turbulent urban flow. Our LBM code is based on the Boltzmann equation derived from a microscopic scale point of view [27]. The LBM discretizes the Boltzmann equation with a discrete velocity set, yielding a numerical method for computing macroscopic distribution functions on a Cartesian grid. The macroscopic aerodynamic quantities, such as pressure and velocity, are obtained as low-order moments of these distribution functions. It can be shown that, for sufficiently small discretization steps in space time, Mach and Knudsen numbers, the LBM solution converges to that of the classical Eulerian/ Navier-Stokes equations. The Boltzmann equation thus governs the fluid particle distribution functions (PDF), f (x, t, ) , which specify the probability of finding a fluid particle at position x and time t with velocity [27,28]. The Boltzmann equation is given by where Ω is the collision operator describing the interaction of the fluid particles. Discretized fluid particle velocities e ijk are introduced to yield a model of reduced computational cost. In this discretized formulation, a fluid particle is only allowed to move from a given lattice point, in a limited number of directions and for specific distances. With these assumptions, Equation (1) converts into a set of discrete Boltzmann equations as where e ijk = c × (i, j, k) and i, j, k ∈ (−1, 0, 1) , c = Δx∕Δt is the lattice speed, with Δx the regular lattice space step and Δt the time step. In the LBM, c is usually set to unity. This means that a fluid particle with velocity c travels over one lattice cell within a discrete time step. The definition of e ijk implies lattice with 27 discretized velocities. This so-called D3Q27 model, with 3 dimensions and 27 velocities, is illustrated in Fig. 1).
A finite difference discretization in space and time over one grid cell yields the lattice Boltzmann equation Equation (3) may be split up into two parts. The first part contains the collision step, for which the fluid particle distribution functions change from the equilibrium state. The equilibrium PDF is only function of the velocity and the pressure. The second part is the propagation step, in which the evolved fluid particle distribution functions are moved to the respective neighboring grid points. The collision step is a complete local operator and this LBM characteristic makes the numerical method highly suitable to take advantage of GPGPU for high-performance computing. The collision and propagation steps are implemented as where f ijk is the post collision particle distribution function. For the collision operator, a relatively new model, based on the so-called cumulant lattice Boltzmann equation, has been used. This cumulant method has some crucial advantages over typical LBM collision operators. The cumulant operator is Galilean invariant, more stable and produces smaller errors and noises. A comprehensive explanation of the cumulant method can be found in the authors' previous work [29]. The pressure and velocities are obtained from the zeroth and first order moments of the fluid particle distribution functions 3 is the speed of sound in the lattice. The present numerical approach is based on the cumulant LBM, introduced in [30]. It is an improvement of the cascaded LBM. In the cumulant LBM, the collisions occur in cumulant space instead of central moments. Cumulants are statistically independent and Galilean invariant variables calculated by mapping the moment space to the frequency space. The cumulant LBM is more stable numerically even for flows associated with very high Reynolds numbers. It also has smaller errors and is generally more accurate than the commonly used multi-relaxation time collision operator.

Particle resuspension model
Particle resuspension, that is the re-entrainment of particles initially at rest on the ground into the flow [31], is here assumed to occur when the friction velocity u * at the wall exceeds a critical velocity u * c [32,33], also called threshold or pick-up velocity. A general difficulty associated with the resuspension of micron particles in a large domain lies in evaluating these two velocities u * and u * c . Assumptions must be made to link the spatial scales of the urban system ranging from the size of the particle to that of a city. First, the friction √ w ∕ f , where w is the wall shear stress calculated from the velocity gradient at the wall and f the viscosity of the fluid, that here is air. In the presented simulations, the grid spacing used to simulate the urban flow is in the order of some meters ( Δx ≈ 1 m) and is thus by several orders of magnitude greater than the particle size ( d p ≈ 100 μm). This in-turn makes it difficult to provide an exact estimation of the friction velocity based on the velocity gradient at the wall. In the presented LBM fluid solver, the instantaneous shear stress tensor is evaluated directly using the non-equilibrium part of the distribution function [34] as where ̂= Δt∕Δx 2 is the non-dimensional kinematic viscosity in LBM units, the air viscosity, e ijk and e ijk the -th, respectively, the -th component of the vector e ijk at the lattice position (i, j, k) and the Kronecker delta. The fluid density f does not explicitly appear in Eq. (8) because it is already incorporated in the density function f ijk . In the presented simulations, the grid spacing used to simulate the urban flow is not fine enough for a precise wall shear stress calculation. For a better estimation of the shear stress tensor at the wall, that is ( ) z=0 , an extrapolation using the first three nodes departing from the wall is used. Note that, due to the bounce back scheme associated with the wall boundary condition, the wall locates at a distance Δx∕2 to the closest fluid node. The calculation of the shear stress tensor field could be further improved next to the wall by using a local grid refinement technique. With the z-axis chosen as direction normal to the ground, the wall shear stress is given by w = √ 2 13 + 2 23 . We stress, that, in the adopted resuspension model, the friction velocity is a locally computed value that is function of time and space and not a single value used for the whole domain. We had originally assumed a log-wind profile [35] near the ground fed by the instantaneous velocity profile. It turns out, this assumption is not valid for flows inside canopies [36,37], as is the case here. Assuming a log-wind profile within the canopy gave incorrect estimations of the friction velocity u * , that were about 50% higher than those obtained with Eq. (8).
A second numerical difficulty lies in estimating the critical velocity u * c . Previous experimental and theoretical works on dry resuspension [31,[38][39][40] have shown that the particle detachment depends on a number of parameters, such as particle and substrate material, particle and substrate roughness, particle shape and size, height of the vegetation, air humidity and temperature at the particle location. A change in any of the above parameters may affect the critical velocity u * c . In those above experiments, which involved sparsely distributed particles initially at rest on a substrate, u * c was defined as the friction velocity at which half of the particles resuspended. The particles were mostly spherical and made of glass, alumina, iron, polymer and had sizes typically ranging from 1 to 100 μ m. Hence, they appear suitable for typical particles on urban ground surfaces. Typically, resuspension occurs for critical velocities in the range of 0.2 m/s < u * c < 1.5 m/s when the ground surface is unvegetated [38,[41][42][43]. Another important finding was, that the critical velocity increases with particle size. For vegetated surfaces, experimental data are much more limited. The sward type, grass height and leaf density play an important role in particle resuspension. Giess et al. [44] studied the resuspension of 1, 10, and 20 μ m silica particles from short and long grassed surfaces in a wind tunnel. In short grass, that is vegetation with a height set to 10 cm, 10%, 27% and 36% of the 10 μ m particles resuspended at u * = 0.27, 0.64 and 1.14 m/s, respectively corresponding to wind speeds of 3, 5 and 7.8 m/s measured at a height of 60 cm above the soil surface. The 50% resuspension rate could not be determined in their work. By extrapolating, we here estimate u * c = 2.35 m/s. In tall grass, that is vegetation with a height set to 25 cm, only values for the 1 μ m particles were reported. The highest resuspension fraction, that could be obtained was 6% at u * = 1.15 m/s. In the following, we hence introduce a correction factor to include the effects of vegetation on resuspension. If the deposited particles are located on a vegetated ground surface, we set = 2. We are well aware, that this value is an approximation. The above experimental results suggest however that the critical velocity u * c doubles for an urban surface with grass.

Particle transport model
After a particle detaches from the surface, its transport is governed by the aerodynamic drag and the gravity. All other forces are neglected. Neither the particle-particle nor particle-surface collisions are considered. If a resuspended particle hits a building or the ground, it sticks to it. This assumption is uncritical with respect to horizontal ground surfaces because almost all resuspended particles do not return onto the ground after detachment. It is however a rather strong assumption for vertical walls. Sommerfeld et al. [45] investigated the rebound of 100 μ m glass beads and quartz particles off substrates made of steel, plexiglass and rubber. The restitution coefficient, that was found to be function of substrate, particle material and impact angle, varied between 0.4 and 0.9. Based on Maxey and Riley [46], further assumptions are made with respect to the particles. They are rigid and their volume do not change in time. They are spherical and no shape effects are considered. The particles are smaller than the Kolmogorov length scale. The particle Reynolds number defined by , so that the drag formulation determined for a small sphere applies. The term u f is fluid velocity at the particle location x p , u p the particle velocity and is the kinematic viscosity of the air. Based on this, the transport of each resuspended particles is governed by Newton's second law. The transport equation for each particle is hence given by where m p is the particle mass, m f the mass of the displaced fluid, D p the drag force acting onto a sphere and g the gravity. The drag is given by where A p is the cross-sectional area of the particle and C D the drag coefficient. C D is function of the particle Reynolds number Re p and is given by [47]. After simplifications, the final form of the particle transport equation becomes p ∕(18 f ) is the particle response time [48]. In Eq. (11), u f is the velocity of fluid at the particle position needed at every time step. In the LBM code the fluid velocity is only known at each grid node of the discretized domain. Hence an interpolation is 1 3 performed to determine the fluid velocity at the exact particle position. In our case, the fluid velocity field is obtained via a tricubic interpolation [49]. More details about this interpolation scheme can be found in the authors' previous work [50]. The new particle velocity is computed via a Taylor series expansion as The position of each particle is updated in time, again via Taylor series expansion, as The authors do not claim here new contributions to the LBM implementation. As pointed out in the introduction, most studies found in the literature focus on the dispersion of lowinertia pollutants from source point or puff. With its applied character, our study differs from others because pollutants originate from a large urban ground surface and are greater in size. Besides, using the LBM for wind simulation and coupling it with a Lagrangian particle-tracking algorithm allows for an almost real-time simulation of particle resuspension in a city part. Both the fluid and the Lagrangian particle solver benefit here from the locality of the operators. Hence, the simulations can be performed with relatively low-cost computational resources. In the present case, a single NVIDIA Tesla P100 was used. A simulation based on multiple GPUs would further reduce the execution time. Because of the increase in frequency of extreme events caused by climate change [51], studies on particle resuspension by strong winds are in fact expected to gain importance in the future. We believe this work will serve as basis for further investigations on particle resuspension in urban systems.

Grid convergence study and LBM validation
The purpose of this section is twofold. First, the LBM flow solver is validated against literature data. Second, a grid convergence study is performed. We point out, that in the fundamental context of turbulent particle-laden flows at high Reynolds, our method has already been validated several times in [29,50]. In the present context however, that is urban fluid mechanics [52], wind velocity profiles obtained in a wind-tunnel [53] and their respective counterparts obtained with a classical Eulerian/Navier-Stokes solver [54] are used as reference data. In the original experiment, an idealized city is exposed to a turbulent boundary-layer flow and was, at that time, performed at the Environmental Wind Tunnel Laboratory of the Meteorological Institute of Hamburg University. The idealized city portion, referred to as the Michel-Stadt (case BL3-3), is illustrated in Fig. 2a. With this experimental set-up, the flow is physically downscaled to 1:225. At full scale, length and width of the city portion are 1320 × 830 m 2 and the building heights are 15, 18 and 24 m. This semi-idealized urban geometry exhibits typical features reminiscent of Northern and Central European cities, such as courtyards, oblique road arrangements, perimeter blocks and open spaces. The urban quarter "Bayerischer Bahnhof", which is later presented, exhibits many of these features as well. Further reading on the original experiment may be found in [21,53,54].

3
Our corresponding full-scale three-dimensional LBM model of the city is shown in Fig. 2b. A no-slip boundary condition is enforced onto the ground and building surfaces. On the two lateral domain sides, a periodic boundary condition is used. The top boundary is set to free slip and a zero gradient is used for the outlet boundary condition. A second order zero gradients is used at the outlet for the unknown particle distribution functions, that is with L x the domain size in the streamwise direction and Δx the grid spacing. On the inlet surface, the boundary condition is enforced as where u x is the streamwise component of the velocity vector u = (u x u y u z ) ⊤ , the power index which is function of the urban roughness, that is the terrain type. For open land is approximately 0.143. For flows over relatively small buildings, as is the case here, = 0.18 is recommended and used throughout this work [55]. The velocity u w = 6.84 m/s is the free stream wind speed and z h = 60 m is the vertical distance from the ground beyond which the velocity profile becomes more or less uniform [55]. Three LBM grids, here termed coarse, medium and fine, with respectively 488×304×44, 590×368× 52 and 732×456× 64 grid points, are tested. The instantaneous and time-averaged velocity fields, passing through the domain midsection and obtained with the fine grid, are shown in Fig. 2c, d. For comparison purposes, we select two vertical segments passing through two points labeled Pt14 and Pt37 in the original experiment (See Fig. 2b). These two points are on the flow centerline of the domain. The time-averaged velocity profiles ũ x are shown in Fig. 3a, b along these two selected positions. Overall, the results from our LBM simulations match relatively well those taken from the literature. Good results are achieved with the fine grid in both the low altitude region ( z∕L z < 0.25 ), where particle resuspension occurs, and in the upper region Here L z is the domain height. In the intermediate region, that is 0.25 < z∕L z < 0.5 , a deviation is observed. This is attributed to the inlet boundary condition, which does not exactly correspond to that in the original flow, especially in terms of inlet velocity modeling. Tolias et al. [54] used a Langevin-type inflow boundary condition. The first point, Pt14, is located close to the channel inlet, upstream of a relatively large open area, that is canopy-free. Near the ground the streamwise velocity is negative. With higher altitude, it eventually becomes positive and gradually increases throughout. The second point, Pt37, is further downstream of the domain in a narrow street canyon, where near-zero values of the streamwise velocity are observed. Here, the air flow within the canyon almost separates from the above flow. Finally, the performance of the code for this scenario is around 1000 million nodes update per second. Time averaging started from time t 0 = 180 s (3 min) to t 1 = 1800 s (30 min). It took 1080 s (18 min) of execution time with the fine grid on the high-performance computing center, which is less than the time difference t 1 − t 0 .

Digital urban city portion
As application, the authors choose the city quarter "Bayerischer Bahnhof" in Leipzig, Germany [56]. With the exception of one 48 m high building, all buildings' heights in this section of the city are in the range of 20 m. The three-dimensional urban data of the selected quarter is illustrated in Fig. 4 29%. The original satellite data of the city have a spatial resolution of Δx = 3 m. These data have been interpolated to Δx = 1.3 m to obtain a finer CFD grid with a total number of nodes set to 498 × 498 × 98. The flow penetration across trees and bushes is neglected. In the presented simulations, these objects hence behave as solid obstacles.
We are well aware that better approximations exist to represent trees and bushes. The porous media model has for example often been used [57,58]. This common correction involves the use of a sink term in the momentum equation, that is a function of the foliage density [59]. With a high leaf area density, the air speed right downstream of a tree is however significantly reduced [60] and eventually approaches that of solid body.
In the present scenarios, dry resuspension is considered and mostly occurs in spring or summer. At this time of the year, little precipitations are expected [61] and trees in urban environments exhibit a relatively dense foliage reducing solar radiation within the city [62]. In this context, it is fair to represent trees and bushes as solid obstacles. A no-slip boundary condition is enforced onto the ground, building and grass surfaces. On the two lateral domain side, a periodic boundary condition is used. The top boundary is set to free slip and a zero gradient is used for the outlet boundary condition. The inlet boundary is discussed below. Wind speed and direction are known to vary seasonally. Figure 5 shows the wind measurement statistics recorded at the Leipzig Institute for Meteorology weather station. This figure, colored with four wind speed classes, is obtained for the time period 2007-2017. Over the years, the south wind direction is predominant, also in the summer, where dry resuspension is most likely to occur. Hence, the south wind flow direction has been chosen and implemented as inlet boundary condition. The inlet velocity profile follows the power law in Eq. (14) with = 0.18, z h = 60 m. The free stream wind speed is varied depending on the scenarios later described.

Turbulent flow field
To quantitatively appreciate the simulated flow field in the chosen city quarter, the velocity magnitude along the vertical line passing through point A (see Fig. 4), is recorded and averaged over time. The profiles are obtained at three wind speeds, that are u w = 7, 10 and 14 m/s, hereafter labeled as moderate, high, and very-high wind speeds, respectively. The term very-high is here employed as a reference to the highest wind speed used in this study. Actual wind speeds, especially those in the upper atmospheric region, are often greater [64]. Table 1 lists the free stream wind speeds u w presently used, along with the equivalent inlet velocities at pedestrian level u x (z=1m) calculated using Eq. (14).  The point A is close to the center of the domain and also not located in or on a building. The profile of the time-averaged wind velocity is shown in Fig. 6 and is compared to the power-law inlet function. The numerically determined profile matches well the expectation despite some differences close to the ground surface, which can be in part attributed to the grid resolution. Furthermore, the wind velocity along the vertical line passing through point A is affected by the upstream buildings which tend to slow down the velocity with respect to the undisturbed profile.
In a second stage, we discuss the turbulent flow fields obtained with u w = 7, 10 and 14 m/s. A snapshot of the instantaneous velocity magnitude in the horizontal plane at one meter above the ground is shown in Fig. 7 (top left and right). These instantaneous velocity fields are taken at time t = 120 s after the initial state. The simulation only took 130 s and hence is almost real-time. The overall flow fields obtained at u w = 7 and 10 m/s look somewhat similar. More turbulent structures are however observed for u w = 10 m/s. The strong ground-level wind flows, that are shown at one meter above the ground in Fig. 7, may endanger the safety of pedestrians, cause discomfort and trigger high particle resuspension. With a ground wind speed ranging from 9 to 13 m/s, a loss of balance when walking or cycling may occur [65]. In a second series of simulations, vegetation is removed. Figure 7 (see bottom left and bottom right) shows that trees considerably reduce the strong groundlevel winds. Inline with other works [66], green infrastructure is generally in favor of public health. Tree placement should however be carefully planned. Inside a street canyon, the  [63] Wind velocity (m/s). presence of tree leads to smaller flow velocities [67], which in-turns reduces the resuspension rate of large-inertia particles. It however also means a modified air exchange with the upper roof flows. Consequently, traffic exhaust emissions may adversely concentrate in a street canyon. We deliberately refrain from further discussing the adverse effect of strong ground-level winds on pedestrian comfort because the primary focus of this work is on particle resuspension and their subsequent dispersion. A peculiarity of the chosen city quarter lies in the frequent occurrence of building structures having a closed circumference and a small yard in the center (Fig. 4). Such "O"like urban structures, i.e. perimeter blocks, cause the formation of air recirculation zones. Oke et al. [11] reported that, in a two-dimensional system, skimming flow occurs whenever H∕W > 0.7 , where H is the building height and W the yard length. In such cases, a stable recirculation zone develops in the yard while the upper flow appears decoupled. Figure 8 shows the flow streamlines for a selected closed building structure with H∕W > 0.85 . A stable re-circulation zone, where the air is trapped inside, can be clearly observed. Even though it is not shown here, this skimming flow was also observed for other closed building structures exhibiting a greater ratio H/W. Under such circumstances, resuspended pollutants could not easily escape, potentially resulting in poor air quality [15].

Particle resuspension
The numerical simulation of the resuspension and subsequent dispersion of the particles for different magnitudes of the wind speed are here discussed. To investigate the dispersion of resuspending particles, four sets of simulations are carried out. In the first set S1, the cloud is composed of 40,000 poly-disperse particles randomly deposited on all ground surfaces. Using a larger number of particles was not found to significantly improve the statistical data presented here. The polydispersity is achieved by using a cloud composed of five particle classes, namely 60 μ m microplastics, 20 μ m metal particles, 80μ m, 90μ m and 100μ m pollen particles. In each class, 8000 particles are used. Size, density, and critical velocity u * c associated with each particle class are given in Table 2. The critical velocity u * c associated with each particle class is function of the substrate. These critical velocities are chosen to quantitatively match the values cited in [38,42,43]. Exact critical velocities for particles on asphalt or vegetated surface are not always available in the literature. We try to remain within the accepted critical velocity range. In real scenarios, deposited micro-plastics would be more concentrated on asphalt surfaces, heavy Fig. 8 Recirculation zone in a perimeter block. The top left schematic is taken from [15] metal particles near construction sites and pollen particles by trees or bushes. Yet, for simplicity and also lack of reported data, the particles are assumed evenly distributed on all ground surfaces present in the digital city model. In this first set S1, the trajectories of the mixed particles resuspending from the ground are investigated at three different wind speeds (Table 1). Simulations performed with a wind speed below 7 m/s did not show any significant resuspension. They are hence not shown here. Figure 9 illustrates the particle trajectories, which are plotted within 3 meters from the ground, where resuspension poses the greatest risk to human health. The particle trajectories are obtained from the simulation time t 0 = 0 s to t 1 = 120 s.
The particles traveling beyond this height are here not considered. In our simulations, the vast majority of resuspended particles did not redeposit and only a small fraction, typically less than 10% were found to hit buildings and trees within 3 m from the ground. As seen in Fig. 9a-c (top), only a small portion of the particles are brought into motion at moderate wind speed, that is u w = 7 m/s. By increasing the wind speed to u w = 10 m/s some high-risk areas with high particle concentrations are observed. Sensitive people should avoid these city areas even under moderate wind conditions. At very-high wind speed, that is u w = 14 m/s, almost all streets turn red and it is highly recommended not to walk or stay outside under such wind conditions. The equivalent wind inlet velocity perceived at pedestrian level (Table 1) is consistent with in-situ measurements on particle resuspension in London [68]. There, the authors reported an increase in the resuspension emission factor for ground wind speeds ranging from 1 to 8 m/s. Vegetation, beside changing the flow field, also reduces the resuspension. To quantify the overall effect of green infrastructure on particle resuspension, we removed all trees and replaced green surfaces with asphalt. Figure 9d, e (bottom) shows the resulting particle trajectories at moderate, high and veryhigh wind speeds. By comparing to reference Fig. 9a-c (top), the general absence of vegetation leads to a dramatic increase in particle resuspension. The ratio of high-risk area to total domain area is calculated using conventional image processing tools. It is the number of red pixels divided by the total number of pixels. The respective ratios of high-risk to total area for u w = 7, 10 and 14 m/s are given in Table 3.
At moderate high wind speed, that is u w = 7 m/s, the high-risk area is low irrespective of the green infrastructure. The ratio of high-risk area to total domain area is typically Table 2 Particle properties and critical friction velocities used in simulation sets S1-4 The cloud in S1 is composed of five particle classes evenly mixed, hence the 20% cloud fraction associated with each particle class lower than 0.2%. At high wind speed, that is u w = 10 m/s, the high-risk area almost doubles after removing the vegetation. The ratio of high-risk area to total domain area increases from 1.1 to 2.3%. At very-high wind speed and without green infrastructure, the number of red pixels in the Fig. 9f is significant and the ratio of high-risk area to total domain area reaches 8.1%. In the presence of green infrastructure, the latter drops to 3.1%. In an attempt to correlate the occurrence of high-risk areas with flow dynamics, we further investigate the scenario associated with Fig. 9e, that is the one without green infrastructure obtained for u w = 10 m/s. In Fig. 10, the particle trajectories within three meters above the ground (d) are laid side by side with the time-averaged friction velocity ũ * calculated on the ground surfaces (a), its instantaneous counterpart u * (b) and the time-averaged wind magnitude ‖ũ‖ at one meter above the ground (c). Figure 10d, in which some selected high-risk areas are   There are more of such outbreaks in the figure. In these high-risk city areas, the local wind orientation is relatively parallel to the street axis. Upstream of a narrow street, the wind within the canopy first strengthens near the entrance and then squeezes to eventually trigger high particle resuspension. While channeling, the wind accelerates carrying the resuspended particles upwards. In a perimeter block, on the contrary, the recirculation vortex associated with the urban canyon similar to that shown in Fig. 8, is too weak to trigger resuspension. In such blocks, a stay is relatively safe. The low resuspension in a perimeter block is in fact well reflected in Fig. 10a, b, where the friction velocity u * is generally low, typically well below 0.1 m/s. We recall, that resuspension of sub-millimeter particles from vegetated and relatively smooth surfaces occurs for critical friction velocities in the range (c) Time-average wind magnitude ũ one meter above the ground.
(d) Particle trajectories within 3 m from the ground. , the time-averaged velocity magnitude 1 m above the ground (c) and the particle trajectories within 3 m above the ground (d).
All data are taken from the simulation set S1 and shown for u w = 10 m/s without green infrastructure. The circles with dark blue contours highlight high-risk areas. Unobstructed open space (large circle, top left) and channeling areas (smaller circles) are propitious for high resuspension. In perimeter blocks, almost no resuspension occurs 0.2 < u * c < 1.5 m/s. In these perimeter blocks, the time-averaged wind velocity magnitude at one meter above the ground is also low, typically about 1 m/s or less. In high-risk areas, the wind velocity magnitude at one meter above the ground often exceeds 8 m/s and the friction velocity fluctuates significantly, suggesting that the near-wall turbulence intensity is also key to resuspension. Mukai et al. [69] investigated the effect of turbulence intensity on the resuspension of particles with diameters ranging from 13 to 20 μ m from a synthetic carpet, somewhat reminiscent of a short grassed surface. By using a flow straightener, the mean air velocity in the channel was kept constant. In their study, the number of detached particles was up by 30% after increasing the turbulence intensity from 9 to 30%.
In the simulation sets S2, S3 and S4, we remain more generic with respect to the particle material, size and density. To do so, the resuspension of three mono-disperse particle clouds, namely 30 μ m microplastics (set S2), 60 μ m heavy metal particles (S3) and 80 μ m pollen particles (S4), is investigated in a separate manner. In each simulation, only one cloud is used. Figure 11 shows the trajectories of the resuspended particles for the three different particle clouds at high wind speed, that is u w = 10 m/s. The ratio of high-risk area to total domain area for the three wind speeds and with green infrastructure is reported in Table 3. At moderate wind speed, that is u w = 7 m/s, the ratio of high-risk area to total domain area for the 30 μ m microplastic is 1.1%. At the same speed wind but, this time with the 60 μ m heavy metal particles, the high-risk area increases to 1.6%. With the 80 μ m pollen particles, the ratio is 2.0%. The increase in high-risk area from S2, to S3, to S4 is not significant and is mostly attributed to the increasing particle size. In fact, the particle density was found to play a negligible role. Similar conclusions can be drawn at high and very-high wind speeds. These numbers again show that, green infrastructures significantly reduce particle resuspension even at high and very-high wind speeds. Such findings are supported by other independent studies, see for instance Tiwari et al. [70], who reported that coniferous trees along traffic lanes are effectively reducing pollutant resuspension. Two tree chains on each side of a road results in the formation of a street canyon with smaller near-ground velocities [67], hence resuspension is reduced.

Conclusion
The contribution of our work is twofold. First, the coupling of a LBM solver with a resuspension model inspired from sediment resuspension in riverine flows, has been performed. Second, we have applied the model to simulate the real-time particle resuspension in a portion of the Leipzig city. So far, the state of the art has focused on mostly LES of pollutant dispersion and, very recently, on real-time simulation of turbulent urban flows. The wall-pollutant interaction in a city parts has generally received scarce attention in terms of modeling. Our work alleviates this shortcoming by providing an attempt at simulating the spread of resuspending pollutant particles with relatively high inertia. Here, the resuspension and dispersion of particles in the size range 10 to 100 μ m has been investigated. Resuspension was modeled using the concept of critical friction velocity. We have assumed a logarithmic wind profile anywhere in the city to link the particle and city scales. A Lagrangian model has been implemented to track the resuspended particles. For the fluid solver, which takes most of the computational cost, the cumulant Lattice Boltzmann Method has been used. The particle trajectories were analyzed within 3 m from the urban ground surface. The results showed that particle resuspension starts becoming significant at high wind speed, that is u w > 10 m/s. At very-high wind speed, that is beyond 14 m/s, and in the absence of green infrastructure, resuspension occurs in nearly all city parts. Three types of particles have been considered, namely microplastics, heavy-metals and pollen. Particle size is the most important criteria. Hence, strong winds tend to trigger the resuspension of many particle types, possibly resulting in carcinogenic and cardiovascular diseases, asthma and skin allergies. Another important aspect is the ambient wind direction, which significantly affects flow and spread of micro-pollutants in an urban environment. Taking advantage of graphics processing unit, simulations at different wind speeds have been performed in almost real time, meaning that the execution time is almost equal to the time span set in the simulation. With the upcoming increase in frequency of extreme events, studies on particle resuspension by strong winds will gain momentum in the future. Further works will include the development of mobile applications with data transfer from high-performance computers.
Consent to publication All authors give their consent for publication of this work.

Ethics approval Not applicable.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.