CET-vol95


 

 DOI: 10.3303/CET2295001 

 
 

 
 
 

 
 

 
 

 
 
 

 
 
 

 
 
 

 

 
 
 
 
 
 

 
 
 

 
 
 

 
 
 

 
 
 

 
 
 

 
 

Paper Received: 15 April 2022; Revised: 26 May 2022; Accepted: 12 June 2022 
Please cite this article as: Tinarelli G.L., Sozzi R., Barbero D., 2022, Implementation of a Simplified Micromixing Model Inside a Lagrangian 
Particle Dispersion Code for the Estimation of Concentration Variances and Peaks, Chemical Engineering Transactions, 95, 1-6 
DOI:10.3303/CET2295001 

CHEMICAL ENGINEERING TRANSACTIONS 

VOL. 95, 2022 

A publication of 

The Italian Association 
of Chemical Engineering 
Online at www.cetjournal.it 

Guest Editors: Selena Sironi, Laura Capelli 

Copyright © 2022, AIDIC Servizi S.r.l. 
ISBN 978-88-95608-94-5; ISSN 2283-9216 

Implementation of a Simplified Micromixing Model Inside a 

Lagrangian Particle Dispersion Code for the Estimation of 

Concentration Variances and Peaks 

Gianni Luigi Tinarellia*, Roberto Sozzib, Daniela Barberoa

a ARIANET S.r.l. Via Benigno Crespi 57, Milan 20159, Italy 
b Independent Scientific Consultant, Italy 

g.tinarelli@aria-net.it

In the frame of the modelling simulation of odour nuisances, the estimation of concentration peaks, representing

values averaged over a relatively short time of the order of the interval between subsequent breathes, is of 

fundamental importance. Dispersion models currently used in this field cannot reconstruct this kind of values at 

relatively high frequency, due to their intrinsic theoretical design that allows to give time- or ensemble average- 

concentrations only. The scope of this work is to describe the implementation of a simplified micromixing model 

inside a standard ensemble average Lagrangian Particle Dispersion Model, with the aim of simulating the field 

of concentration variances together with concentration averages. A simplified micromixing model represents a 

way to describe the interaction between the part of the emitted plume and the rest of the atmospheric flow,

derived through bulk entrainment relationships. This simplified view allows the description of the first two

moments of the concentration distribution which is however sufficient to describe a peak-to-mean relationship 

making some hypotheses about the form of the distribution. Some preliminary results of the application of this

method inside the SPRAY Lagrangian Particle Dispersion Model are shown, comparing both the instantaneous

concentration and the peak-to-mean ratio together with their spatial behaviour derived in some controlled

conditions with those obtained from the application of other schemes currently included in the code. 

1. Introduction

Atmospheric dispersion models are widely used tools to give a spatial and temporal description of the 

concentrations of emitted substances that could potentially generate different problems to the population. In

most of their applications, typically related to the comparison with air quality regulatory rules, atmospheric

dispersion models are used to reconstruct mean concentrations representing for example the average over time

intervals of the order of 1 hour. There exist some specific applications, as the dispersion of explosive or

flammable, toxic or odorous emissions, where the simple reconstruction of the mean value is not sufficient due

to the intrinsic characteristics of the consequences generated by the involved substances. In particular, the 

sensation of olfactory nuisance occurs in an individual during the normal respiratory activity happening at a 

relatively high frequency, typically every less than 10 seconds. This implies that the reconstruction of the impact 

generated by the odorous emissions should in principle take into account this peculiarity and a reconstruction 

of the short-term concentrations, which are more assimilated to instantaneous events, is potentially required 

(Capelli et al., 2013). This means that models should include a way to predict also concentration fluctuations or

a statistical indicator to describe them. The availability of both the average concentrations and concentration 

fluctuations as model output, allows to make a more complete and correct estimation of the impact that odorous 

emissions can have on the population. Not many of the typical models used in the framework of odour impact 

have a direct mechanism inside to fulfil this request. A useful review of the efforts made by the scientific

community in the last seventy years to include the reconstruction of concentration fluctuations inside dispersion

models can be found in Cassiani et al. (2020). A particular attention was devoted to Lagrangian one-particle

dispersion models, derived from the fundamental work of Thomson (1987), and widely used in the frame of air 

quality to reconstruct ensemble average concentrations over complex terrain. Different methods have been

1



proposed for these models to couple the calculation of concentration fluctuations with the standard calculation 

of average concentrations. Among them, a parameterization of the concentration variance transport and 

dissipation equation (Ferrero et al., 2017) and a simplified form of a micromixing model allowing to take into 

account only particles coming from the emissions (Cassiani, 2013) can be cited. Both these methods have been 

incorporated into the Lagrangian Particle Dispersion Model SPRAY (Tinarelli et al., 2000), in order to calculate 

concentration fluctuations close to the ground, allowing to estimate a peak-to-mean concentration ratio to be 

applied in simulations of odorous emissions. In this paper, the implementation of the simplified micromixing 

method, named Volume Particle Approach (VPA), is explained and some results in controlled conditions are 

shown, also comparing this method with the Variance Transport Equation (VTE) method already present in the 

code. Being the implementation mainly targeted to the reconstruction of odorous emissions, results from the 

two methods are also compared with a more trivial peak-to-mean parameterization widely used in this field, to 

understand the main differences that can arise. 

2. Theoretical description of the SPRAY model 

In the Lagrangian Particle SPRAY code, the dispersion of an airborne pollutant is simulated following the motion 

of a large number of fictitious particles. Each particle moves through a ‘‘transport’’ component of the velocity 

provided by an external meteorological driver and a stochastic or ‘‘turbulent’’ component obtained by solving a 

3-D form of a Langevin equation defined in its general form, following Thomson (1987) as: 

𝑑𝑢𝑖 = 𝑎𝑖 (𝒙, 𝒖, 𝑡)𝑑𝑡 + 𝑏𝑖𝑗 (𝒙, 𝒖, 𝑡)𝑑𝑊𝑗 (𝑡) 

where x represents the vector of the particle position, and u the Lagrangian velocity vector,  
𝑎𝑖 (𝒙, 𝒖, 𝑡)𝑑𝑡 is a deterministic term depending on Eulerian probability density function (PDF) of the turbulent 
velocity, determined from the Fokker–Planck equation and 
𝑏𝑖𝑗 (𝒙, 𝒖, 𝑡)𝑑𝑊𝑗 (𝑡) represents a stochastic term where 𝑑𝑊𝑗 (𝑡) is the incremental Wiener process. Each particle 

carries a mass of pollutants that in case of odorous emission can be considered in terms of “Odour Units” and 

ensemble average concentrations can be computed considering the total mass of all the particles contained into 

control volumes of a 3D grid. Although Lagrangian Particle Dispersion Models are used in many situations to 

describe the dispersion in complex conditions, they cannot in principle directly address the calculation of 

concentration fluctuations useful for a correct estimation of instantaneous concentrations.  

Currently, the simplified VTE method to compute concentration variances is implemented into the SPRAY code 

as follows. An instantaneous concentration c can be supposed as the sum of an (time or ensemble) average 

value C plus a fluctuation c’. In general, the conservation equation for the concentration variance c’2 can be 

written (Stull, 1989, Milliez and Carissimo, 2008) in the following form: 

𝜕𝑐′2̅̅ ̅̅

𝜕𝑡
+ 𝑈𝑗

𝜕𝑐′2̅̅ ̅̅

𝜕𝑥𝑗
= −2𝜈𝐶 (

𝜕𝑐′

𝜕𝑥𝑗
)

2̅̅ ̅̅ ̅̅ ̅̅ ̅

−
𝜕

𝜕𝑥𝑗
(𝑢𝑗

′𝑐′2̅̅ ̅̅ ̅̅ ̅) − 2𝑢𝑗
′𝑐′̅̅ ̅̅ ̅ (

𝜕𝐶

𝜕𝑥𝑗
)

2

 

where the two terms on the left represent the storage and advection, the first term on the right represents the 

viscous dissipation, the second one represents the divergence of the turbulent flow and the third one the source 

term for the variance, proportional to the average concentration gradient. Following Hsieh (2007) and Milliez 

and Carissimo (2008) an algebraic closure for the molecular dissipation can be rewritten as: 

−2𝜈𝐶 (
𝜕𝑐′

𝜕𝑥𝑗
)

2

= 𝜀𝑐 =
𝑐′2̅̅ ̅̅

𝑡𝑑
 

where td represents a dissipation time scale characteristic of the concentration variance decay. For the third 

order term representing the turbulent divergence of the concentration variance a K closure can be adopted, 

obtaining: 

 𝑢𝑗
′𝑐′̅̅ ̅̅ ̅ = −𝐾𝑗

𝜕𝐶

𝜕𝑥𝑗
           where: 

𝐾𝑗 = 𝜎𝑢𝑗
2 𝑇𝐿𝑗 

Oettl and Ferrero (2017) evaluated the relative importance of the various terms inside the conservation equation 

for the concentration variance through some considerations derived from available measurements. As a result, 

they estimated that the advection term on the left and the diffusion term on the right can be reasonably 

neglected, obtaining an equation independent from the mean flow field in the form: 

2



𝜕𝑐′2̅̅ ̅̅

𝜕𝑡
= 2𝜎𝑢𝑗

2 𝑇𝐿𝑗 (
𝜕𝐶

𝜕𝑥𝑗
)

2

−
𝑐′2̅̅ ̅̅

𝑡𝑑
 

This equation has the following analytical solution: 

𝑐′2̅̅ ̅̅ = 2𝜎𝑢𝑗
2 𝑇𝐿𝑗 𝑡𝑑 (

𝜕𝐶

𝜕𝑥𝑗
)

2

[1 − 𝑒 (−𝑡 𝑡𝑑⁄ )] 

It is suggested (Ferrero et al., 2020) that the approximation represented by the last relationship can be 

considered acceptable for the calculation of the concentration variance, once the turbulent variables and the 

average concentration field are known. Using the average concentration field directly computed by the code and 

the turbulent variables considered for the particle movement, a time varying concentration variance field is 

calculated with the previous relationship on the same grid considered for the average concentrations. It is 

supposed that td = TLw. Once the concentration variance is computed, the peak concentration is then estimated 

considering the 98 percentile of a two parameters Weibull distribution having the known average concentration 

and concentration variance. 

2.1 Simplified Micromixing VPA model 

A micromixing model (or PDF model) simulates the intrinsically continuous Planetary Boundary Layer as a 

geometric space in which a very large number of air particles, each fully detectable, are uniformly distributed. 

Each of them is completely characterized at a generic instant t by:    

• a position in space X(t), 

• a velocity fluctuation u(X,t) with respect to a mean (Eulerian) field of motion U(X,t), 

• a concentration of the interested pollutant  C(X,t). 

At each instant t prior to an initial instant t0 all particles (initially uniformly distributed in space) possess a 

concentration C(t<t0) = 0. At a given initial time t0, some of these particles will transit through the source (e.g. a 

point source) and will acquire mass from it and therefore an initial concentration C0(t0), while all the others will 

continue to keep zero concentration.  From the instant t0 onwards, the model will begin to simulate the dispersion 

of all the particles (both those with non-zero concentration and those with zero concentration), considering their 

different stochastic trajectories and their mutual interaction. This interaction is constituted by a mass exchange 

of the pollutant between a generic particle and the adjacent particles, induced by molecular diffusivity and driven 

by the turbulence that characterizes locally the PBL. Therefore, the model will simulate the trajectory of all 

particles as a normal Lagrangian one-particle model with the additional practical problem related to the huge 

number of particles whose trajectory and mass exchange must be simulated.  

In a micromixing model, the pollutant exchange among close particles is modelled, for simplicity, through a bulk 

law describing, for each particle, the mass exchange with the adjacent external environment, simulated as a 

continuous fluid characterized by an average concentration C. In practice, for the p-th particle this exchange is 

described by the following micromixing relation that simulates the action of molecular diffusivity: 

𝑑𝐶𝑝(𝑿𝑝, 𝑡)

𝑑𝑡
= −

𝐶𝑝(𝑿𝑝, 𝑡) − 𝐶(𝑿𝑝, 𝑡)
̅̅ ̅̅ ̅̅ ̅̅ ̅̅ ̅

𝜏𝑚
 (1) 

where 𝜏𝑚 represents the so-called micromixing time scale. Some particles will decrease in concentration (those 

passing through the source) while others (those that constitute the surrounding air) will increase it. 

Assuming to divide the entire computational domain in cells, at the end of the time step there will be Nk particles 

in the k-th cell, each with its own concentration. The average characteristic concentration of the cell can be 

computed as: 

 𝐶𝑘̅̅ ̅(𝑡) =
1

𝑁𝑘
∑ 𝐶𝑝(𝑡)

𝑁𝑘
𝑝_1                                          

while the second moment can be computed as: 

𝐶𝑘
2̅̅̅̅ (𝑡) =

1

𝑁𝑘
∑ (𝐶𝑝(𝑡))

2
𝑁𝑘

𝑝_1
 

and the concentration variance can be computed as: 

(𝜎𝐶
2)𝑘 = 𝐶𝑘

2̅̅̅̅ (𝑡) − (𝐶𝑘̅̅ ̅(𝑡))
2
 (2) 

3



To overcome the problem of simulating a huge number of particles, both those emitted by the considered 

sources and those referring to the background, Cassiani (2013) proposed the VPA approach. This approach, 

allows to consider only the particles emitted by the source, as during the normal execution of the SPRAY code. 

To this aim, the volume Vp is added to each particle as an additional property. Considering the appropriate initial 

conditions for the particle volume, related mainly to the initial source size, it is possible to define the evolution 

of the particle concentration from its initial value, as described in equation (1), during a model discretization time 

step. Since the particle mass must be conserved, this means that also a change in volume must be considered 

as: 

𝑉𝑝(𝑡 + ∆𝑡) = 𝑉𝑝(𝑡) (
𝐶𝑝 (𝑡)

𝐶𝑃(𝑡) + ∆𝐶𝑝
) 

In the current implementation of the VPA model inside SPRAY, the micromixing time scale 𝜏𝑚 is defined, as 

described in Dixon and Tomlin (2007), with the following simplified formula: 

𝜏𝑚 = 0.75 𝑘 𝜀⁄  

where k represents the turbulent kinetic energy and 𝜀 represents its dissipation rate. Both terms can be 
calculated in each point using the variables that define the turbulence parameters inside the code. Once the 

average concentrations are computed and concentration variances are estimated using (2), then the peak 

concentration is calculated as the 98 percentile of a two parameters Weibull distribution. 

3. Results 

3.1 Peak-to-mean in controlled conditions 

To obtain some preliminary results, a simulation in controlled conditions has been performed in order to compare 

the results between the VTE and VPA methods. A 12 x 12 km2 computational domain with flat terrain and a 

surface roughness of 0.1 m has been considered. The meteorological conditions consider a constant mean wind 

of 1 m/s blowing towards North and unstable convective stability conditions characterized by a Monin-Obuhkov 

length of -10 m, a vertical mixed layer of 1500 m and a convective velocity scale of 0.2 m/s. In these conditions, 

a neutrally buoyant point source with stack height = 10 m, stack diameter = 0.5 m, stationary emitting flow = 

1700 eu/s has been simulated, where eu represents generic “emission units” that can be interpreted as odour 

units in case of odour emissions or mass units in case of pollutant mass emissions.   Figure 1a represents the 

average concentration field close to the ground given by the standard dispersion model, while Figure 1b and 

Figure 1c reproduce peak concentrations using the VTE method and the VPA method respectively, as described 

in Chapter 2. For both methods, peak concentration fields show a pattern evidently larger than the one of 

average values, with values systematically higher in each simulated point. The VTE method shows larger 

concentration values in the proximity to the source with a maximum of 5 eu/m3, compared to a maximum of the 

average value at the same point of 1.57 eu/m3 (corresponding to a peak-to-mean value for the maximum 

average concentration of 3.18). The VPA method shows instead a larger pattern extending to a greater distance 

from the source, even if characterized by smaller concentrations close to the emission with a maximum of 1.96 

eu/m3 (corresponding to a peak-to mean value for the maximum average concentrations of 1.25 only). 

 

 
Figure 1: ground level average concentrations (a), peak concentrations using VTE method (b) and peak 

concentrations using VPA method (c) 

 

To have a clearer view of the spatial behaviour of the peak values and their relationship with the average 

concentrations, two sections along the x axis of the ground level average and peak concentrations have been 

a b

 

c 

4



extracted corresponding to the two horizontal dashed and continuous lines depicted in Figure 1 a,b,c, located 1 

and 5 km downwind the source respectively. Figures 2a and 2b illustrate the behaviour of average and peak 

concentrations at the 1 km and 5 km horizontal sections respectively. Along section 1, closer to the source, VPA 

method shows similar values with respect to VTE with slightly larger values close to the plume centreline where 

a peak-to-mean ratio is about to 1.75 for VPA and 1.55 for VTE. Moving away from the source position, the VPA 

method shows a clearer signal, with larger peak values compared to VTE methods and a peak-to-mean ratio 

close to the plume centreline of about 2.4 for VPA and 1.54 for VTE. 

 

 
Figure 2: ground level average and peak concentrations at the two sections along x axis 1 km (a) and 5 km (b) 

downwind the source position. 

 

A preliminary comparison of the peak-to-mean ratios obtained by the VPA model against crosswind theoretical 

values 250 m downwind the source position obtained by a semi-empirical Gaussian Plume model validated with 

experiments by Lǿfstrom at al. (1994) in unstable conditions is illustrated in Figure 3, showing promising results. 

 

 
Figure 3: comparison of crosswind peak-to-mean ratios 250 m downwind source in unstable conditions obtained 

by the semi-empirical Gaussian model validated in Lofstrom et al. (1994) and results from the VPA model 

3.2 Comparison with constant peak-to-mean ratio 

It is a common procedure, in the frame of odour applications, to calculate peak concentrations by applying a 

simple peak-to-mean ratio, kept uniformly constant throughout the domain and over time. This way of operating 

is meant to represent a sort of precautional procedure, giving safe peak values everywhere. As an example, 

although there is not a legislation at national level, some Italian guidelines are applied at a regional level, 

suggesting the use of a constant 2.3 peak-to-mean ratio. It is hence interesting to understand the differences 

that can arise using the two VTE and VPA methods introduced here and a simplified procedure in the calculation 

of peak concentrations. Figures 4a and 4b illustrate the absolute differences between peak values computed 

using VTE and VPA methods respectively and those computed applying a 2.3 uniform peak-to-mean ratio. In 

the figure, regions represented by colours from green to red are those containing a peak-to-mean ratio larger 

than 2.3, while regions in blue are those with peak-to-mean smaller. Both VTE and VPA methods calculate a 

non-uniform distribution of the peak-to-mean ratio both with a general similar behaviour showing smaller than 

2.3 values downwind the source close to the plume centreline and larger than 2.3 values close to the plume 

external border. The VTE method gives larger peak concentrations close to the source, where the VPA method 

doesn’t calculate such big values. The VPA method shows instead peak-to-mean ratios significantly higher than 

the constant 2.3 factor also at large distances from the source, where the VTE method doesn’t manifest such 

big values. 

a b 

5



 
Figure 4: absolute differences between peak concentrations at ground computed using the VTE method (a) and 

VPA method (b) and the application of a 2.3 constant peak-to-mean ratio. 

4. Conclusions 

A simplified micromixing model (VPA) was introduced inside the Lagrangian Particle Dispersion Model SPRAY 

with the aim of calculating concentration variances and estimating peak concentrations, specifically required to 

reconstruct the impact of odorous emissions. This method allows to take into consideration only the particles 

related to the emissions, increasing the computational efficiency with respect to a standard micromixing method 

that requires to follow also particles representing the atmospheric background. The micromixing model has been 

compared in a controlled condition with another existing method already present in the code implementing a 

simplified form for the transport equation of the concentration variance (VTE). Both methods show a general 

behaviour of giving non uniform peak-to-mean ratios with larger peak concentrations close to the plume border 

and smaller around the plume centreline. Results show that the VTE method determines higher peaks close to 

the source, with peak-to-mean ratios that can be significantly larger than a factor 3, while the VPA method shows 

smaller peak-to-mean ratios close to the source but extending over regions further away from the source, even 

for more than 5 km. These preliminary results must be confirmed through a comparison with available measured 

data in a real environment and this validation procedure will be implemented in the next future.  

References 

Capelli L., Sironi S., Del Rosso R., Guillot J.M., 2013, Measuring odours in the environment vs. dispersion 

modelling: a review - Atmospheric Environment, 79, 731–743. 

Cassiani M., 2013, The Volumetric Particle Approach for Concentration Fluctuations and Chemical Reactions 

in Lagrangian Particle and Particle-grid Models - Boundary-Layer Meteorology 146, 207–233. 

Cassiani M., Bertagni M.B., Marro M., Salizzoni P., 2020, Concentration Fluctuations from Localized 

Atmospheric Releases - Boundary-Layer Meteorology, 177, 461–510. 

Dixon N.S., Tomlin A.S., 2007, A Lagrangian stochastic model for predicting concentration fluctuations in urban 

areas - Atmospheric Environment, 41, 8114-8127. 

Ferrero E., Mortarini L., Purghè F., 2017, A Simple Parametrization for the Concentration Variance Dissipation 

in a Lagrangian Single-Particle Model - Boundary-Layer Meteorology, 163, 91–101. 

Ferrero E., Monor A., Mortarini L., Oettl D., 2020, Concentration fluctuations and odor dispersion in Lagrangian 

Models – Atmosphere, 11, 27; doi:10.3390/atmos11010027. 

Hsieh K.K., Lien F.S., Yee E., 2007, Numerical modelling of passive scalar dispersion in an urban canopy layer 

- J. Wind Eng. Ind. Aerodyn., 95, 1611-1636. 

Lǿfstrom P., Jǿrgensen H., Lyck E., Mikkelsen T.,1994, A concentration fluctuation model for decision-makers 

based on joint tracer and lidar measurements from a non-buoyant elevated plume – Trans. Ecology and 

Environ., 3, 571-579. 

Milliez M., B. Carissimo (2008): Computational fluid dynamical modelling of concentration fluctuations in an 

idealized urban area – Boundary Layer Meteorology, 127, 241-259. 

Oettl D., Ferrero E., 2017, A simple model to assess odour hour for regulatory purposes - Atmospheric 

Environment, G155, 162-173. 

Thomson D., 1987, Criteria for the selection of stochastic models of particle trajectories in turbulent flows - J. 

Fluid Mech., 180, 529–556. 

Stull R., 1988, An introduction to Boundary Layer Meteorology - Kluwer Academic Publishers, Dordrecht. 

Tinarelli G., Anfossi D., Trini Castelli S., Bider M., Ferrero E., 2000, A new high performance version of the 

Lagrangian particle dispersion model SPRAY, some case studies, in Air pollution modeling and its 

application XIII, 499-507. 

a 
b 

6


	Implementation of a Simplified Micromixing Model Inside a Lagrangian Particle Dispersion Code for the Estimation of Concentration Variances and Peaks